求解斐波那契数列模$p$意义下最短循环节

如题，毕克老师给我们出的noip(NOIplus)模拟赛的(Day1T1)

首先我们知道斐波那契数列的特征根

[phi_1=frac{1+sqrt{5}}{2} ]

[phi_2=frac{1-sqrt{5}}{2} ]

于是

[F_n=frac{phi_1^n-phi_2^n}{sqrt{5}} ]

对于不为(5)的质数(p)

若(5)是模(p)意义下的二次剩余

设最小循环节长度为(m)
那么(phi_1,phi_2)是模(p)的完全剩余系中的元素
根据费马小定理

[phi_1^{p-1} equiv 1 mod p​ ]

[phi_2^{p-1} equiv 1 mod p​ ]

[F_{p-1} equiv frac{phi_1^{p-1}-phi_2^{p-1}}{sqrt{5}} equiv 0 mod p ]

[F_{p} equiv frac{phi_1^{p}-phi_2^{p}}{sqrt{5}} equiv frac{phi_1-phi_2}{sqrt{5}} equiv 1 mod p ]

[m | p-1 ]

若(5)是模(p)意义下的二次非剩余

设最小循环节长度为(m)
根据欧拉判别准则

[5^{frac{p-1}{2}} equiv -1 mod p ]

[phi_1^p equiv (frac{1+sqrt{5}}{2})^p equiv (frac{1}{2})^p (1+sqrt{5})^p equiv (frac{1}{2})^p (1+sqrt{5}^p) equiv frac{1}{2} (1+5^{frac{p-1}{2}} sqrt{5}) equiv phi_2 mod p ]

同理可得

[phi_2^p equiv phi_1 mod p ]

[F_{2p+1} equiv frac{phi_1^{2p+1}-phi_2^{2p+1}}{sqrt{5}} equiv frac{phi_1^{2p}phi_1-phi_2^{2p}phi_2}{sqrt{5}} equiv frac{phi_2^{2}phi_1-phi_1^{2}phi_2}{sqrt{5}} equiv phi_1phi_2frac{phi_2-phi_1}{sqrt{5}} equiv 1 mod p ]

[F_{2p+2} equiv frac{phi_1^2phi_2^2-phi_2^2phi_1^2}{sqrt{5}} equiv 0 mod p ]

[F_{2p+3} equiv F_{2p+1}+ F_{2p+2}equiv 1 mod p ]

[m|2p+2 ]

对于质数的幂(p^k)

