vjudge

description

你有一个\(n \times n\)的矩阵\(A\)，下标从\(0\)开始，其中\(A_{i,j}=di + q^i\)。
给你\(d,q,n,k,s,t\)，求\((A^k)_{s,t}\)的值模\(10^9+7\)。
\(n\le10^4,d,q,k\le10^9,0\le s,t <n\)

sol

先写出一个极其丑陋的柿子：
\[(A^k)_{s,t}=\sum_{x_1}\sum_{x_2}...\sum_{x_{k-1}}(ds+q^{x_1})(dx_1+q^{x_2})...(dx_{k-1}+q^t)\]
等于你从\(s\)出发走\(k\)次走到\(t\)，那么中间就有\(k-1\)个位置可以任选，所以才会有前面的\(k-1\)个\(\sum\)。
把后面的括号暴力展开，一共有\(2^k\)项。然而还是做不了。
考虑每个\(x_i\)。它出现在答案中的形式只有\(4\)种：\(1,dx_i,q_{x_i},dx_iq^{x_i}\)。这取决于这个\(x_i\)出现的两个相邻的括号内取的是哪一个。
预处理出\(sd=\sum_{i=0}^{n-1}di,sq=\sum_{i=0}^{n-1}q^i,ss=\sum_{i=0}^{n-1}diq^i\)。
设\(f_{i,0/1}\)表示考虑了前\(i\)个括号，第\(i\)个括号取的是前一项/后一项时，只计算前\(i-1\)个\(x_i\)时的答案。显然有转移式：
\[f_{i+1,0}=f_{i,0}\times sd+f_{i,1}\times ss\\f_{i+1,1}=f_{i,0}\times n+f_{i,1}\times sq\]
初值\(f_{1,0}=ds,f_{1,1}=1\)。最后一次转移特判一下（因为最后一个\(x\)必须强制为\(t\)）。
这样就可以做\(O(k)\)的\(dp\)了。考虑到每次的转移都是一样的，所以直接矩乘就好了。

code

为什么没有人告诉我，定义在\(\mbox{class}\)里面的变量的初值不是\(0\)。
害我浪费了宝贵的\(\mbox{20min}\)才发现这个错误\(\mbox{qwq}\)。

#include<cstdio>
#include<algorithm>
#include<cstring>
#include<vector>
using namespace std;
int gi(){
    int x=0,w=1;char ch=getchar();
    while ((ch<'0'||ch>'9')&&ch!='-') ch=getchar();
    if (ch=='-') w=0,ch=getchar();
    while (ch>='0'&&ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar();
    return w?x:-x;
}
const int N = 10005;
const int mod = 1e9+7;
int fastpow(int a,int b){
    int s=1;
    while (b) {if (b&1) s=1ll*s*a%mod;a=1ll*a*a%mod;b>>=1;}
    return s;
}
struct matrix{
    int a[2][2];
    matrix(){memset(a,0,sizeof(a));}
    matrix(int q,int w,int d,int b){
        a[0][0]=q;a[0][1]=w;a[1][0]=d;a[1][1]=b;
    }
    int* operator [](int x){return a[x];}
    matrix operator * (matrix &b){
        matrix c;
        for (int i=0;i<2;++i)
            for (int j=0;j<2;++j)
                for (int k=0;k<2;++k)
                    c[i][k]=(1ll*a[i][j]*b[j][k]+c[i][k])%mod;
        return c;
    }
}S,T,A,B;
class MatrixPower{
public:
    int m,sd,sq,ss;vector<int>ans;
    vector<int> getElements(int d,int q,int n,int k,vector<int>r,vector<int>c){
        m=r.size();ans.resize(m);
        if (k==1){
            for (int i=0;i<m;++i) ans[i]=(1ll*d*r[i]+fastpow(q,c[i]))%mod;
            return ans;
        }
        k-=2;sd=sq=ss=0;
        for (int i=0;i<n;++i) sd=(sd+1ll*d*i)%mod;
        for (int i=0;i<n;++i) sq=(sq+fastpow(q,i))%mod;
        for (int i=0;i<n;++i) ss=(ss+1ll*d*i%mod*fastpow(q,i))%mod;
        T=matrix(sd,n,ss,sq);S=matrix(1,0,0,1);
        while (k) {if (k&1) S=S*T;T=T*T;k>>=1;}
        for (int i=0;i<m;++i){
            A=matrix(1ll*d*r[i]%mod,1,0,0);
            B=matrix(sd,1ll*n*fastpow(q,c[i])%mod,ss,1ll*sq*fastpow(q,c[i])%mod);
            A=A*S*B;
            ans[i]=(1ll*A[0][0]+A[0][1]+A[1][0]+A[1][1])%mod;
        }
        return ans;
    }
};