斐波那契循环节 笔记

发表于 | 更新于 | 阅读次数：
本文字数： 2.2k 字 | 阅读时长 ≈ 10 分钟

这篇是比较具体的吧…很多证明网络上别的博客是没有的，只有一个结论。还有一篇博客用了很多复杂的东西里证明，虽然我能勉强看懂，但是很多人珂能看不懂。这篇也许算一个适中的了，容易理解些，但是证明少些。

如果您有更好的证明，或者能补全我没有的证明，也很感谢了。

正片开始

关于斐波那契数列，大家应该都不陌生。$f_1=1,f_2=1,f_n=f_{n-1}+f_{n-2}$。斐波那契循环节，就是求斐波那契数列膜一个数的循环节长度。假设这个膜数为$m$。那么以下的式子，如果没有说明，都是膜$m$意义下的。

一些简单的变换

设斐波那契数列膜$m$的循环节长度为$l(m)$。形式化地，$m$为最小的满足$f_i=0,f_{i+1}=1$的正整数。而且，如果存在$n$满足$f_n=0,f_{n+1}=1$，那么$n$是$m$的正整数倍。我们管这个叫珂朵莉的定义式。

我们把$m$质因数分解成$a_1^{p_1}a_2^{p_2}a_3…a_k^{p_k}$的形式。对于一个$a^p$，满足$l(a^p)$=$l(a)\times a^{p-1}$。然后$l(m)$，是每个$l(a^p)$的最小公倍数的因数。我们管这个叫珂朵莉第一定理。

那么接下来的部分，我们都是在求膜一个质数的循环节长度。

膜质数$p$的循环节长度

特殊说明来了：以下的式子，如果没有第二个特殊说明，都是膜$p$意义下的。

我们知道斐波那契数列有通项公式，你可以用特征方程，或者百度，必应和谷歌等方式证明这个公式：
设$\sqrt{5}=q$，则$f_i=\dfrac{1}{q}[(\dfrac{1+q}{2})^i-(\dfrac{1-q}{2})^i]$。

分类讨论$5$模$p$的二次剩余情况。

0. 特判

如果$p=2,3,5$，直接判掉，因为后面的式子仿佛都要$p>5$。顺便说一句，此时$l(p)$的值分别是$3,8,20$。

1. $5$膜$p$有二次剩余 （p=1或4 膜5）

此时，$q$是存在的，满足$q^2=5 (mod p)$。然后，根据费马小定理，循环节的长度就一定是$p-1$的倍数，具体是哪个，就暴力验证一下即珂。因为因数最多$\sqrt{p}$个，一次验证（矩阵快速幂）是$O(log)$的，顶多$O(\sqrt{p}logp)$。一般这也就够了…

2. $5$膜$p$没有二次剩余 (p=2或3 膜5)

由勒让德记号（自行百度）的特征得知，$5^{\dfrac{p-1}{2}}=-1$。

考虑$({1+q})^p$这个式子。由于$p$是质数，所以$C_p^i (0<i<p)=\dfrac{p(p-1)(p-2)…(p-i+1)}{i!}$中，$i!$和$p$肯定是互质的，所以最后的结果肯定等于$p$乘以一大堆式子。又因为这个组合数肯定是整数，所以后面一大堆式子也就是整数了。换句话说，$C_p^i$肯定是$p$的倍数。

然后我们用二项式定理拆开$(1+q)^p$这个式子，中间有一些项，膜$p$肯定都是$0$，只剩下$1+q^p=1+5^{\dfrac{p-1}{2}}q$。

然后因为$5^{\dfrac{p-1}{2}}=-1$，所以$(1+q)^p=1-q$。

又因为$p>5$，所以显然满足$2^{(p-1)}=1$

所以$(\dfrac{1+q}{2})^p=(\dfrac{1-q}{2})$。

带入$i=2p+2$到$f_i$，直接上通项公式。把$(\dfrac{1+q}{2})^{2p+2}$代换一下，变成$(\dfrac{1+q}{2})^2(\dfrac{1-q}{2})^2$。同理，后面的$(\dfrac{1-q}{2})^{2p+2}$变成$(\dfrac{1-q}{2})^2(\dfrac{1+q}{2})^2$。相等，然后把它们相减，为$0$。

然后我们带入$i=2p+3$进去，自己推一下，它等于$1$。

由珂朵莉的定义式知，此时循环节长度是$2p+2$的倍数。

这上面三个东西，我们得到一些结论：

  1. p=2,3,5,l(p)=3,8,20
  2. p=1,4 (mod 5)，l(p)|p-1
  3. p=2,3 (mod 5)，l(p)|2p+2 （珂朵莉：这个暴力算就好了）

我们管这个叫珂朵莉第二定理。

总结

求循环节的步骤：

质因数分解

求每个质因子的循环节，然后用珂朵莉第一定理合并答案

对于每个质因子，用珂朵莉第二定理求出循环节长度

然后你就求出了循环节。求出了循环节之后，您就珂以做一些奇妙的题目，比如这个：

求$f_i\%p$，其中$i<=10^{30000000}$，$p<=1e9$。（洛谷4000）

（先求循环节，然后直接算）

这个的代码（用了$vector$，比较好理解，但是常数巨大）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
#include <bits/stdc++.h>
using namespace std;
namespace Flandre_Scarlet
{
    #define N 30000007
    #define int long long 
    #define F(i,l,r) for(int i=l;i<=r;++i)
    #define D(i,r,l) for(int i=r;i>=l;--i)
    #define Fs(i,l,r,c) for(int i=l;i<=r;c)
    #define Ds(i,r,l,c) for(int i=r;i>=l;c)
    #define Tra(i,u) for(int i=G.Start(u),__v=G.To(i);~i;i=G.Next(i),__v=G.To(i))
    #define MEM(x,a) memset(x,a,sizeof(x))
    #define FK(x) MEM(x,0)