设模(p)意义下的最小循环节长度为(m)，模(p^k)意义下的最小循环节长度为(m')

[F_m equiv frac{phi_1^m-phi_2^m}{sqrt{5}} equiv 0 mod p ]

[phi_1^m equiv phi_2^m mod p ]

[F_{m+1} equiv F_1 mod p ]

[phi_1^{m+1}-phi_2^{m+1}-phi_1+phi_2equiv phi_1(phi_1^m-1)-phi_2(phi_2^m - 1) equiv ( phi_1-phi_2)(phi_1^m-1) equiv 0 mod p ]

[phi_1^m equiv phi_2^m equiv 1 mod p ]

[{(phi_1^m)}^{p^{k-1}} equiv {(phi_2^m)}^{p^{k-1}} equiv 1 mod p^k ]

注:这里引用了一个定理

若(a equiv 1 mod p),则(a^{p^k} equiv 1 mod p^{k+1})

可以用数学归纳法证明
(k=1)时，令(a=np+1)有

[a^p equiv (np+1)^p equiv sum_{i=0}^{p}{C_p^i* (np)^p} equiv 1+p*(np)+sum_{i=2}^{p}{C_p^i* (np)^p} equiv 1 mod p^2 ]

若当(k=m)时成立，令(a^{p^m}=np^{m+1}+1)有

[a^{p^{m+1}} equiv (a^{p^m})^p equiv (np^{m+1}+1)^p equiv 1+p*np^{m+1}+sum_{i=2}^{p}{C_p^i *(np^{m+1})^p} equiv 1 mod p^{m+2} ]

得证

所以

[F_{mp^{k-1}} equiv 0 mod p^k ]

[F_{mp^{k-1}+1} equiv 1 mod p^k ]

[m'|mp^{k-1} ]

注：
就我所知对于目前已知的所有情况,都有$ m'=mp^{k-1} $
然而我所能查到的论文都说数学上还没有证明
如果哪位大佬知道请赐教……

对于合数(p_1^{k_1}*p_2^{k_2}*...*p_c^{k_c})

设模(p_1^{k_1}*p_2^{k_2}*...*p_c^{k_c})意义下的最小循环节长度为(m)

[F_m equiv 0 mod p_1^{k_1}*p_2^{k_2}*...*p_c^{k_c} ]

[F_{m+1} equiv 1 mod p_1^{k_1}*p_2^{k_2}*...*p_c^{k_c} ]

由中国剩余定理

[F_m equiv 0 mod p_1^{k_1} ]

[F_{m+1} equiv 1 mod p_1^{k_1} ]

[F_m equiv 0 mod p_2^{k_1} ]

[F_{m+1} equiv 1 mod p_2^{k_1} ]

[... ]

[F_m equiv 0 mod p_c^{k_c} ]

[F_{m+1} equiv 1 mod p_c^{k_c} ]

所以在模(p_i^{k_i})意义下分别求解
答案取(lcm)即可

优化&总结&代码

关于判断(5)是否是模(p)((p)为奇素数且(p eq 5))意义下的二次剩余
可以直接用欧拉判别准则
但是需要快速幂，带一个(log)
我们还可以考虑优化这个算法
由二次互反律
可以得到

[(frac{5}{p})=(frac{p}{5})(-1)^{frac{(5-1)(p-1)}{4}}=(frac{p}{5}) ]

因为

[(frac{1}{5}) = (frac{4}{5}) =1 ]

[(frac{2}{5}) = (frac{3}{5}) =-1 ]

所以
当且仅当(p equiv 1 mod 5)或(p equiv 4 mod 5)时(5)是模(p)意义下的二次剩余
当且仅当(p equiv 2 mod 5)或(p equiv 3 mod 5)时(5)是模(p)意义下的二次非剩余

(2)和(5)比较特殊
同样由二次互反律可以注意到(2)的条件是相反的
而(5)则不在以上所有讨论范围内，没有理论上的证明
然而不是有样例吗……
所以小于等于(5)的我都直接返回答案了(实际上(1e6)内都可以暴力跑出来)
对于其它质数求解时直接枚举(p-1)或(2p+2)的约数判断即可

然后这道题就完美地做完了

#include<bits/stdc++.h>

using namespace std;

#define gc c=getchar()
#define r(x) read(x)
#define ll long long

template<typename T>
inline void read(T&x){
    x=0;T k=1;char gc;
    while(!isdigit(c)){if(c=='-')k=-1;gc;}
    while(isdigit(c)){x=x*10+c-'0';gc;}x*=k;
}

const int N=5e4+7;

int tot;
int pri[N];
bool mark[N];

inline void init(){
	for(int i=2;i<N;++i){
		if(!mark[i])pri[++tot]=i;
		for(int j=1,tmp;j<=tot&&(tmp=i*pri[j])<N;++j){
			mark[tmp]=1;
			if(i%pri[j]==0)break;
		}
	}
}

int p;

inline int add(int a,int b){
	a+=b;
	if(a>=p)a-=p;
	return a;
}

inline int mul(int a,int b){
	return (ll)a*b%p;
}

inline int spr(int a){
	return mul(a,a);
}

ll gcd(ll a,ll b){
	return b?gcd(b,a%b):a;
}

ll lcm(ll a,ll b){
	return a/gcd(a,b)*b;
}

map<int,int>F;

int fib(int n){
	if(n<=1)return n;
	int &ans=F[n];
	if(ans)return ans;
	if(n&1)return ans=add(spr(fib((n+1)>>1)),spr(fib((n-1)>>1)));
	int tmp=fib(n>>1);
	return ans=mul(tmp,add(mul(2,fib((n>>1)-1)),tmp));
}

vector<int> d;

inline void div(int x){
	int t=sqrt(x);
	d.clear();
	for(int i=2;i<t;++i){
		if(x%i==0){
			d.push_back(i);
			d.push_back(x/i);
		}
	}
	if(t*t==x)d.push_back(t);
	sort(d.begin(), d.end());
}

inline int check(int x){
	div(x);
	F.clear();
	for(int i=0;i<d.size();++i){
		if(fib(d[i])==0&&fib(d[i]+1)==1)return d[i];
	}
//	assert(fib(x)==0&&fib(x+1)==1);
	return x;

}

inline int query(int x){
	if(x==2)return 3;
	if(x==3)return 8;
	if(x==5)return 20;
	p=x;
	if(x%5==1||x%5==4)return check(x-1);
	if(x%5==2||x%5==3)return check(2*x+2);
//	assert(0);
}

inline ll solve(int n){
	ll ans=1;
	for(int i=1;i<=tot;++i){
		if(pri[i]>n)break;
		if(n%pri[i]==0){
			n/=pri[i];
			ll t=1;
			while(n%pri[i]==0){
				n/=pri[i];
				t*=pri[i];
			}
			ans=lcm(ans,t*query(pri[i]));
		}
	}
	if(n>1){
		ans=lcm(ans,query(n));
	}
	return ans;
}

int main(){
    init();
    int T;
    for(r(T);T;--T){
    	int n;r(n);
    	printf("%lld
",solve(n));
    }
}

后记&其它做法

毕老师：“其实我只是想考察大家的找规律技巧，noip不会考证明的”
我：“其实我也只是自己感兴趣而已”

如果感兴趣可以看看这个网站以及它推荐的其它网站
OEIS也还行
我觉得在网上花些时间认真看论文还是有收获的,遗憾的是几乎没有中文资料,所以自己来写一写
学信息学竞赛就是这点好,有足够的时间和充足的资料(虽然大部分是英语的)可以自己去研究,而不是被别人牵着走
说实话我觉得凭兴趣和爱好去学习是一件很幸福的事,即使绕弯路,碰壁都会有收获,都会感到快乐

另外这道题还可以用BSGS水过去
虽然没有仔细看证明,有结论是答案不会超过(6*p)
然后大概就像这道题一样搞搞就好了