[模板] 常系数线性递推

Posted 2022-01-30 ubospica

常系数线性递推

给定向量 \(A_0 = (a_1, a_2, \dotsc, a_k)\), 和向量 \(H = (h_1, h_2, \dotsc, h_k)\), 同时

\[ a_n = \sum_i=1^k a_n-i h_i \]

求 \(a_n\).

算法

我们只需求出 \(A_n = (a_n, a_n+1, \dotsc, a_n+k-1)\) 即可.

设 \(f(\lambda)\) 表示转移方程的特征多项式, 有

\[ f(\lambda) = \lambda^k - \sum_i=0^k-1 h_k-i \lambda^i \]

设 \(g(\lambda) \equiv \lambda^n-1 \pmodf(\lambda)\), 那么有

\[ a_n = \sum_i=0^k-1 g_i a_i+1 \]

求 \(g(\lambda)\) 如果直接取模则复杂度为 \(O(k^2\log n)\), 多项式取模则为 \(O(k\log k \log n)\).

代码

bzoj4161: Shlw loves matrixI

#include<cstdio>
#include<iostream>
#include<cmath>
#include<cstring>
#include<algorithm>
#include<queue>
#include<map>
#include<set>
using namespace std;
#define rep(i,l,r) for(register int i=(l);i<=(r);++i)
#define repdo(i,l,r) for(register int i=(l);i>=(r);--i)
#define il inline
typedef long long ll;
typedef double db;

//---------------------------------------
const int nsz=4050;
const ll nmod=1e9+7;
int n,k;

ll v1[nsz],v2[nsz];//a_n = \sum_i=1^k v1[i]v2[n-i]
ll f[nsz],g[nsz],h[nsz],c[nsz];

ll qp(ll a,ll b)
    ll res=0;
    for(;b;b>>=1,a=a*a%nmod)if(b%1)res=res*a%nmod;
    return res;


void mul(ll *a,ll *b,ll *res)//res=a*b%f
    rep(i,0,k*2)c[i]=0;
    rep(i,0,k-1)rep(j,0,k-1)c[i+j]=(c[i+j]+a[i]*b[j])%nmod;
    repdo(i,k*2-2,k)
        if(c[i])
            rep(j,0,k)
                c[i-k+j]=(c[i-k+j]-c[i]*f[j])%nmod;
            
        
    
    rep(i,0,k-1)res[i]=c[i];


void qp(ll *a,ll b,ll *c)
    c[0]=1;
    for(;b;b>>=1,mul(a,a,a))
    if(b&1)
    mul(c,a,c); 


ll getrecurrence()
    if(n<=k)return v2[n]%nmod;
    if(k==1)return qp(v1[1],n)*v2[0]%nmod;
    ++n; 
    rep(i,0,k-1)f[i]=-v1[k-i];
    f[k]=1;
    g[1]=1;
    qp(g,n-1,h);
    ll res=0; 
    rep(i,0,k-1)
        res=(res+h[i]*v2[i])%nmod;
    
    return res;
 

int main()
    ios::sync_with_stdio(0),cin.tie(0);
    cin>>n>>k;
    rep(i,1,k)cin>>v1[i];
    rep(i,0,k-1)cin>>v2[i];
    cout<<(getrecurrence()+nmod)%nmod<<'\n';
    return 0;