#836536 (C++17) No.2993 冪乗乗 mod 冪乗

提出ソース

結果

問題	No.2993 冪乗乗 mod 冪乗
ユーザー	👑 p-adic
提出日時	2023-02-02 11:41:00
言語	C++17 (gcc 12.3.0 + boost 1.83.0)
結果	AC
実行時間	465 ms / 10,000 ms
コード長	12,449 bytes
コンパイル時間	3,899 ms
コンパイル使用メモリ	226,300 KB
実行使用メモリ	6,820 KB
最終ジャッジ日時	2024-12-19 16:19:30
合計ジャッジ時間	13,254 ms
ジャッジサーバーID （参考情報）	judge2 / judge1

このコードへのチャレンジ
（要ログイン）

テストケース

テストケース表示

入力	結果	実行時間実行使用メモリ
testcase_00	AC	277 ms 6,816 KB
testcase_01	AC	2 ms 6,816 KB
testcase_02	AC	2 ms 6,820 KB
testcase_03	AC	2 ms 6,820 KB
testcase_04	AC	2 ms 6,820 KB
testcase_05	AC	2 ms 6,816 KB
testcase_06	AC	2 ms 6,816 KB
testcase_07	AC	6 ms 6,816 KB
testcase_08	AC	3 ms 6,816 KB
testcase_09	AC	3 ms 6,816 KB
testcase_10	AC	3 ms 6,816 KB
testcase_11	AC	465 ms 6,820 KB
testcase_12	AC	327 ms 6,816 KB
testcase_13	AC	341 ms 6,816 KB
testcase_14	AC	236 ms 6,816 KB
testcase_15	AC	187 ms 6,816 KB
testcase_16	AC	416 ms 6,816 KB
testcase_17	AC	390 ms 6,816 KB
testcase_18	AC	428 ms 6,820 KB
testcase_19	AC	315 ms 6,816 KB
testcase_20	AC	370 ms 6,816 KB
testcase_21	AC	413 ms 6,816 KB
testcase_22	AC	312 ms 6,816 KB
testcase_23	AC	458 ms 6,816 KB
testcase_24	AC	317 ms 6,816 KB
testcase_25	AC	315 ms 6,816 KB
testcase_26	AC	411 ms 6,820 KB
testcase_27	AC	248 ms 6,820 KB
testcase_28	AC	177 ms 6,816 KB
testcase_29	AC	97 ms 6,820 KB
testcase_30	AC	101 ms 6,816 KB

権限があれば一括ダウンロードができます

ソースコード

raw source code

// 入力制約チェック
#pragma GCC optimize ( "O3" )
#pragma GCC optimize( "unroll-loops" )
#pragma GCC target ( "sse4.2,fma,avx2,popcnt,lzcnt,bmi2" )
#include <bits/stdc++.h>
using namespace std;

using uint = unsigned int;
using ll = long long;

#define MAIN main
#define TYPE_OF( VAR ) remove_const<remove_reference<decltype( VAR )>::type >::type
#define UNTIE ios_base::sync_with_stdio( false ); cin.tie( nullptr ) 
#define CEXPR( LL , BOUND , VALUE ) constexpr const LL BOUND = VALUE 
#define CIN( LL , A ) LL A; cin >> A 
#define ASSERT( A , MIN , MAX ) assert( ( MIN ) <= A && A <= ( MAX ) ) 
#define CIN_ASSERT( A , MIN , MAX ) CIN( TYPE_OF( MAX ) , A ); ASSERT( A , MIN , MAX ) 
#define FOR( VAR , INITIAL , FINAL_PLUS_ONE ) for( TYPE_OF( FINAL_PLUS_ONE ) VAR = INITIAL ; VAR < FINAL_PLUS_ONE ; VAR ++ ) 
#define FOREQ( VAR , INITIAL , FINAL ) for( TYPE_OF( FINAL ) VAR = INITIAL ; VAR <= FINAL ; VAR ++ ) 
#define REPEAT( HOW_MANY_TIMES ) FOR( VARIABLE_FOR_REPEAT , 0 , HOW_MANY_TIMES ) 
#define QUIT return 0 
#define COUT( ANSWER ) cout << ( ANSWER ) << "\n"; 
#define RETURN( ANSWER ) COUT( ANSWER ); QUIT 

