矩阵乘+坑 [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 long using 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(){;}