    class Matrix//这是一个矩阵，这个实现自己写就好了，不需要阅读
    {
        #define NN 5
        private:
            int a[NN][NN];
        public:
            int n;
            Matrix()
            {
                memset(a,0,sizeof(a));
                n=0;
            }
            Matrix(int _n)
            {
                memset(a,0,sizeof(a));
                n=_n;
            }
            Matrix(int _n,int _x)
            {
                n=_n;
                for(int i=0;i<NN;++i)
                {
                    for(int j=0;j<NN;++j)
                    {
                        a[i][j]=_x;
                    }
                }
            }
    
            int* operator[](int i)
            {
                return *(a+i);
            }
    
            void Set(int x)
            {
                for(int i=0;i<NN;++i)
                {
                    for(int j=0;j<NN;++j)
                    {
                        a[i][j]=x;
                    }
                }
            }
            void Identity()
            {
                memset(a,0,sizeof(a));
                for(int i=0;i<NN;++i)
                {
                    a[i][i]=1;
                }
            }
            #undef NN //5
    };
    Matrix Times(Matrix x,Matrix y,int mod) //两个矩阵%mod的积
    {
        Matrix ans(x.n,0);
        int n=ans.n;
        for(int i=1;i<=n;++i)
        {
            for(int j=1;j<=n;++j)
            {
                if (x[i][j])
                {
                    for(int k=1;k<=n;++k)
                    {
                        ans[i][k]+=(x[i][j]*y[j][k])%mod;
                        ans[i][k]%=mod;
                    }
                }
            }
        }
        return ans;
    }
    Matrix Qpow(Matrix x,int p,int mod) //x^p %mod
    {
        Matrix ans(x.n,1);
        ans.Identity();
        while(p)
        {
            if (p&1) ans=Times(ans,x,mod);
            x=Times(x,x,mod),p>>=1;
        }
        return ans;
    } //以上为矩阵
    int Fib(int x,int mod) //求f[x] % mod
    {
        Matrix Trans(3,0);
        Trans[1][1]=Trans[1][2]=Trans[2][1]=1;

        if (x==0) return 0;
        if (x==1 or x==2) return 1;
        Matrix Ans=Qpow(Trans,x-1,mod);
        return Ans[1][1];
    }

    int p;char a[N];
    void Input()
    {
        scanf("%s%lld",a,&p);
    }

    struct node{int fac,pow,sum;}; //这是一个因式，它等于fac^pow，sum记录fac^pow的值，方便后面求fac^(pow-1) （珂朵莉第一定理的式子）
    vector<node> Factor(int x) // 把x分解成很多个node
    {
        vector<node> ans;ans.clear();
        for(int i=2;i*i<=x;++i)
        {
            if (x%i==0)
            {
                node cur;cur.fac=i;cur.pow=0;cur.sum=1;
                while(x%i==0)
                {
                    cur.sum*=i;
                    x/=i;
                    ++cur.pow;
                }
                ans.push_back(cur);
            }
        }
        if (x>1)
        {
            ans.push_back((node){x,1,x});
        }
        return ans;
    }
    int Get1(int p) //珂朵莉第二定理求f[i]%p的值
    {
        if (p==2) return 3;
        if (p==3) return 8;
        if (p==5) return 20; // 特判

        int len; //l(p)是len的因数
        if (p%5==1 or p%5==4) len=p-1; //有二次剩余的情况
        else len=2*(p+1); //没有

        // 你在这里直接return len也不是没有问题，因为这题你只要求斐波那契数列的值。至于你求的是最小的循环节，还是它的某一个倍数，这不重要

        vector<int> fac; //保存因数
        for(int i=1;i*i<=len;++i)
        {
            if (len%i==0)
            {
                fac.push_back(i);
                if (i*i!=len) fac.push_back(len/i); 
            }
        } 
        sort(fac.begin(),fac.end()); //排个序，以便求最小的
        F(i,0,(int)fac.size()-1)
        {int v=fac[i];
            if (Fib(v,p)==0 and Fib(v-1,p)==1)
            {
                return v; //找到第一个就返回吧
            }
        }
        return len;//这句没什么用...
    }
    int lcm(int a,int b){return a*b/__gcd(a,b);}
    int GetLen(int x) //用珂朵莉第一定理合并
    {
        vector<node> fac=Factor(x);
        int ans=1;
        F(i,0,(int)fac.size()-1)
        {node u=fac[i];
            int cur=Get1(u.fac)*u.sum/u.fac;
            ans=lcm(ans,cur);
        }
        return ans; //严格上来说，l(x)应该是ans的因数，而不一定是ans。但是，和上面类似的道理，我懒的写了，直接返回ans。
    }
    void Soviet()
    {
        int mod=GetLen(p); //求循环节
        int aa=0,la=strlen(a);
        F(i,0,la-1)
        {
            aa=(aa*10+a[i]-'0')%mod;
        } //取一下膜
        printf("%lld\n",Fib(aa,p)); //然后直接做即珂
    }

    #define Flan void
    Flan IsMyWife()
    {
        Input();
        Soviet();
    }
    #undef int //long long 
}
int main(){
    Flandre_Scarlet::IsMyWife();
    getchar();getchar();
    return 0;
}

w