template <typename INT1 , typename INT2> inline constexpr INT1& Residue( INT1& n , const INT2& M ) noexcept { return n < 0 ? ( ( ( ( ++n ) *= -1 ) %= M ) *= -1 ) += M - 1 : n %= M; }
template <typename INT1 , typename INT2> inline constexpr INT1 Residue( INT1&& n , const INT2& M ) noexcept { return move( Residue( n , M ) ); }

#define POWER( ANSWER , ARGUMENT , EXPONENT )				\
  static_assert( ! is_same<TYPE_OF( ARGUMENT ),int>::value && ! is_same<TYPE_OF( ARGUMENT ),uint>::value ); \
  TYPE_OF( ARGUMENT ) ANSWER{ 1 };					\
  {									\
    TYPE_OF( ARGUMENT ) ARGUMENT_FOR_SQUARE_FOR_POWER = ( ARGUMENT );	\
    TYPE_OF( EXPONENT ) EXPONENT_FOR_SQUARE_FOR_POWER = ( EXPONENT );	\
    while( EXPONENT_FOR_SQUARE_FOR_POWER != 0 ){			\
      if( EXPONENT_FOR_SQUARE_FOR_POWER % 2 == 1 ){			\
	ANSWER *= ARGUMENT_FOR_SQUARE_FOR_POWER;			\
      }									\
      ARGUMENT_FOR_SQUARE_FOR_POWER *= ARGUMENT_FOR_SQUARE_FOR_POWER;	\
      EXPONENT_FOR_SQUARE_FOR_POWER /= 2;				\
    }									\
  }									\

#define POWER_MOD( ANSWER , ARGUMENT , EXPONENT , MODULO )		\
  ll ANSWER{ 1 };							\
  {									\
    ll ARGUMENT_FOR_SQUARE_FOR_POWER = ( ( ARGUMENT ) % ( MODULO ) ) % ( MODULO ); \
    ARGUMENT_FOR_SQUARE_FOR_POWER < 0 ? ARGUMENT_FOR_SQUARE_FOR_POWER += ( MODULO ) : ARGUMENT_FOR_SQUARE_FOR_POWER; \
    TYPE_OF( EXPONENT ) EXPONENT_FOR_SQUARE_FOR_POWER = ( EXPONENT );	\
    while( EXPONENT_FOR_SQUARE_FOR_POWER != 0 ){			\
      if( EXPONENT_FOR_SQUARE_FOR_POWER % 2 == 1 ){			\
	ANSWER = ( ANSWER * ARGUMENT_FOR_SQUARE_FOR_POWER ) % ( MODULO ); \
      }									\
      ARGUMENT_FOR_SQUARE_FOR_POWER = ( ARGUMENT_FOR_SQUARE_FOR_POWER * ARGUMENT_FOR_SQUARE_FOR_POWER ) % ( MODULO ); \
      EXPONENT_FOR_SQUARE_FOR_POWER /= 2;				\
    }									\
  }									\

