UVALive3704(Cellular Automaton)

• 题目链接

• 题意
  一个细胞自动机含$n$个格子，每个格子的取值为$0{\sim}m-1$。给定距离$d$，则每次操作后每个格子的值将变为到它的距离不超过$d$的所有格子在操作之前的值之和除以$m$的余数，其中i和j的距离为$min\{|i-j|,n-|i-j|\}$。给出$n，m，d，k$和自动机各格子的初始值，计算k次操作以后各格子的值。

• 分析
  把每 $t$ 次操作以后的格子值写成列向量$v_t$，则$v_{t+1}$的每个数都是$v_t$的线性组合，因此可以用矩阵快速幂来求解。例如设$v_{t+1}=Av_t$ ，当$d=1$时有：

$$A=\left( \begin{matrix} 1& 1 &0& \cdots &0&1\\ 1 & 1&1& 0&\cdots & 0\\ 0&1 & 1&1&\cdots & 0\\ \vdots & \cdots & \vdots & \vdots&\ddots&\vdots\\ 1&\cdots&0&0&1&1 \end{matrix} \right)$$
然后，由于 $n\le500$, $O(n^3logk)$的复杂度会超时。而矩阵 $A$是循环矩阵，循环矩阵相乘的结构还是循环矩阵，因此在计算矩阵乘法时只需计算结果矩阵的第一行，其余由第一行推出并保存。矩阵相乘的复杂度降为$O(n^2)$,总时间复杂度降为 $O(n^2logk)$。
- 代码

/*
    矩阵快速幂
    循环矩阵相乘的结果还是循环矩阵
        循环矩阵优化O(N^3) -> O(N^2)
*/
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=510;
ll n,m,d,k,v[maxn];
typedef vector<ll> vec;
typedef vector<vec>mat;

mat mul(mat &A,mat &B,ll m)
{
    mat C(A.size(),vec(B[0].size()));
    for(int k=0;k<B[0].size();++k)
        for(int j=0;j<B[0].size();++j)
            C[0][k]=(C[0][k]+A[0][j]*B[j][k])%m;
    for(int i=1;i<n;++i)
        for(int j=0;j<n;++j)
        {
            C[i][j]=C[i-1][(j-1+n)%n];
        }
    return C;
}

mat qpow(mat A,ll n,ll m)
{
    mat B(A.size(),vec(A.size()));
    for(int i=0;i<A.size();++i)
        B[i][i]=1;
    while(n>0){
        if(n&1) B=mul(B,A,m);
        A=mul(A,A,m);
        n>>=1;
    }
    return B;
}

int main()
{
    while(~scanf("%lld%lld%lld%lld",&n,&m,&d,&k))
    {
        for(int i=0;i<n;++i) scanf("%lld",&v[i]);
        ll ans[maxn]={0};
        mat A(n,vec(n));
        for(int i=0;i<n;++i) A[i][i]=1;
        for(int i=0;i<n;++i)
        {
            for(int j=1;j<=d;++j)
                A[i][(i+j+n)%n]=A[i][(i-j+n)%n]=1;
        }
        A=qpow(A,k,m);
        for(int i=0;i<n;++i)
        {
            ans[i]=0;
            for(int j=0;j<n;++j)
            {
                ans[i]=(ans[i]+v[j]*A[i][j]%m)%m;
            }
            printf("%lld",ans[i]%m);
            if(i==n-1) printf("\n");
            else printf(" ");
        }
    }
    return 0;
}