BZOJ 2432 [Noi2011]兔农 矩乘+数论

Posted 2020-10-01 TS_Hugh

这道题的暴力分还是很良心嘛~~~~~

直接刚的话我发现本蒟蒻只会暴力，矩乘根本写不出来，然后让我们找一下规律，我们发现如果我们把这个序列在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;
}