// 通常の二分探索（単調関数-目的値が一意実数解を持つ場合にそれを超えない最大の整数を返す）
#define BS( ANSWER , MINIMUM , MAXIMUM , EXPRESSION , TARGET )		\
  ll ANSWER = MAXIMUM;							\
  {									\
    ll VARIABLE_FOR_BINARY_SEARCH_L = MINIMUM;				\
    ll VARIABLE_FOR_BINARY_SEARCH_U = ANSWER;				\
    ll VARIABLE_FOR_DIFFERENCE_FOR_BINARY_SEARCH = ( TARGET ) - ( EXPRESSION ); \
    if( VARIABLE_FOR_DIFFERENCE_FOR_BINARY_SEARCH == 0 ){		\
      VARIABLE_FOR_BINARY_SEARCH_L = ANSWER;				\
    } else {								\
      ANSWER = ( VARIABLE_FOR_BINARY_SEARCH_L + VARIABLE_FOR_BINARY_SEARCH_U ) / 2; \
    }									\
    while( VARIABLE_FOR_BINARY_SEARCH_L != ANSWER ){			\
      VARIABLE_FOR_DIFFERENCE_FOR_BINARY_SEARCH = ( TARGET ) - ( EXPRESSION ); \
      if( VARIABLE_FOR_DIFFERENCE_FOR_BINARY_SEARCH == 0 ){		\
	break;								\
      } else {								\
	if( VARIABLE_FOR_DIFFERENCE_FOR_BINARY_SEARCH > 0 ){		\
	  VARIABLE_FOR_BINARY_SEARCH_L = ANSWER;			\
	} else {							\
	  VARIABLE_FOR_BINARY_SEARCH_U = ANSWER;			\
	}								\
	ANSWER = ( VARIABLE_FOR_BINARY_SEARCH_L + VARIABLE_FOR_BINARY_SEARCH_U ) / 2; \
      }									\
    }									\
  }									\

template <typename INT , INT val_limit  , int length_max>
class PrimeEnumeration
{

public:
  INT m_val[length_max];
  int m_length;
  inline constexpr PrimeEnumeration();

};

template <typename INT , INT val_limit , int length_max>
inline constexpr PrimeEnumeration<INT,val_limit,length_max>::PrimeEnumeration() : m_val() , m_length( 0 )
{

  bool is_comp[val_limit] = {};

  for( INT i = 2 ; i < val_limit ; i++ ){

    if( is_comp[i] == false ){

      INT j = i;

      while( ( j += i ) < val_limit ){

	is_comp[j] = true;

      }

      m_val[m_length++] = i;

    }

  }

}

template <typename INT> inline int log( const ll& b , INT n )
{
  int answer = -1;
  while( n > 0 ){
    n /= b;
    answer++;
  }
  return answer;
}

template <typename INT>
INT ChineseRemainderTheorem( const INT& b_0 , const INT& c_0 , const INT& b_1 , const INT& c_1 , ll& lcm )
{
  
  INT a[2][2];
  INT b[2] = { b_0 , b_1 };
  int i_0 = ( b_0 >= b_1 ? 0 : 1 );
  int i_1 = 1 - i_0;

  for( uint i = 0 ; i < 2 ; i++ ){

    INT ( & ai )[2] = a[i];

    for( uint j = 0 ; j < 2 ; j++ ){
      ai[j] = ( i == j ? 1 : 0 );
    }

  }
  
  INT q;

  while( b[i_1] != 0 ){

    INT& b_i_0 = b[i_0];
    INT& b_i_1 = b[i_1];
    INT ( &a_i_0 )[2] = a[i_0];
    INT ( &a_i_1 )[2] = a[i_1];
    q = b_i_0 / b_i_1;
    a_i_0[i_0] -= q * a_i_1[i_0];
    a_i_0[i_1] -= q * a_i_1[i_1];
    b_i_0 %= b_i_1;
    swap( i_0 , i_1 );

  }

  INT& gcd = b[i_0];
  INT c = c_0 % gcd;

  if( c_1 % gcd != c ){

    return -1;

  }
  
  lcm = ( b_0 / gcd ) * b_1;
  INT ( &a_i_0 )[2] = a[i_0];
  INT& a_i_00 = a_i_0[0];
  a_i_00 *= ( c_1 - c ) / gcd;
  Residue( a_i_00 , lcm );
  INT& a_i_01 = a_i_0[1];
  a_i_01 *= ( c_0 - c ) / gcd;
  a_i_01 = ( a_i_01 >= 0 ? a_i_01 % lcm : lcm - ( - a_i_01 - 1 ) % lcm - 1 );
  return ( c + a_i_00 * b_0 + a_i_01 * b_1 ) % lcm;

}

