矩阵乘+坑 [NOI2011] 兔农

     

问题 C: 兔农

时间限制: 1 Sec  内存限制: 256 MB

题目描述

农夫栋栋近年收入不景气，正在他发愁如何能多赚点钱时，他听到隔壁的小 朋友在讨论兔子繁殖的问题。 问题是这样的：第一个月初有一对刚出生的小兔子，经过两个月长大后，这 对兔子从第三个月开始，每个月初生一对小兔子。新出生的小兔子生长两个月后 又能每个月生出一对小兔子。问第 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 的余数即可。

输入

输入一行，包含三个正整数 n, k, p。

输出

输出一行，包含一个整数，表示栋栋第 n 个月的兔子对数除 p 的余数。

样例输入

6 7 100

样例输出

7

提示

    好一道国赛的题，觉得巨坑很好。。

    菲波数%k所得到的数排成排，在需要-1（也就是%k==1）时，重启一行，就会发现能搞出一个循环，也不一定，我下面说。。

    有两个结论，我不会证，但对这道题很有帮助，

     一：如果有循环节，最多找6*k个数，就一定能循环一次，

    二：设x为这一行行首的值，x*fib[i]%k==这一行第i位的值。

    然后发现，x*fib[i]%k==1就是行尾，而fib[i]即为x的乘法逆元，所以能求出fib[i]的值，

    并且如结论一所提，在这里找时，记录某个fib[i]第一次出现的位置，因此求出逆元后，就变相得知了这行的长度。当把循环节求出后，就好说了。。

    注：有可能出现每个fib[i]%k都无逆元的情况，就是说这一行无限延伸下去了。。直接求就行了。

    当然乘的过程中用矩阵乘优化，其实那个矩阵怎么构造的我也不太懂

    

#include<cstdio>#include<cstdlib>#include<cstring>#include<algorithm>#include<iostream>#define ll long longusing namespace std;ll read(){ll sum=0,f=1;char x=getchar();while(x<'0'||x>'9'){if(x=='-')f=-1;x=getchar();}while(x>='0'&&x<='9'){sum=sum*10+x-'0';x=getchar();}return sum*f;}ll f[6000405],vis[1000405],ine[1000405],don[1000405],len[1000405];ll n,p,k;bool flag=0;struct hh{ll n,m,f[3][3];hh(){}hh(int x,int y){n=x,m=y;memset(f,0,sizeof(f));}} mul,d,ans,res[1000505];inline hh operator * (hh a,hh b){hh c(a.n,b.m);for(int i=0;i<a.n;i++)   for(int j=0;j<b.m;j++)      for(int k=0;k<a.m;k++)          (c.f[i][j]+=(a.f[i][k]*b.f[k][j])%p)%=p;return c;}inline ll ny(ll a,ll b,ll c){if(a==0)return -1;else if((c%a)==0)return c/a;ll t=ny(b%a,a,((-c%a)+a)%a);if(t==-1)return -1;return (t*b+c)/a;}inline hh quick(hh a,ll x){hh sum(a.n,a.m);for(int i=0;i<a.n;i++)   sum.f[i][i]=1;for(;x;x>>=1,a=a*a)    if(x&1)sum=sum*a;return sum;}int yjn(){//freopen("noi2011_rabbit.in","r",stdin);//freopen("noi2011_rabbit.out","w",stdout);n=read();k=read();p=read();f[1]=f[2]=1;for(int i=3;i<=6*k;i++){f[i]=(f[i-1]+f[i-2])%k;if(!vis[f[i]])vis[f[i]]=i;if(f[i]==1&&f[i-1]==1)break;}mul.n=mul.m=d.n=d.m=3;ans.n=1;ans.m=3;ans.f[0][0]=ans.f[0][2]=1;mul.f[0][1]=mul.f[1][0]=mul.f[1][1]=mul.f[2][2]=1;d.f[0][0]=d.f[1][1]=d.f[2][2]=1;d.f[2][1]=-1;for(ll i=1;n;){if(!ine[i])ine[i]=ny(i,k,1);if(ine[i]==-1){    ans=ans*quick(mul,n);n=0;}else{if(!don[i]||flag){don[i]=1;if(!vis[ine[i]]){ans=ans*quick(mul,n);n=0;}else{len[i]=vis[ine[i]];if(n>=len[i]){n-=len[i];res[i]=quick(mul,len[i])*d;ans=ans*res[i];(i*=f[len[i]-1])%=k;}else{ans=ans*quick(mul,n);n=0;}}}else{ll cnt=0;hh l(3,3);l.f[0][0]=l.f[1][1]=l.f[2][2]=1;for(ll j=i*f[len[i]-1]%k;j!=i;(j*=f[len[j]-1])%=k)    cnt+=len[j],l=l*res[j];cnt+=len[i],l=res[i]*l;ans=ans*quick(l,n/cnt);n%=cnt;flag=1;}}}cout<<(ans.f[0][1]+p)%p;}int qty=yjn();int main(){;}