指数循环节:
a
n
m
o
d
p
=
a
n
%
ϕ
(
p
)
+
ϕ
(
p
)
m
o
d
p
(
n
≥
ϕ
(
p
)
)
a^nmod\ p=a^{n\%\phi(p)+\phi(p)}mod\ p \ (n\geq\phi(p))
anmod p=an%ϕ(p)+ϕ(p)mod p (n≥ϕ(p))
斐波那契模p循环节:
对于一个正整数n,我们求Fib数模n的循环节的长度的方法如下
(1)把n素因子分解,即
n
=
p
1
a
1
p
2
a
2
.
.
.
p
k
a
k
n=p_1^{a_1}p_2^{a_2}...p_k^{a_k}
n=p1a1p2a2...pkak
(2)分别计算Fib数模每个的循环节长度,假设长度分别是
x
1
,
x
2
,
.
.
.
,
x
k
x_1,x_2,...,x_k
x1,x2,...,xk
(3)那么Fib模n的循环节长度
l
c
m
(
x
1
,
x
2
,
.
.
.
,
x
k
)
lcm(x_1,x_2,...,x_k)
lcm(x1,x2,...,xk)
Fib数模 p m p^m pm的最小循环节长度等于 G ( p ) ∗ p m − 1 G(p)*p^{m-1} G(p)∗pm−1,其中 G ( p ) G(p) G(p)表示Fib数模素数p的最小循环节长度。
对于求 G ( p ) G(p) G(p)我们利用如下定理:
如果5是模p的二次剩余(即是 5 p − 1 2 % = 1 5^{\frac{p-1}{2}}\%=1 52p−1%=1),那么循环节的的长度是的 p − 1 p-1 p−1因子,否则,循环节的长度是 2 ∗ ( p − 1 ) 2*(p-1) 2∗(p−1)的因子。
对于小于等于5的素数,我们直接特殊判断, G ( 2 ) = 3 , G ( 3 ) = 8 , G ( 5 ) = 20 G(2)=3,G(3)=8,G(5)=20 G(2)=3,G(3)=8,G(5)=20。
那么我们可以先求出所有的因子,然后用矩阵快速幂来一个一个判断.
#include <iostream>
#include <string.h>
#include <algorithm>
#include <stdio.h>
#include <math.h>
using namespace std;
typedef unsigned long long LL;
const int M = 2;
struct Matrix
{
LL m[M][M];
};
Matrix A;
Matrix I = {1,0,0,1};
Matrix multi(Matrix a,Matrix b,LL MOD)
{
Matrix c;
for(int i=0; i<M; i++)
{
for(int j=0; j<M; j++)
{
c.m[i][j] = 0;
for(int k=0; k<M; k++)
c.m[i][j] = (c.m[i][j]%MOD + (a.m[i][k]%MOD)*(b.m[k][j]%MOD)%MOD)%MOD;
c.m[i][j] %= MOD;
}
}
return c;
}
Matrix power(Matrix a,LL k,LL MOD)
{
Matrix ans = I,p = a;
while(k)
{
if(k & 1)
{
ans = multi(ans,p,MOD);
k--;
}
k >>= 1;
p = multi(p,p,MOD);
}
return ans;
}
LL gcd(LL a,LL b)
{
return b? gcd(b,a%b):a;
}
const int N = 400005;
const int NN = 5005;
LL num[NN],pri[NN];
LL fac[NN];
int cnt,c;
bool prime[N];
int p[N];
int k;
void isprime()
{
k = 0;
memset(prime,true,sizeof(prime));
for(int i=2; i<N; i++)
{
if(prime[i])
{
p[k++] = i;
for(int j=i+i; j<N; j+=i)
prime[j] = false;
}
}
}
LL quick_mod(LL a,LL b,LL m)
{
LL ans = 1;
a %= m;
while(b)
{
if(b & 1)
{
ans = ans * a % m;
b--;
}
b >>= 1;
a = a * a % m;
}
return ans;
}
LL legendre(LL a,LL p)
{
if(quick_mod(a,(p-1)>>1,p)==1) return 1;
else return -1;
}
void Solve(LL n,LL pri[],LL num[])
{
cnt = 0;
LL t = (LL)sqrt(1.0*n);
for(int i=0; p[i]<=t; i++)
{
if(n%p[i]==0)
{
int a = 0;
pri[cnt] = p[i];
while(n%p[i]==0)
{
a++;
n /= p[i];
}
num[cnt] = a;
cnt++;
}
}
if(n > 1)
{
pri[cnt] = n;
num[cnt] = 1;
cnt++;
}
}
void Work(LL n)
{
c = 0;
LL t = (LL)sqrt(1.0*n);
for(int i=1; i<=t; i++)
{
if(n % i == 0)
{
if(i * i == n) fac[c++] = i;
else
{
fac[c++] = i;
fac[c++] = n / i;
}
}
}
}
LL find_loop(LL n)
{
Solve(n,pri,num);
LL ans=1;
for(int i=0; i<cnt; i++)
{
LL record=1;
if(pri[i]==2)
record=3;
else if(pri[i]==3)
record=8;
else if(pri[i]==5)
record=20;
else
{
if(legendre(5,pri[i])==1)
Work(pri[i]-1);
else
Work(2*(pri[i]+1));
sort(fac,fac+c);
for(int k=0; k<c; k++)
{
Matrix a = power(A,fac[k]-1,pri[i]);
LL x = (a.m[0][0]%pri[i]+a.m[0][1]%pri[i])%pri[i];
LL y = (a.m[1][0]%pri[i]+a.m[1][1]%pri[i])%pri[i];
if(x==1 && y==0)
{
record = fac[k];
break;
}
}
}
for(int k=1; k<num[i]; k++)
record *= pri[i];
ans = ans/gcd(ans,record)*record;
}
return ans;
}
void Init()
{
A.m[0][0] = 1;
A.m[0][1] = 1;
A.m[1][0] = 1;
A.m[1][1] = 0;
}
int main()
{
LL n;
Init();
isprime();
while(cin>>n)
cout<<find_loop(n)<<endl;
return 0;
}