結果

問題 No.950 行列累乗
ユーザー tempura_pptempura_pp
提出日時 2019-12-10 08:13:04
言語 C++14
(gcc 12.3.0 + boost 1.83.0)
結果
WA  
(最新)
AC  
(最初)
実行時間 -
コード長 6,064 bytes
コンパイル時間 2,263 ms
コンパイル使用メモリ 144,520 KB
実行使用メモリ 12,580 KB
最終ジャッジ日時 2023-09-09 10:42:07
合計ジャッジ時間 5,570 ms
ジャッジサーバーID
（参考情報） 		judge13 / judge15
このコードへのチャレンジ
（要ログイン）

テストケース

テストケース表示
入力 結果 実行時間
実行使用メモリ
testcase_00 WA -
testcase_01 AC 2 ms
4,376 KB
testcase_02 AC 1 ms
4,376 KB
testcase_03 AC 2 ms
4,376 KB
testcase_04 AC 2 ms
4,380 KB
testcase_05 AC 3 ms
4,376 KB
testcase_06 AC 2 ms
4,380 KB
testcase_07 AC 2 ms
4,376 KB
testcase_08 AC 2 ms
4,376 KB
testcase_09 AC 1 ms
4,380 KB
testcase_10 AC 1 ms
4,376 KB
testcase_11 AC 2 ms
4,376 KB
testcase_12 AC 2 ms
4,376 KB
testcase_13 AC 2 ms
4,380 KB
testcase_14 AC 2 ms
4,376 KB
testcase_15 AC 2 ms
4,380 KB
testcase_16 AC 2 ms
4,380 KB
testcase_17 AC 1 ms
4,376 KB
testcase_18 AC 1 ms
4,380 KB
testcase_19 AC 2 ms
4,376 KB
testcase_20 AC 2 ms
4,376 KB
testcase_21 AC 3 ms
4,380 KB
testcase_22 AC 5 ms
4,376 KB
testcase_23 AC 4 ms
4,376 KB
testcase_24 AC 3 ms
4,376 KB
testcase_25 AC 5 ms
4,376 KB
testcase_26 AC 4 ms
4,380 KB
testcase_27 AC 2 ms
4,380 KB
testcase_28 AC 7 ms
4,380 KB
testcase_29 AC 12 ms
4,380 KB
testcase_30 AC 3 ms
4,376 KB
testcase_31 AC 77 ms
9,480 KB
testcase_32 AC 129 ms
12,304 KB
testcase_33 AC 98 ms
10,264 KB
testcase_34 AC 146 ms
12,196 KB
testcase_35 AC 92 ms
10,340 KB
testcase_36 AC 2 ms
4,376 KB
testcase_37 AC 2 ms
4,376 KB
testcase_38 AC 3 ms
4,384 KB
testcase_39 AC 2 ms
4,376 KB
testcase_40 AC 4 ms
4,376 KB
testcase_41 AC 4 ms
4,380 KB
testcase_42 AC 3 ms
4,376 KB
testcase_43 AC 4 ms
4,380 KB
testcase_44 AC 4 ms
4,376 KB
testcase_45 AC 3 ms
4,376 KB
testcase_46 AC 93 ms
9,888 KB
testcase_47 AC 3 ms
4,376 KB
testcase_48 AC 38 ms
6,284 KB
testcase_49 AC 61 ms
6,208 KB
testcase_50 AC 148 ms
12,580 KB
testcase_51 AC 232 ms
12,472 KB
testcase_52 AC 8 ms
4,376 KB
testcase_53 AC 4 ms
4,376 KB
testcase_54 AC 1 ms
4,376 KB
testcase_55 AC 2 ms
4,376 KB
testcase_56 AC 3 ms
4,376 KB
testcase_57 AC 1 ms
4,376 KB
testcase_58 AC 2 ms
4,380 KB
testcase_59 AC 2 ms
4,376 KB
testcase_60 AC 3 ms
4,376 KB
権限があれば一括ダウンロードができます

ソースコード

diff # 

#include<iostream>
#include<string>
#include<algorithm>
#include<vector>
#include<iomanip>
#include<math.h>
#include<complex>
#include<queue>
#include<deque>
#include<stack>
#include<map>
#include<set>
#include<bitset>
#include<functional>
#include<assert.h>
#include<numeric>
using namespace std;
#define REP(i,m,n) for(int i=(int)(m) ; i < (int) (n) ; ++i )
#define rep(i,n) REP(i,0,n)
using ll = long long;
const int inf=1e9+7;
const ll longinf=1LL<<60 ;
const ll mod=1e9+7 ;


