本文版权归ljh2000和博客园共有，欢迎转载，但须保留此声明，并给出原文链接，谢谢合作。
本文作者：ljh2000
作者博客：http://www.cnblogs.com/ljh2000-jump/
转载请注明出处，侵权必究，保留最终解释权！
Description
农夫栋栋近年收入不景气，正在他发愁如何能多赚点钱时，他听到隔壁的小朋友在讨论兔子繁殖的问题。
问题是这样的：第一个月初有一对刚出生的小兔子，经过两个月长大后，这对兔子从第三个月开始，每个月初生一对小兔子。新出生的小兔子生长两个月后又能每个月生出一对小兔子。问第n个月有多少只兔子？
聪明的你可能已经发现，第n个月的兔子数正好是第n个Fibonacci(斐波那契)数。栋栋不懂什么是Fibonacci数，但他也发现了规律：第i+2个月的兔子数等于第i个月的兔子数加上第i+1个月的兔子数。前几个月的兔子数依次为：
1 1 2 3 5 8 13 21 34 …
栋栋发现越到后面兔子数增长的越快，期待养兔子一定能赚大钱，于是栋栋在第一个月初买了一对小兔子开始饲养。
每天，栋栋都要给兔子们喂食，兔子们吃食时非常特别，总是每k对兔子围成一圈，最后剩下的不足k对的围成一圈，由于兔子特别害怕孤独，从第三个月开始，如果吃食时围成某一个圈的只有一对兔子，这对兔子就会很快死掉。
我们假设死去的总是刚出生的兔子，那么每个月的兔子数仍然是可以计算的。例如，当k=7时，前几个月的兔子数依次为：
1 1 2 3 5 7 12 19 31 49 80 …
给定n，你能帮助栋栋计算第n个月他有多少对兔子么？由于答案可能非常大，你只需要告诉栋栋第n个月的兔子对数除p的余数即可。
Input
输入一行，包含三个正整数n, k, p。
Output
输出一行，包含一个整数，表示栋栋第n个月的兔子对数除p的余数。
Sample Input
6 7 100
Sample Output
7
HINT
1<=N<=10^18
Source
5,5,3,0,
3,3,6,2,0,
2,2,4,6,3,2,5,0,5,5,3,0,
3,3,6,2,0,
是不是很有规律？
会发现每行第一个数是上一行最后一个非0数（这不是废话吗，0+任何数=自己啊...），而这一行的第i个数=第一个数×斐波那契数列的第i个数取模即可得到。
这就在暗示我们或许可以做了！令这一行的长度为len，那么第一个数×fib[len]=1(mod k)，我们可以通过已知得到fib[len] 的值，假如我们预处理出斐波那契数列在模k意义下的值，就可以做出每个模第一次出现的位置，就从而得到长度。所以只要对第一个数求个逆元就可以了，但是模数可能不是质数...那就上exgcd！当然如果没有逆元，意味着可以直接快速幂结束了...
但是这似乎还是会T，观察到上述序列很有可能出现循环节！就比如样例中的这个序列，我们没有必要反复计算，可以考虑一个循环节视为一个整体，对这个整体矩乘快速幂，剩下的零头再单独做，就可以完美解决这道题了。
1 //It is made by ljh2000
2 #include <iostream> 3 #include <cstdlib> 4 #include <cstring> 5 #include <cstdio> 6 #include <cmath> 7 #include <algorithm> 8 #include <ctime> 9 #include <vector> 10 #include <queue> 11 #include <map> 12 #include <set> 13 using namespace std; 14 typedef long long LL; 15 const int inf = (1<<30); 16 const int MAXN = 1000011; 17 const int MAXL = 6000011; 18 LL n,fir; 19 int k,MOD; 20 int fib[MAXL],vis[MAXN],len[MAXN]; 21 LL F[MAXN];//每个开头的数对应的使得其最后一个数为1的斐波那契数 22 bool hav[MAXN];//这个k意义下的模数作为开头是否出现过 23 struct matrix{ 24 int n,m; 25 LL s[3][3]; 26 }ans,ini,dan,mul,jian,I,old[MAXN]; 27 28 inline matrix operator * (matrix a,matrix b){ 29 matrix tmp=ini; tmp.n=a.n; tmp.m=b.m; 30 for(int i=0;i<a.n;i++) 31 for(int j=0;j<b.m;j++) 32 for(int l=0;l<a.m;l++) 33 tmp.s[i][j]+=a.s[i][l]*b.s[l][j],tmp.s[i][j]%=MOD; 34 return tmp; 35 } 36 37 inline int getint() 38 { 39 int w=0,q=0; char c=getchar(); 40 while((c<'0' || c>'9') && c!='-') c=getchar(); if(c=='-') q=1,c=getchar(); 41 while (c>='0' && c<='9') w=w*10+c-'0', c=getchar(); return q ? -w : w; 42 } 43 inline LL gcd(LL x,LL y){ if(y==0) return x; return gcd(y,x%y); } 44 inline matrix fast_pow(matrix a,LL y){ matrix tmp=dan; while(y>0) { if(y&1) tmp=tmp*a; a=a*a; y>>=1; } return tmp; } 45 inline LL exgcd(LL aa,LL bb,LL &x,LL &y){ 46 if(bb==0) { 47 x=1; y=0; 48 return aa; 49 } 50 LL cun=exgcd(bb,aa%bb,x,y),lin=x; 51 x=y; y=lin-(aa/bb)*y; 52 return cun; 53 } 54 inline void out(){ if(n!=0) ans=ans*fast_pow(mul,n); ans.s[0][1]%=MOD; ans.s[0][1]+=MOD; ans.s[0][1]%=MOD; printf("%lld",ans.s[0][1]); exit(0); } 55 inline void work(){ 56 scanf("%lld",&n); k=getint(); MOD=getint(); fib[1]=fib[2]=1; 57 for(int i=3;;i++) {//斐波那契数列的模k意义下的循环节不超过6*k 58 fib[i]=fib[i-1]+fib[i-2]; fib[i]%=k; 59 if(!vis[fib[i]]) vis[fib[i]]=i; 60 if(fib[i-1]==1 && fib[i]==1) break; 61 } 62 for(int i=0;i<3;i++) for(int j=0;j<3;j++) ini.s[i][j]=0; for(int i=0;i<3;i++) dan.s[i][i]=1; jian=dan; ans.s[0][0]=ans.s[0][2]=1; 63 ans.n=1; ans.m=3; mul.n=mul.m=jian.n=jian.m=dan.n=dan.m=3; bool FFF=false; 64 mul.s[0][1]=mul.s[1][0]=mul.s[1][1]=mul.s[2][2]=1; jian.s[2][1]=-1; fir=1; LL xx,yy,gg; 65 for(;n;) { 66 if(F[fir]==0) { 67 gg=gcd(fir,k); if(gg!=1) F[fir]=-1; 68 else { exgcd(fir,k,xx,yy); xx%=k; xx+=k; xx%=k;/*exgcd!!!*/ F[fir]=xx; } 69 } 70 if(F[fir]==-1) out(); 71 if(!hav[fir] || FFF) {//未出现过 72 if(!vis[F[fir]]) out();//不存在这种模数 73 len[fir]=vis[F[fir]];//出现的第一个位置 74 if(n>=len[fir]) { 75 old[fir]=fast_pow(mul,len[fir]); 76 n-=len[fir]; 77 } 78 else out(); 79 old[fir]=old[fir]*jian; 80 ans=ans*old[fir]; 81 hav[fir]=1; fir*=fib[len[fir]-1]; fir%=k; 82 } 83 else{//出现循环节 84 LL ff=fir; LL cnt=0/*!!!*/;//统计计算的次数 85 I=dan; 86 for(ff=fir*fib[len[fir]-1]%k;ff!=fir;(ff*=fib[len[ff]-1])%=k){ 87 I=I*old[ff]; cnt+=len[ff]; 88 } 89 cnt+=len[ff]; I=old[ff]*I;/*矩乘不满足交换律！！！*///I=I*old[ff]; 90 ans=ans*fast_pow(I,n/cnt); 91 n=n%cnt; FFF=true;//标记 92 } 93 } 94 out(); 95 } 96 97 int main() 98 { 99 work(); 100 return 0; 101 }