void ComputeLocalAnswer( int& factor_num , vector<ll>& factor , vector<ll>& factor_power , vector<ll>& local_answer , bool& unsolvable , int& B , const int& N , const ll& M , ll prime_curr )
{
  int prime_curr_minus = prime_curr - 1;
  int order = prime_curr_minus;
  FOR( j , 0 , factor_num ){
    const int& factor_j = factor[j];
    while( order % factor_j == 0 ){
      order /= factor_j;
    }
  }
  order = prime_curr_minus / order;
  factor_num++;
  factor.push_back( prime_curr );
  int exponent = 1;
  B /= prime_curr;
  while( B % prime_curr == 0 ){
    exponent++;
    B /= prime_curr;
  }
  exponent *= N;
  POWER( factor_power_curr , prime_curr , exponent );
  factor_power.push_back( factor_power_curr );
  int valuation = 0;
  int case_num;
  ll M_shifted;
  if( N == 0 ){
    POWER_MOD( M_power , M , order , prime_curr );
    M_shifted= Residue( 1 - M_power , prime_curr );
    case_num = 0;
  } else if( M == 1 ){
    M_shifted = 0;
    case_num = 0;
  } else if( order == 1 ){
    M_shifted = 1 - M;
    while( M_shifted % prime_curr == 0 ){
      M_shifted /= prime_curr;
      valuation++;
    }
    Residue( M_shifted , factor_power_curr );
    case_num = 1;
  } else if( order == 2 ){
    ll M_plus = 1 + M;
    while( M_plus % prime_curr == 0 ){
      M_plus /= prime_curr;
      valuation++;
    }
    M_plus %= factor_power_curr;
    M_shifted = 1 - M;
    while( M_shifted % prime_curr == 0 ){
      M_shifted /= prime_curr;
      valuation++;
    }
    ( Residue( M_shifted , factor_power_curr ) *= M_plus ) %= factor_power_curr;
    case_num = 1;
  } else if( prime_curr <= 11 ){
    int extra = exponent << 1;
    if( extra % prime_curr_minus == 0 ){
      extra /= prime_curr_minus;
    } else {
      ++( extra /= prime_curr_minus );
    }
    POWER( factor_power_extra , prime_curr , extra );
    factor_power_extra *= factor_power_curr;
    POWER_MOD( M_power , M , order , factor_power_extra );
    M_shifted = 1 - M_power;
    while( M_shifted % prime_curr == 0 && valuation < extra ){
      M_shifted /= prime_curr;
      valuation++;
    }
    Residue( M_shifted , factor_power_curr );
    case_num = valuation < extra ? 1 : 2;
  } else {
    POWER_MOD( M_power , M , order , factor_power_curr );
    M_shifted = Residue( 1 - M_power , factor_power_curr );
    case_num = 2;
  }
  if( valuation == 0 && M_shifted % prime_curr != 0 ){
    unsolvable = true;
    return;
  }
  int Ki;
  if( case_num == 0 ){
    Ki = 0;
  } else if( case_num == 1 ){
    Ki = ( exponent / valuation ) + 1;
    Ki += log( prime_curr , Ki );
  } else {
    Ki = exponent + 1;
    Ki += log( prime_curr , Ki );
    if( Ki >= prime_curr ){
      Ki = prime_curr_minus;
    }
  }
  local_answer.push_back( 0 );
  ll& local_answer_curr = local_answer.back();
  ll M_shifted_power = 1;
  vector<ll> inverse{};
  inverse.push_back( 0 );
  inverse.push_back( 1 );
  int inverse_size = 2;
  FOREQ( k , 1 , Ki ){
    int k_coprime = k;
    int valuation_temp = valuation * k;
    while( k_coprime % prime_curr == 0 ){
      k_coprime /= prime_curr;
      valuation_temp--;
    }
    k_coprime %= factor_power_curr;
    while( k_coprime >= inverse_size ){
      if( inverse_size % prime_curr == 0 ){
	inverse.push_back( 0 );
      } else {
	// gcd( prime_curr , inverse_size ) == 1かつ
	// inverse_size != 1より
	// inverse_size_sub != 0である。
	const ll inverse_size_sub = factor_power_curr % inverse_size;
	const ll& inverse_size_sub_inv = inverse[inverse_size_sub];
	if( inverse_size_sub_inv == 0 ){
	  // inverse_size_sub != 0より
	  // inverse_size_sub_minus < inverse_sizeである。
	  const ll inverse_size_sub_minus = inverse_size - inverse_size_sub;
	  // gcd( prime_curr , inverse_size ) == 1かつ
	  // inverse_size == inverse_size_sub + inverse_size_sub_minusかつ
	  // gcd( prime_curr , inverse_size_sub ) != 1より
	  // gcd( prime_curr , inverse_size_sub_minus ) == 1である。
	  const ll& inverse_size_sub_minus_inv = inverse[ inverse_size_sub_minus ];
	  inverse.push_back( ( inverse_size_sub_minus_inv * ( ( factor_power_curr / inverse_size ) + 1 ) ) % factor_power_curr );
	} else {
	  inverse.push_back( factor_power_curr - ( ( inverse_size_sub_inv * ( factor_power_curr / inverse_size ) ) % factor_power_curr ) );
	}
      }
      inverse_size++;
    }
    ll& k_coprime_inverse = inverse[k_coprime];
    POWER_MOD( prime_curr_power , prime_curr , valuation_temp , factor_power_curr );
    local_answer_curr += ( ( ( M_shifted_power *= M_shifted ) %= factor_power_curr ) * ( ( k_coprime_inverse * prime_curr_power ) % factor_power_curr ) ) % factor_power_curr;
  }
  if( order < inverse_size ){
    ( local_answer_curr *= inverse[order] ) %= factor_power_curr;
  } else {
    POWER_MOD( order_inverse , order , factor_power_curr - factor_power_curr / prime_curr - 1 , factor_power_curr );
    ( local_answer_curr *= order_inverse ) %= factor_power_curr;
  }
  ( local_answer_curr *= -1 ) < 0 ? local_answer_curr += factor_power_curr : local_answer_curr;
  return;
}

