首页 > 代码库 > 【BZOJ 2432】 [Noi2011]兔农 矩乘+数论
【BZOJ 2432】 [Noi2011]兔农 矩乘+数论
这道题的暴力分还是很良心嘛~~~~~
直接刚的话我发现本蒟蒻只会暴力,矩乘根本写不出来,然后让我们找一下规律,我们发现如果我们把这个序列在mod k的意义下摆出,并且在此过程中把值为1的的数减一,我们发现他可以成为一段一段的被0(我们在此只关注减1变为0的点)区间,我们继续分析我们分析出来了这样的性质:如果存在这样的点,那么他右边的点一定是两个重复的数,而且往后是fibonacci数列(重头开始)乘第一个数,那么他之后再出现这样的0,的充要条件是其后存在一个fibonacci数是这段数第一个数的逆元。
我们先介绍一个结论(本蒟蒻并不会证):斐波那契数列模k后一定是0,1,1开头的纯循环,而且这个循环节的长度≤6k。我们开一个数组vis[i]表示第一个在mod k意义下值为i的fibonacci数的位置,这样我们求出某个数的逆元就知道以这个数为首项的数段的长度了(这样我们就可以跳啦)。
现在我说一下到了现在我们做了什么,我们发现在mod k 的意义下,这个数列会是一个由0(再次强调我们在此只关注减1变为0的点)隔开散区间(一整个也算)开头,之后要么成为循环,要么不再有0。先解释成为循环:因为我们在mod k的意义下因此最多不过k个数段就可以找到循环。我们再解释一下不在有0,这个就是当这个数段的第一项在mod k的意义下没有逆元,或者有逆元但是找不到vis[]。
现在我们用k再乘一个不大的数的时间复杂度找出来我们要处理的假fibonacci在哪些位置减一,以及他到底是存在循环节还是到后来没有了0,现在我们就可以想到一些矩阵,进行操作,但是这些矩阵一定要满足可乘。对于循环节我们暴力处理两边和一个循环节,对于最后没有0,我们就在后面直接fibonacci。
坑:I.via[1]!=1!!!!我们要找的是他后面第一个1,而我们前两个数是受法律保护的。
II.对于处理一段一段的,最后一段可能顶到n,不满
|||.一定要注意矩阵乘没有交换律
#include <cstdio> #include <cmath> #include <cstring> #include <iostream> using namespace std; typedef long long LL; const LL N=1000010; LL n,k,p,vis[N],f[N*6],Ola,Stop,L,R,pos[N]; LL a[4][4],b[4][4],temp_ab[4][4],F[4],temp_F[4],S[4][4]; bool huzhi[N]; LL GCD(LL x,LL y){ return x==0?y:GCD(y%x,x); } inline LL Min(LL x,LL y){ return x<y?x:y; } inline void get_Ola(){ LL lim=(LL)sqrt(k+0.5); LL x=k;Ola=k; for(LL i=2;i<=lim;i++) if(x%i==0){ Ola=Ola/i*(i-1); while(x%i==0)x/=i; } if(x!=1){ Ola=Ola/x*(x-1); } } //*********************Multi**********************// inline void Multi_One(){ memset(temp_F,0,sizeof(temp_F)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) temp_F[i]=(temp_F[i]+F[j]*a[j][i]%p+p)%p; memcpy(F,temp_F,sizeof(F)); } inline void Multi_Two(){ memset(temp_ab,0,sizeof(temp_ab)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) for(int l=1;l<=3;l++) temp_ab[i][j]=(temp_ab[i][j]+a[i][l]*a[l][j]%p+p)%p; memcpy(a,temp_ab,sizeof(a)); } inline void Multi_Three(){ memset(temp_F,0,sizeof(temp_F)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) temp_F[i]=((temp_F[i]+F[j]*b[j][i]%p)%p+p)%p; memcpy(F,temp_F,sizeof(F)); } inline void Multi_Four(){ memset(temp_ab,0,sizeof(temp_ab)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) for(int l=1;l<=3;l++) temp_ab[i][j]=(temp_ab[i][j]+S[i][l]*a[l][j]%p+p)%p; memcpy(S,temp_ab,sizeof(S)); } inline void Multi_Five(){ memset(temp_ab,0,sizeof(temp_ab)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) for(int l=1;l<=3;l++) temp_ab[i][j]=(temp_ab[i][j]+a[i][l]*a[l][j]%p)%p; memcpy(a,temp_ab,sizeof(a)); } inline void Multi_Six(){ memset(temp_ab,0,sizeof(temp_ab)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) for(int l=1;l<=3;l++) temp_ab[i][j]=(temp_ab[i][j]+S[i][l]*b[l][j]%p+p)%p; memcpy(S,temp_ab,sizeof(S)); } inline void Multi_Seven(){ memset(temp_F,0,sizeof(temp_F)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) temp_F[i]=(temp_F[i]+F[j]*S[j][i]%p+p)%p; memcpy(F,temp_F,sizeof(F)); } inline void Multi_E(){ memset(temp_ab,0,sizeof(temp_ab)); for(int i=1;i<=3;i++) for(int j=1;j<=3;j++) for(int l=1;l<=3;l++) temp_ab[i][j]=(temp_ab[i][j]+S[i][l]*S[l][j]%p+p)%p; memcpy(S,temp_ab,sizeof(S)); } //*******************GCD&&Ola*********************// inline LL Pow(LL x,LL y,LL P){ LL ans=1; while(y){ if(y&1)ans=ans*x%P; y>>=1,x=x*x%P; } return ans; } inline void POW(int step){ while(step){ if(step&1)Multi_Seven(); step>>=1,Multi_E(); } } //**********************POW***********************// inline void Stop_Forever(){ LL now=1,step=1; while(1){ if(step>=n){ Stop=n; break; } if(pos[now]){ L=pos[now],R=step-1; break; } pos[now]=step; if(huzhi[now]==0){ Stop=step-1;break; } LL Now=Pow(now,Ola-1,k); if(vis[Now]==0){ Stop=step-1;break; } step+=vis[Now]; now=f[vis[Now]-1]*now%k; } } //*********************Judge*********************// inline void F_H(){ f[1]=1,f[2]=1; for(int i=3;i<=6*k;i++) f[i]=(f[i-1]+f[i-2])%k,vis[f[i]]=vis[f[i]]==0?i:vis[f[i]]; for(int i=1;i<k;i++) if(GCD(i,k)==1)huzhi[i]=1; } inline void Pre(){ b[1][1]=1; b[2][2]=1; b[3][2]=-1,b[3][3]=1; F[1]=0,F[2]=F[3]=1; } inline void Init(){ scanf("%lld%lld%lld",&n,&k,&p); F_H(); get_Ola(); Stop_Forever(); Pre(); } //***********************Pre**********************// inline void GO(LL step){ memset(a,0,sizeof(a)); a[1][2]=1; a[2][1]=1,a[2][2]=1; a[3][3]=1; while(step){ if(step&1)Multi_One(); step>>=1,Multi_Two(); } } inline void GO_ON(int step){ memset(a,0,sizeof(a)); a[1][2]=1; a[2][1]=1,a[2][2]=1; a[3][3]=1; while(step){ if(step&1)Multi_Four(); step>>=1,Multi_Five(); } } inline void Get_Boss(LL &now,LL &step){ S[1][1]=1; S[2][2]=1; S[3][3]=1; while(step<R){ LL Now=Pow(now,Ola-1,k); LL y=vis[Now]; if(!step)y--; GO_ON(y); if(!step)y++; step+=y; Multi_Six(); now=f[vis[Now]-1]*now%k; } } //****************PPRREE******************// inline void Work_Stop(){ LL now=1,step=0; while(step<Stop){ LL Now=Pow(now,Ola-1,k); LL y=Min(n-step,vis[Now]); if(!step)y--; GO(y); if(!step)y++; step+=y; if(y==vis[Now])Multi_Three(); now=f[vis[Now]-1]*now%k; } GO(n-step); printf("%lld",F[2]); } inline void Work_Forever(){ LL now=1,step=0; while(step<L-1){ LL Now=Pow(now,Ola-1,k); LL y=vis[Now]; if(!step)y--; GO(y); if(!step)y++; step+=y; Multi_Three(); now=f[vis[Now]-1]*now%k; } Get_Boss(now,step); POW((n-L+1)/(R-L+1)); step=(n-L+1)/(R-L+1)*(R-L+1)+L-1; while(step<n){ LL Now=Pow(now,Ola-1,k); LL y=Min(n-step,vis[Now]); GO(y); step+=y; if(y==vis[Now])Multi_Three(); now=f[vis[Now]-1]*now%k; } printf("%lld",F[2]); } //*********************War*****************// int main(){ Init(); if(Stop)Work_Stop(); else Work_Forever(); return 0; }
【BZOJ 2432】 [Noi2011]兔农 矩乘+数论