给出矩阵A，求矩阵\(A+A^2+...+A^k\)各个元素\(mod\ yyb\)的值，\(n\leq 30,k\leq 10^9,yyb\leq 10^4\)。

解

法一：分治

显然是数列题，故数列最浅显的减法是分治，寻找其中重复计算的部分，故可以

\[A+A^2+...+A^{k/2}+A^{k/2+1}+...+A^k=\]

\[A+A^2+...+A^{k/2}+A^{k/2}(A+...+A^{k/2})+(k\&1)A^k=\]

\[(A^{k/2}+1)(A+...+A^{k/2})+(k\&1)A^k\]

对式子主体进行分治，其他的部分快速幂即可。

参考代码：

#include 
    
     #include 
     
      #include 
      
       #define il inline#define ri registerusing namespace std;int n,yyb;struct matrix{    int jz[30][30];    il void clear(){        memset(jz,0,sizeof(jz));    }    il void unit(){        clear();ri int i;        for(i=0;i
       
            il matrix operator^(free y){        matrix x(*this),ans;ans.unit();        while(y){            if(y&1)ans=ans*x;            x=x*x,y>>=1;        }return ans;    }}state,unit,zero;il matrix efen(int);int main(){    int k,i;    scanf("%d%d%d",&n,&k,&yyb);    unit.unit(),state.read();    efen(k).print();    return 0;}il matrix efen(int y){    if(!y)return unit;    if(y==1)return state;    return ((state^(y>>1))+unit)*efen(y>>1)        +((y&1)?(state^y):zero);}
       
      
     
    

法二：矩阵快速幂

显然数列的题目，经常会存在递推方程，于是矩阵快速幂会在其中大有用武之地，于是设\(f_i=A+A^2+...+A^i\),不难有递推方程\(f_i=f_{i-1}+A^i\)，于是我们可以同时转移\(f_i,A^i\),故状态矩阵为

\[\begin{bmatrix}A_i&f_{i-1}\end{bmatrix}\]

转移矩阵为

\[\begin{bmatrix}A&I\\0&I\end{bmatrix}\]

以此矩阵中套矩阵仿照套路即可。

参考代码：

#include 
    
     #include 
     
      #include 
      
       #define il inline#define ri registerusing namespace std;int n,yyb;struct matrix1{    int jz[30][30];    il void read(){        ri int i,j;        for(i=0;i
       
            il matrix1 operator^(free y){        matrix1 ans,x(*this);ans.unit();        while(y){            if(y&1)ans=ans*x;            x=x*x,y>>=1;        }return ans;    }}A;struct matrix2{    matrix1 jz[2][2];    il void clear(){        memset(jz,0,sizeof(jz));    }    il void unit(){        clear();ri int i;        for(i=0;i<2;++i)jz[i][i].unit();    }    il matrix2 operator*(matrix2 x){        matrix2 y;y.clear();        ri int i,j,k;        for(i=0;i<2;++i)            for(j=0;j<2;++j)                for(k=0;k<2;++k)                    y.jz[i][j]=y.jz[i][j]+jz[i][k]*x.jz[k][j];        return y;    }template
        
             il matrix2 operator^(free y){        matrix2 ans,x(*this);ans.unit();        while(y){            if(y&1)ans=ans*x;            x=x*x,y>>=1;        }return ans;    }}state,tran;int main(){    int k;    scanf("%d%d%d",&n,&k,&yyb);    A.read(),state.jz[0][0]=A;    tran.jz[0][0]=A,tran.jz[0][1].unit();    tran.jz[1][1].unit(),state=state*(tran^k);    state.jz[0][1].print();    return 0;}
        
       
      
     
    

法三：倍增前缀和优化

数列求和可以利用倍增优化，主要思想是维护\(A^{2^i}\)，这个显然可通过\(A^{2^{i+1}}=(A^{2^i})^2\)来暴力递推，再维护一个倍增的数列，\(f_i=A^1+A^2+...+A^{2^i}\)，不难得知其转移方程为\(f_i=f_{i-1}A^{2^{i-1}}+f_{i-1}\)，而这个的维护也可以暴力维护，于是对于我们的求和，以\(A+A^2+...+A^{15}\)为例，有

\[A+A^2+...+A^{15}=f_3+A^9+...+A^{15}=\]

\[f_3+A^8(A+...+A^7)=f_3+A^8(f_2+A^5+A^6+A^7)=\]

\[f_3+A^8[f_2+A^4(A+A^2+A^3)]=\]

\[f_3+A^8[f_2+A^4(f_1+A^2A^1)]=f_3+A^8[f_2+A^4(f_1+A^2f_0)]\]

于是按照这个先维护好一段倍增前缀和，再对后式提公因式，继续利用倍增前缀和优化，不停地继续，即可得到答案，当然你也可以递归处理，以下参考程序用的是非递归的方式。

参考代码：

#include 
    
     #include 
     
      #include 
      
       #define il inline#define ri registerusing namespace std;int n,yyb;struct matrix{    int jz[30][30];    il void clear(){        memset(jz,0,sizeof(jz));    }    il void unit(){        clear();ri int i;        for(i=0;i
       
        =0;--i)        if(k>>i&1){            ans=ans+sum[i]*jilu;            jilu=jilu*p[i];        }ans.print();    return 0;}
       
      
     
    

总上所诉，不难得知数列题目的一般求法，矩阵快速幂和分治，而特殊地对于数列前缀和的问题，我们可以利用倍增前缀和优化，但无论如何，数列问题绝不只一条道路通往罗马。