結果

問題 No.950 行列累乗
コンテスト
ユーザー chocorusk
提出日時 2019-12-10 07:59:38
言語 C++14
(gcc 13.3.0 + boost 1.87.0)
結果
WA  
実行時間 -
コード長 5,920 bytes
コンパイル時間 2,270 ms
コンパイル使用メモリ 143,780 KB
実行使用メモリ 12,544 KB
最終ジャッジ日時 2024-06-23 20:03:51
合計ジャッジ時間 4,979 ms
ジャッジサーバーID
（参考情報） 		judge3 / judge2
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
sample AC * 4
other AC * 51 WA * 6
権限があれば一括ダウンロードができます

ソースコード

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];
	}
	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){
	ll a= A[0][0], b=A[0][1], c=A[1][0], d = A[1][1];
	mat P = {{b, b}, {p-a,d}};
	ll inv = powmod(b*(a+d)%p, p-2, p);
	mat Pi = {{d*inv%mod,(p-b)*inv%mod},{a*inv%mod, b*inv%mod}};
	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])-4*(A[0][0]*A[1][1]%p-A[0][1]*A[1][0]%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