bzoj#4161-Shlw loves matrixI【常係數線性齊次遞推】
正題
題目連結:https://darkbzoj.tk/problem/4161
題目大意
給出序列\(a\)，和\(h\)的\(0\sim k-1\)項，滿足
\[h_n=\sum_{i=1}^na_ih_{n-i} \]求\(h_n\)。
\(1\leq n\leq 10^9,1\leq k\leq 2000\)
解題思路
顯然\(k\leq 2000\)我們就不能矩陣乘法了，然後就去研究了一下常係數線性齊次遞推，而且這題\(k^2\)能過所以不用\(O(n\log n)\)的多項式取模。
首先對於一個\(k\times k\)的矩陣\(A\)的特徵多項式\(f(x)\)的定義是
\[f(x)=\sum_{i=0}^kdet(A-k\lambda)x^i \]而我們知道對於這種遞推式的特徵多項式就是
然後有個Cayley-Hamilton定理就是說對於\(A\)的特徵多項式\(f(x)\)滿足\(f(A)=0\)。
如果按照這種特殊的特徵多項式去理解的話挺好懂的，但是實際上的就不會證了。
那麼我們如果要求一個轉移矩陣\(A^n=A^n\mod f(A)\)，我們可以求出一個多項式\(r(x)=x^n\mod f(x)\)，然後直接把\(A\)帶入這個多項式就好了。
這個東西要用到快速冪+多項式取模，\(O(k^2)\)的話實現起來很簡單，考慮一個\(x^{k+i}(i\geq 0)\mod f(x)\)會發生什麼，因為\(f(x)[k]=1\)
然後我們得到了一個多項式\(r(x)\)，考慮怎麼計算答案，只需要計算\(A^{n-k}\vec{h}\)的最後一項，所以我們實際上求的\(r(x)=x^{n-k}\mod f(x)\)，記\(r(x)\)的\(i\)次係數為\(b_i\)，那麼我們有
\[r(A)\vec{h}=\sum_{i=0}^{k-1}b_iA^{i}\vec{h} \]所以我們算出\(h\)的\(k\sim 2k-1\)項然後就有
\[h_n=\sum_{i=0}^{k-1}r(x)[i]h_{k+i} \]時間複雜度：\(O(k^2\log n)\)
值得一提的是如果轉移矩陣\(A\)不是一個遞推式的矩陣，我們也可以用拉格朗日插值插出它的特徵多項式從而降低平均的複雜度。
code
#include<cstdio>
#include<cstring>
#include<algorithm>
#define ll long long
using namespace std;
const ll N=4100,P=1e9+7;
ll n,k,a[N],p[N],h[N],b[N],f[N],tmp[N];
void Mul(ll *a,ll *b,ll *c){
memset(tmp,0,sizeof(tmp));
for(ll i=0;i<k;i++)
for(ll j=0;j<k;j++)
(tmp[i+j]+=a[i]*b[j]%P)%=P;
for(ll i=2*k-2;i>=k;i--)
for(ll j=0;j<k;j++)
(tmp[i+j-k]+=P-tmp[i]*p[j]%P)%=P;
for(ll i=0;i<k;i++)c[i]=tmp[i];
return;
}
signed main()
{
scanf("%lld%lld",&n,&k);
for(ll i=1;i<=k;i++)
scanf("%lld",&a[i]),a[i]=(a[i]+P)%P,p[k-i]=P-a[i];
for(ll i=0;i<k;i++)
scanf("%lld",&h[i]),h[i]=(h[i]+P)%P;
for(ll i=k;i<(k<<1);i++)
for(ll j=1;j<=k;j++)
(h[i]+=h[i-j]*a[j]%P)%=P;
if(n<(k<<1))return printf("%lld\n",h[n])&0;
b[1]=p[k]=f[0]=1;
ll B=n-k,ans=0;
while(B){
if(B&1)Mul(f,b,f);
Mul(b,b,b);B>>=1;
}
for(ll i=0;i<k;i++)
(ans+=f[i]*h[i+k]%P)%=P;
printf("%lld\n",ans);
return 0;
}