bzoj3231[Sdoi2008]递归数列

Description

一个由自然数组成的数列按下式定义：
对于i <= k：ai = bi
对于i > k: ai = c1ai-1 + c2ai-2 + … + ckai-k
其中bj和 cj （1<=j<=k）是给定的自然数。写一个程序，给定自然数m <= n, 计算am + am+1 + am+2 + … + an, 并输出它除以给定自然数p的余数的值。

Input

由四行组成。
第一行是一个自然数k。
第二行包含k个自然数b1, b2,…,bk。
第三行包含k个自然数c1, c2,…,ck。
第四行包含三个自然数m, n, p。

Output

仅包含一行：一个正整数，表示(am + am+1 + am+2 + … + an) mod p的值。

HINT

对于100%的测试数据：
1<= k<=15
1 <= m <= n <= 10^18
对于所有测试数据：
0<= b1, b2,… bk, c1, c2,…, ck<=10^9
1 <= p <= 10^8

矩阵快速幂裸题.
我们只需要在构造矩阵时再维护一个前缀和就好，这样两遍快速幂，再用结果做个差即可
构造矩阵为

$$\begin{bmatrix} C1& C2&C3&...&CK \\ 1 & 0&0&...&0\\ 0&1 &0&...&0\\ ...&...&...&...&...\\ 1&0&0&...&0 \end{bmatrix}\quad \times \begin{bmatrix} Ak\\ Ak-1\\ Ak-2\\ ...\\ A1\\ sum_{k-1}\\ \end{bmatrix}\quad$$
注意特判掉 $n,m<K$的情况.
代码

    #include <cstdio>
#include <iostream>
#include <cstring>
#define maxn 18
#define mem(a,b) memset(a,b,sizeof(a))
using namespace std;
inline void read(int& x)
{   char c=getchar();x=0;int y=1;
    while(c<'0'||c>'9'){if(c=='-') y=-1;c=getchar();}
    while(c>='0'&&c<='9') x=x*10+c-'0',c=getchar();
    x*=y;
}
typedef long long ll;
ll sum[maxn],n,m;
int C[maxn],K,mod;
struct matrix
{   int s[maxn][maxn];
    matrix(){mem(s,0);}
    inline void init(){mem(s,0);for(int i=0;i<K;++i) this->s[i][i]=1;}
    inline int* operator [](int x){return s[x];}
    inline friend matrix operator *(matrix x,matrix y)
    {   matrix z;
        for(int i=0;i<K;++i)
            for(int j=0;j<K;++j)
                for(int k=0;k<K;++k)
                    z[i][j]=(1ll*z[i][j]+1ll*x[i][k]*y[k][j]%mod)%mod;
        return z;
    }
}A,ans,F,__backup,__rebackup;
int main()
{   read(K);
    for(int i=1;i<=K;++i) read(F[K-i][0]);
    for(int i=1;i<=K;++i) read(C[i]);
    scanf("%lld%lld%d",&m,&n,&mod);
    for(int i=1;i<=K;++i) F[K-i][0]%=mod,sum[i]=(sum[i-1]+F[K-i][0])%mod;
    sum[K]=(sum[K]-F[0][0]+mod)%mod;++K;
    F[K-1][0]=sum[K-1];
    for(int i=0;i<K;++i) if(i!=K-1) A[0][i]=C[i+1];
    for(int i=1,k=K-1;i<k;++i) A[i][i-1]=1;
    A[K-1][0]=1,A[K-1][K-1]=1;__backup=A;
    ll y=m-K+1;ans.init();
    if(y>=0){for(;y;y>>=1,A=A*A) if(y&1) ans=ans*A;__rebackup=ans*F;}
    else __rebackup[K-1][0]=sum[m-1];
    y=n-K+1;A=__backup;ans.init();
    if(y>=0){for(;y;y>>=1,A=A*A) if(y&1) ans=ans*A;F=ans*F;}
    else F[K-1][0]=sum[n];
    printf("%d",((ll)F[K-1][0]+F[0][0]-__rebackup[K-1][0]+mod)%mod);
    return 0;
}