很容易看出来这道题是求模n意义下fib数列的最小循环节
对于fib数列的最小循环节的求法,我们可以这样:
1、令n=p1^m1 * p2^m2 * p3^m3……
2、分别计算fib数列在模p1^m1,p2^m2……意义下的最小循环节
3、模n意义下的最小循环节为2步骤各循环节的LCM
首先步骤三是非常容易证明的,CRT直接可以看出来
我们考虑如何计算p1^m1的循环节
这里有一个定理:G(p1^m1)=G(p1)*p1^(m-1)
我们只需要计算G(p)就可以了
我们知道对于Fib数列,满足f(i)=1.0/sqrt(5) *( (1+sqrt(5))/2)^i - ((1-sqrt(5))/2)^i)
那么在模p意义下存在以下两种情况:
1、5是模p意义下的二次剩余,即满足欧拉判别 5^((p-1)/2)=1
那么原式中的sqrt(5)在模意义下可以用整数x替代
由费马小定理很容易知道F(p-1)=0
又知道
((1+x)/2)^p=(1+x)/2
((1-x)/2)^p=(1-x)/2
两式做差等于x
所以F(p)=1
因此我们知道p-1是Fib数列的一个循环节,那么其最小循环节一定是p-1的约数
枚举+矩阵乘法判断即可
2、5不是模p意义下的二次剩余,则存在5^((p-1)/2)=-1
那么由二项式定理得
(1+sqrt(5))^p=1+(sqrt(5))^p
(1-sqrt(5))^p=1-(sqrt(5))^p
(因为其系数的组合数里含有p这个因子的时候值为0)
两式做差得2*sqrt(5)^p
进而得F(p)=(sqrt(5)/2)^(p-1)
我们知道在模p意义下2^(p-1)的逆元为1
所以F(p)=(sqrt(5))^(p-1)=5^((p-1)/2)=-1
同理,我们可以得到F(p+1)=0
设Fib数列在模意义下第一个F(i)=0且i不等于0的位置为i
设F(i-1)=a 易证(a,p)=1
则一定有i是p+1的约数
且F(p)=a^j*F(p-i*j)
显然F(p-i*j)=a,所以F(p)=a^(j+1)
进而知道F(2*p+1)=a^(2*j+2)=(-1)^2=1
F(2*p+2)=0
所以2*(p+1)是Fib数列的一个循环节,那么其最小循环节一定是2*(p+1)的约数
枚举+矩阵乘法判断即可
#include<cstdio> #include<cstring> #include<cstdlib> #include<iostream> #include<algorithm> using namespace std; const int maxn=20010; typedef long long LL; LL n,mod,now; LL a,b,ans; LL st[maxn],top=0; struct Matrix{ LL a[2][2]; void clear(){memset(a,0,sizeof(a));} void print(){ for(int i=0;i<2;++i){ for(int j=0;j<2;++j){ printf("%lld ",a[i][j]); }printf(" "); }printf(" "); } }A; LL mul(LL a,LL b){ LL s=0; while(b){ if(b&1)s=(s+a)%mod; a=(a<<1)%mod;b>>=1; }return s; } LL GCD(LL a,LL b){return b==0?a:GCD(b,a%b);} LL LCM(LL a,LL b){return a*b/GCD(a,b);} LL Get_pow(LL v,LL p){ LL tmp=1; while(p){ if(p&1)tmp=mul(tmp,v); v=mul(v,v);p>>=1; }return tmp; } Matrix operator *(const Matrix &A,const Matrix &B){ Matrix C;C.clear(); for(int i=0;i<2;++i){ for(int j=0;j<2;++j){ for(int k=0;k<2;++k){ C.a[i][j]=C.a[i][j]+mul(A.a[i][k],B.a[k][j]); if(C.a[i][j]>=mod)C.a[i][j]-=mod; } } }return C; } Matrix pow_mod(Matrix v,LL p){ Matrix tmp;tmp.clear(); tmp.a[0][0]=tmp.a[1][1]=1; while(p){ if(p&1)tmp=tmp*v; v=v*v;p>>=1; }return tmp; } void F(LL p){ A.a[0][0]=1;A.a[0][1]=1; A.a[1][0]=1;A.a[1][1]=0; A=pow_mod(A,p-1); b=(A.a[0][0]+A.a[0][1])%mod; a=(A.a[1][0]+A.a[1][1])%mod; } LL Get_loop(LL k){ if(k==2)return 3; else if(k==3)return 8; else if(k==5)return 20; LL p;mod=k;top=0; if(Get_pow(5LL,(mod-1)>>1)==1)p=mod-1; else p=2*(mod+1); for(LL i=1;i*i<=p;++i){ if(p%i==0){ LL x=i,y=p/i; F(x); if(a==0&&b==1)return x; if(y!=x)st[++top]=y; } } while(top){ F(st[top]); if(a==0&&b==1)return st[top]; top--; }return 0; } int main(){ while(scanf("%lld",&n)==1){ now=n;ans=1; for(LL i=2;i*i<=now;++i){ if(now%i==0){ LL S=1; while(now%i==0)now/=i,S*=i; S/=i;S=S*Get_loop(i); ans=LCM(ans,S); } } if(now>1){ LL S=Get_loop(now); ans=LCM(ans,S); } if(ans%2==0)ans/=2; printf("%lld ",ans); }return 0; }