ll powmod(ll n,ll k, ll p){
	ll ret=1;
	while(k){
		if(k&1)ret=ret*n%p;
		n=n*n%p;
		k>>=1;
	}
	return ret;
}

using mat = vector<vector<ll>>;

mat mul(mat A, mat B, ll p){
	mat C(2,vector<ll>(2));
	rep(i,2)rep(j,2)rep(k,2){
		C[i][j]+=A[i][k]*B[k][j]%p;
	}
	rep(i,2)rep(j,2)C[i][j]%=p;
	return C;
}

mat powmat(mat A, ll n, ll p){
	mat ret(2,vector<ll>(2));
	ret[0][0]=ret[1][1]=1;
	while(n){
		if(n&1)ret=mul(ret,A,p);
		A = mul(A,A,p);
		n>>=1;
	}
	return ret;
}

ll calc(mat A, mat B, ll T, ll p){
	ll sq = 1;
	while(sq*sq<T)++sq;
	map<mat,int> mp;
	mat C = powmat(A,sq-1,p);
	mat Inv = powmat(A,T-1,p);
	for(int i=sq-1;i>=0;i--){
		mp[C]=i;
		C = mul(C,Inv,p);
	}
	ll ret = longinf;
	C = powmat(A,T,p);
	Inv = powmat(A,T-sq,p);
	rep(i,sq){
		if(T-sq*i<0)break;
		mat g = mul(B,C,p);
		if(mp.count(g))ret= min(ret,mp[g]+sq*i);
		C = mul(C,Inv,p);
	}
	if(ret>T)return -1;
	return ret;
}


ll calc2(mat A, mat B, ll T, ll p){
	if(A[0][1]==0){
		swap(A[0][1],A[1][0]);
		swap(B[1][0],B[0][1]);
	}
	ll a= A[0][0], b=A[0][1], c=A[1][0], d = A[1][1];
	mat P = {{b, b}, {-a,d}};
	ll inv = powmod(b*(a+d)%p, p-2, p);
	mat Pi = {{d*inv%p,-b*inv%p},{a*inv%p, b*inv%p}};
	if(b==0){
		P = {{1,0},{0,1}};
		Pi = {{1,0},{0,1}};
	}
	auto D = mul(P,Pi,p);
	B = mul(Pi,mul(B,P,p),p);
	if(B[0][0]||B[0][1]||B[1][0])return -1;
	ll sq = 1;
	while(sq*sq<T)++sq;
	map<ll,int> mp;
	ll x=(a+d)%p;
	for(int i=sq-1;i>=0;i--){
		mp[powmod(x,i,p)]=i;
	}
	ll ret = longinf;
	rep(i,sq){
		if(T-sq*i<0)break;
		ll g = B[1][1]*powmod(x,T-sq*i,p)%p;
		if(mp.count(g))ret= min(ret,mp[g]+sq*i);
	}
	if(ret>T)return -1;
	return ret;
}

map<int,int> fact(ll p, bool f){
	ll tmp=p-1;
	map<int,int>ret;
	for(ll i=2;i*i<=p-1;++i){
		if(tmp%i==0){
			while(tmp%i==0){
				ret[i]++;
				tmp/=i;
			}
		}
	}
	if(tmp>1)ret[tmp]++;
	tmp = p+f;
	for(ll i=2;i*i<=p+1;++i){
		if(tmp%i==0){
			while(tmp%i==0){
				ret[i]++;
				tmp/=i;
			}
		}
	}
	if(tmp>1)ret[tmp]++;
	return ret;
}


long long extGcd(long long a, long long b, long long &p, long long &q) {  
    if (b == 0) { p = 1; q = 0; return a; }  
    long long d = extGcd(b, a%b, q, p);  
    q -= a/b * p;  
    return d;  
}

// 中国剰余定理
// リターン値を (r, m) とすると解は x ≡ r (mod. m)
// 解なしの場合は (0, -1) をリターン
pair<long long, long long> ChineseRem(long long b1, long long m1, long long b2, long long m2) {
  long long p, q;
  long long d = extGcd(m1, m2, p, q); // p is inv of m1/d (mod. m2/d)
  if ((b2 - b1) % d != 0) return make_pair(0, -1);
  long long m = m1 * (m2/d); // lcm of (m1, m2)
  long long tmp = (b2 - b1) / d * p % (m2/d);
  long long r = (b1 + m1 * tmp)% m;
  if(r<0)r+=m;
  return make_pair(r, m);
}

