很容易看出来这道题是求模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;
}