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

提出ソース

結果

問題	No.2993 冪乗乗 mod 冪乗
ユーザー	👑 p-adic
提出日時	2023-02-02 10:53:51
言語	C++17 (gcc 12.3.0 + boost 1.83.0)
結果	AC
実行時間	409 ms / 10,000 ms
コード長	12,395 bytes
コンパイル時間	3,886 ms
コンパイル使用メモリ	227,240 KB
実行使用メモリ	6,820 KB
最終ジャッジ日時	2024-12-19 16:19:29
合計ジャッジ時間	12,395 ms
ジャッジサーバーID （参考情報）	judge3 / judge4

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

テストケース

テストケース表示

入力	結果	実行時間実行使用メモリ
testcase_00	AC	242 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,820 KB
testcase_07	AC	6 ms 6,816 KB
testcase_08	AC	3 ms 6,816 KB
testcase_09	AC	2 ms 6,820 KB
testcase_10	AC	2 ms 6,820 KB
testcase_11	AC	397 ms 6,816 KB
testcase_12	AC	276 ms 6,816 KB
testcase_13	AC	281 ms 6,816 KB
testcase_14	AC	194 ms 6,816 KB
testcase_15	AC	163 ms 6,816 KB
testcase_16	AC	359 ms 6,816 KB
testcase_17	AC	343 ms 6,820 KB
testcase_18	AC	391 ms 6,816 KB
testcase_19	AC	283 ms 6,820 KB
testcase_20	AC	349 ms 6,816 KB
testcase_21	AC	358 ms 6,820 KB
testcase_22	AC	270 ms 6,816 KB
testcase_23	AC	409 ms 6,820 KB
testcase_24	AC	284 ms 6,816 KB
testcase_25	AC	283 ms 6,816 KB
testcase_26	AC	373 ms 6,820 KB
testcase_27	AC	221 ms 6,816 KB
testcase_28	AC	158 ms 6,820 KB
testcase_29	AC	92 ms 6,820 KB
testcase_30	AC	92 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( int , N );
    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

結果

テストケース

ソースコード