int main(){
	ll p;
	cin>>p;
	mat A(2,vector<ll>(2)), B(2,vector<ll>(2));
	rep(i,2)rep(j,2)cin>>A[i][j];
	rep(i,2)rep(j,2)cin>>B[i][j];
	if(A==B){
		cout<<1<<endl;return 0;
	}
	ll D = (A[0][0]+A[1][1])%p*(A[0][0]+A[1][1])%p-4*(A[0][0]*A[1][1]%p-A[0][1]*A[1][0]%p)%p;
	D%=p;
	if(D<0)D+=p;
	if(D==0){
		if((A[0][0]+A[1][1])%p==0){
			if(A==B){
				cout<<1<<endl;return 0;
			}
			else if(mul(A,A,p)==B){
				cout<<2<<endl;return 0;
			}
			else {
				cout<<-1<<endl;return 0;
			}
		}
		else {
			auto fac = fact(p,false);
			mat E = powmat(A,0,p);
			ll T = p*(p-1);
			for(auto& e:fac){
				while(e.second && powmat(A,T/e.first,p)==E){
					--e.second;
					T/=e.first;
				}
			}
			vector<ll> v, r;
			for(auto e:fac){
				if(e.second==0)continue;
				ll q = 1;
				rep(i,e.second)q*=e.first;
				ll res = calc(powmat(A,T/q,p),powmat(B,T/q,p),q,p);
				if(res==-1){
					cout<<-1<<endl;return 0;
				}
				v.push_back(q);r.push_back(res);
			}
			pair<ll,ll> ans = {0,1};
			rep(i,v.size()){
				ans = ChineseRem(ans.first,ans.second, r[i],v[i]);
				if(ans.second==-1){
					cout<<-1<<endl;return 0;
				}
			}
			if(ans.first==0)ans.first=ans.second;
			if(powmat(A,ans.first,p)!=B){
				cout<<-1<<endl;return 0;
			}
			cout<<ans.first<<endl;return 0;
		}
	}
	else {
		if(A[0][0]*A[1][1]%p==A[1][0]*A[0][1]%p){
			
			auto fac = fact(p,true);
			ll T = (p+1)*(p-1);
			ll d = (A[0][0]+A[1][1])%p;
			for(auto& e:fac){
				while(e.second && powmod(d,T/e.first,p)==1){
					--e.second;
					T/=e.first;
				}
			}
			vector<ll> v, r;
			for(auto e:fac){
				if(e.second==0)continue;
				ll q = 1;
				rep(i,e.second)q*=e.first;
				ll res = calc2(powmat(A,T/q,p),powmat(B,T/q,p),q,p);
				if(res==-1){
					cout<<-1<<endl;return 0;
				}
				v.push_back(q);r.push_back(res);
			}
			pair<ll,ll> ans = {0,1};
			rep(i,v.size()){
				ans = ChineseRem(ans.first,ans.second, r[i],v[i]);
				if(ans.second==-1){
					cout<<-1<<endl;return 0;
				}
			}
			if(ans.first==0)ans.first=ans.second;
			if(powmat(A,ans.first,p)!=B){
				cout<<-1<<endl;return 0;
			}
			cout<<ans.first<<endl;return 0;
		}
		else {
			auto fac = fact(p,true);
			mat E = powmat(A,0,p);
			ll T = (p+1)*(p-1);
			for(auto& e:fac){
				while(e.second && powmat(A,T/e.first,p)==E){
					--e.second;
					T/=e.first;
				}
			}
			vector<ll> v, r;
			for(auto e:fac){
				if(e.second==0)continue;
				ll q = 1;
				rep(i,e.second)q*=e.first;
				ll res = calc(powmat(A,T/q,p),powmat(B,T/q,p),q,p);
				if(res==-1){
					cout<<-1<<endl;return 0;
				}
				v.push_back(q);r.push_back(res);
			}
			pair<ll,ll> ans = {0,1};
			rep(i,v.size()){
				ans = ChineseRem(ans.first,ans.second, r[i],v[i]);
				if(ans.second==-1){
					cout<<-1<<endl;return 0;
				}
			}
			if(ans.first==0)ans.first=ans.second;
			if(powmat(A,ans.first,p)!=B){
				cout<<-1<<endl;return 0;
			}
			cout<<ans.first<<endl;return 0;
		}
	}
	return 0;
}
0