int MAIN()
{
  UNTIE;
  // 1000以下の素数の最大値+1の前計算値
  CEXPR( int , prime_lim , 998 );
  // 1000以下の素数の総数の前計算値
  CEXPR( int , length_max , 168 );
  constexpr PrimeEnumeration<int,prime_lim,length_max> prime{};
  CEXPR( int , bound_T , 100000 );
  CIN_ASSERT( T , 1 , bound_T );
  CEXPR( int , bound_B , 1000000 );
  CEXPR( ll , bound_M , 1000000000000000000 );
  REPEAT( T ){
    CIN_ASSERT( B , 2 , bound_B );
    CIN_ASSERT( N , 0 , log( B , bound_B ) );
    CIN_ASSERT( M , 0 , bound_M );
    int factor_num = 0;
    vector<ll> factor{};
    vector<ll> factor_power{};
    vector<ll> local_answer{};
    bool unsolvable = false;
    FOR( i , 0 , length_max ){
      const int& prime_i = prime.m_val[i];
      if( B % prime_i == 0 ){
	ComputeLocalAnswer( factor_num , factor , factor_power , local_answer , unsolvable , B , N , M , prime_i );
      } else if( B / prime_i < prime_i ){
	break;
      }
      if( unsolvable ){
	break;
      }
    }
    if( B != 1 && ! unsolvable ){
      int prime_last = B;
      ComputeLocalAnswer( factor_num , factor , factor_power , local_answer , unsolvable , B , N , M , prime_last );
    }
    if( unsolvable ){
      COUT( -1 );
    } else {
      ll base = 1;
      ll answer = 0;
      FOR( i , 0 , factor_num ){
	ll lcm;
	answer = ChineseRemainderTheorem( base , answer , factor_power[i] , local_answer[i] , lcm );
	base = lcm;
      }
      COUT( answer );
    }
  }
  QUIT;
}

yukicoder

結果

テストケース

ソースコード