POJ 3150 Cellular Automaton（矩阵快速幂）

Description
给出一个n个元素组成的环状序列，现在对这个序列进行k次同样的操作，每次操作，从第一个元素开始到最后一个元素结束，将每个元素的值更新为这个元素周围d范围内所有元素(2*d+1个)的和模m后的值，输出操作完后的序列
Input
第一行为四个整数n,m,d,k,第二行为n个整数表示这个环状序列
Output
输出操作完后的序列
Sample Input
sample input #1
5 3 1 1
1 2 2 1 2
sample input #2
5 3 1 10
1 2 2 1 2
Sample Output
sample output #1
2 2 2 2 1
sample output #2
2 0 0 2 2
Solution
显然是要构造转移矩阵，首先将原序列序列放入一个1*n矩阵
B=
下面来看转移矩阵的构造，以n=7,d=2为例，转移矩阵如下

轻易得出这个转移矩阵的每一行都是由上一行得到，即有转移方程A[i][j]=A[(i-1+n)%n][(j-1+n)%n]，这种性质十分有用，对于最初转移矩阵的构造遗迹矩阵快速幂时每次矩阵的乘法，我们都只需算出第一行的数据，之后递推求出后面n-1行即可，这样可以在矩阵乘法部分把时间复杂度由O(n^3)降到O(n^2)，下面来证明，两个具有这种性质的矩阵相乘之后仍然具有这种性质，

所以这种性质是具有积性的，故构造完转移矩阵A后用矩阵快速幂算出A^k然后用B*A^k即为答案，总时间复杂度为O(n^2*log k)
注意这道题要开的二维数组非常大，在函数中没法开，只能开全局，所以函数中要用指针传递数组
Code

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
#define maxn 520
typedef long long ll;
int n,m,d,k;
ll num[maxn];
ll temp[maxn][maxn],ans[maxn][maxn],f[maxn][maxn];
void mod_mul(ll a[][maxn],ll b[][maxn])//矩阵乘法,结果放在b数组中 
{
    memset(temp,0,sizeof(temp));//初始化 
    for(int k=0;k<n;k++)//算出第一行 
        if(a[0][k])//优化 
            for(int j=0;j<n;j++)
            {
                temp[0][j]+=a[0][k]*b[k][j];
                temp[0][j]%=m;
            }
    for(int i=1;i<n;i++)//递推算出剩下n-1行 
        for(int j=0;j<n;j++)
            temp[i][j]=temp[i-1][(j-1+n)%n];
    for(int i=0;i<n;i++)//复制结果 
        for(int j=0;j<n;j++)
            b[i][j]=temp[i][j];
}
void mod_pow(ll a[][maxn])//矩阵快速幂 
{
    //初始化一个单位矩阵 
    memset(ans,0,sizeof(ans));
    for(int i=0;i<n;i++) ans[i][i]=1;
    while(k)
    {
        if(k&1) mod_mul(a,ans);
        mod_mul(a,a);
        k>>=1;
    }
}
int main()
{
    scanf("%d%d%d%d",&n,&m,&d,&k);
    for(int i=0;i<n;i++)
    {
        scanf("%lld",&num[i]);
        f[0][i]=0;
    }
    f[0][0]=1;
    for(int i=1;i<=d;i++)//构造转移矩阵的第一行 
        f[0][i]=f[0][n-i]=1;
    for(int i=1;i<n;i++)//递推算出剩下n-1行 
        for(int j=0;j<n;j++)
            f[i][j]=f[i-1][(j-1+n)%n];
    mod_pow(f);
    ll end[maxn];//储存结果的数组 
    memset(end,0,sizeof(end));//初始化 
    for(int i=0;i<n;i++)//初始序列组成的1*n矩阵与转移矩阵相乘即为答案 
        for(int j=0;j<n;j++)
        {
            end[i]+=ans[i][j]*num[j];
            end[i]%=m;
        }
    for(int i=0;i<n;i++)
        printf("%lld%c",end[i],i==n-1?'\n':' ');    
    return 0;
}