[斐波那契數列 求解二次剩餘&二次模方程 BSGS] COGS 2114 [CodeChef FN]斐波那契數
阿新 • • 發佈:2019-02-17
我們可以考慮斐波那契的通項公式 換元可以轉化成一個二次方程 根據求根公式求出後 再用BSGS求出指數
注意分奇偶考慮
#include<cstdio> #include<cstdlib> #include<algorithm> #include<cstring> #include<cmath> #define cl(x) memset(x,0,sizeof(x)) using namespace std; typedef long long ll; inline char nc(){ static char buf[100000],*p1=buf,*p2=buf; if (p1==p2) { p2=(p1=buf)+fread(buf,1,100000,stdin); if (p1==p2) return EOF; } return *p1++; } inline void read(ll &x){ char c=nc(),b=1; for (;!(c>='0' && c<='9');c=nc()) if (c=='-') b=-1; for (x=0;c>='0' && c<='9';x=x*10+c-'0',c=nc()); x*=b; } namespace BSGS{ inline ll Pow(ll a,ll b,ll p){ ll ret=1; for (;b;(a*=a)%=p,b>>=1) if (b&1) (ret*=a)%=p; return ret; } struct hashmap{ #define MOD 1000007 ll num; ll key[1000005],value[1000005],next[1000005]; ll head[MOD]; inline void insert(ll k,ll v){ ll ad=k%MOD; key[++num]=k; value[num]=v; next[num]=head[ad]; head[ad]=num; } inline ll operator[](ll k){ ll ad=k%MOD; for (int p=head[ad];p;p=next[p]) if (key[p]==k) return value[p]; return -1; } inline void clear(){ cl(head); num=0; } }M; ll lastA=-1,lastP=-1; inline ll Solve(ll A,ll B,ll P,ll phi){ A%=P; if (!A && !B) return 1; if (!A) return 1LL<<60; ll m=sqrt((double)P)+1; if (lastA!=A || lastP!=P){ M.clear(); ll tmp=1; M.insert(1,0); for (int i=1;i<m;i++){ (tmp*=A)%=P; if (M[tmp]==-1) M.insert(tmp,i); } lastA=A; lastP=P; } ll inv=1,T=Pow(A,phi-m,P); for (int k=0;k<m;k++){ ll i=M[B*inv%P]; if (i!=-1) return k*m+i; (inv*=T)%=P; } return 1LL<<60; } } namespace ROOT{ ll P; inline ll Pow(ll a,ll b){ ll ret=1; for (;b;b>>=1,a=a*a%P) if (b&1) ret=ret*a%P; return ret; } inline ll legendre(ll a){ return Pow(a,(P-1)>>1); } struct abcd{ ll a,b,w; //a+b*sqrt(w) abcd(ll a=0,ll b=0,ll w=0):a(a),b(b),w(w) { } friend abcd operator *(abcd A,abcd B){ return abcd((A.a*B.a%P+A.b*B.b%P*A.w%P)%P,(A.a*B.b%P+A.b*B.a%P)%P,A.w); } }; inline abcd Pow(abcd a,int b){ abcd ret=abcd(1,0,a.w); for (;b;b>>=1,a=a*a) if (b&1) ret=ret*a; return ret; } inline ll calc(ll n){ //if (P==2) return 1; n%=P; if (!n) return 0; if (legendre(n)==P-1) return -1; ll a,w; while (1){ a=rand()%P; w=((a*a-n)%P+P)%P; if (legendre(w)==P-1) break; } return Pow(abcd(a,1,w),(P+1)>>1).a; } inline pair<ll,ll> Root(ll n,ll p){ P=p; ll ans=calc(n); if (ans==-1) return make_pair(-1LL,-1LL); return make_pair(ans,P-ans); } } ll P,A,SQRT5,INV2; inline ll Pre(ll p){ P=p; INV2=(P+1)>>1; SQRT5=ROOT::Root(5,p).first; A=(1+SQRT5)*INV2%P; } inline ll Pow(ll a,ll b){ ll ret=1; for (;b;b>>=1,a=a*a%P) if (b&1) ret=ret*a%P; return ret; } inline ll Inv(ll a){ return Pow(a,P-2); } inline ll SolveEven(ll F){ ll C=F*SQRT5%P; pair<ll,ll> tem=ROOT::Root(C*C+4,P); if (tem.first==-1) return 1LL<<60; ll x1=(C+tem.first)*INV2%P,x2=(C+tem.second)*INV2%P; ll ret=1LL<<60; tem=ROOT::Root(x1,P); if (tem.first!=-1) ret=min(ret,BSGS::Solve(A,tem.first,P,P-1)<<1); if (tem.second!=tem.first) ret=min(ret,BSGS::Solve(A,tem.second,P,P-1)<<1); if (x1!=x2){ tem=ROOT::Root(x2,P); if (tem.first!=-1) ret=min(ret,BSGS::Solve(A,tem.first,P,P-1)<<1); if (tem.second!=tem.first) ret=min(ret,BSGS::Solve(A,tem.second,P,P-1)<<1); } return ret; } inline ll SolveOdd(ll F){ ll C=F*SQRT5%P; pair<ll,ll> tem=ROOT::Root(C*C+P-4,P); if (tem.first==-1) return 1LL<<60; ll x1=(C+tem.first)*INV2%P,x2=(C+tem.second)*INV2%P; ll ret=1LL<<60; tem=ROOT::Root(x1*A,P); if (tem.first!=-1) ret=min(ret,(BSGS::Solve(A,tem.first,P,P-1)<<1)-1); if (tem.second!=tem.first) ret=min(ret,(BSGS::Solve(A,tem.second,P,P-1)<<1)-1); if (x1!=x2){ tem=ROOT::Root(x2*A,P); if (tem.first!=-1) ret=min(ret,(BSGS::Solve(A,tem.first,P,P-1)<<1)-1); if (tem.second!=tem.first) ret=min(ret,(BSGS::Solve(A,tem.second,P,P-1)<<1)-1); } return ret; } int main(){ ll T,f,p,tem,Ans; freopen("fn.in","r",stdin); freopen("fn.out","w",stdout); read(T); while (T--){ read(f); read(p); Pre(p); Ans=1LL<<60; tem=SolveEven(f); Ans=min(Ans,tem); tem=SolveOdd(f); Ans=min(Ans,tem); if (Ans==1LL<<60) printf("-1\n"); else printf("%lld\n",Ans); } return 0; }