题目描述

给定数列 {hn}前k项，其后每一项满足

hn = a1h(n-1) + a2h(n-2) + ... + ak*h(n-k)

其中 a1,a2...ak 为给定数列。请计算 h(n)，并将结果对 1000000007 取模输出。

解法

一个显然的思路就是矩阵乘法，但这样的话显然超时。

实际上，我们还可以继续对这个矩阵乘法进行优化。

首先，由于这是常系数齐次线性递推式，简单来说就是：

\[f_i=\sum_{j=1}^k a_i*f_{i-j}
\]

然后，我们需要引进特征多项式这个概念。

对于一个矩阵\(A\)，它的特征多项式是\(f(x)=|Ix-A|\)

把行列式展开之后得，\(f(x)=|Ix-A|=x^K-\sum_{i=1}^K a_i*x^{K-i}\)

由Cayley-hamilton定理，那么我们知道\(f(A)=0\)

然后就能知道一个关键的式子：

\[A^K=\sum_{i=1}^K a_i*A^{K-i}
\]

然后由于\(A^i\)都能表示成\(A^1,A^2,....,A^K\)的线性组合，

所以现在矩阵乘法直接使用\(O(K^2)\)的卷积即可。

当得到了\(A^n=\sum_{i=1}^K c_i*A^i\)之后，我们乘上题目给的向量h。

那么就会有\(\sum_{i=1}^K c_i*A^i*h_K\)，即\(Ans=\sum_{i=1}^K c_i*h_{K+i}\)。

复杂度就被优化为\(O(K^2 log n)\)

Code

#include<bits/stdc++.h>

#define ll long long

#define fo(i,x,y) for(int i=x;i<=y;i++)

#define fd(i,x,y) for(int i=x;i>=y;i--)

using namespace std;

const int maxn=4007,mo=1e9+7;

int n,K;

int a[maxn],h[maxn],f[maxn];

int Ans;

void Init(){

	scanf("%d%d",&n,&K);

	fo(i,1,K) scanf("%d",&a[i]);

	memcpy(f,a,sizeof a);

	fo(i,1,K) scanf("%d",&h[i]);

}

#define PLUS(x,y) (x)=((x)+(y))%mo

void mult(int *a,int *b){

	static int c[maxn];

	memset(c,0,sizeof c);

	fo(i,1,K)

		fo(j,1,K)

			PLUS(c[i+j],1ll*a[i]*b[j]);

	fd(i,2*K,K+1)

		fo(j,1,K)

			PLUS(c[i-j],1ll*c[i]*f[j]);

	memcpy(a,c,sizeof c);

}

void qPower(int x){

	bool bz=false;

	static int b[maxn];

	b[1]=1;

	while (x){

		if (x&1){

			if (bz) mult(a,b);

			else{

				bz=true;

				memcpy(a,b,sizeof b);

			}

		}

		mult(b,b);

		x>>=1;

	}

}

void Solve(){

	n++;

	fo(i,K+1,2*K) fo(j,1,K) PLUS(h[i],1ll*h[i-j]*a[j]);

	if (n<=2*K){

		Ans=h[n];

		return;

	}

	qPower(n-K);

	Ans=0;

	fo(i,1,K) PLUS(Ans,1ll*a[i]*h[i+K]);

}

void Print(){

	Ans=(Ans+mo)%mo;

	printf("%d\n",Ans);

}

int main(){

	freopen("1.in","r",stdin);

	freopen("1.out","w",stdout);

	Init();

	Solve();

	Print();

	return 0;

}

ix

【BZOJ4161】Shlw loves matrixI

题目描述

解法

Code