[SDOI2008]递归数列 题解

矩乘板子题

\(k \le 15\)\(a_{i}= {\textstyle \sum_{j=i}^{1}}c_{j}*a_{i-j} \ \ \ \ (i>k)\) (这一看就很矩乘) 考虑矩阵加速。

题目让求一段区间的和,可以转化为两前缀和相减的形式,同时,让矩阵边长加上 \(1\) , 多加的 \(1\) 用来储存前缀和

状态矩阵 \(f_{i}\) 存储:

\[\begin{bmatrix} a_{i}&a_{i+1} &a_{i+2} &··· &a_{i+k} &S_{i-1} \end{bmatrix}\]

我们想让他转移到 \(f_{i+1}\) 即为:

\[\begin{bmatrix} a_{i+1}&a_{i+2} &a_{i+3} &··· &a_{i+k+1} &S_{i} \end{bmatrix}\]

很容易得到,转移矩阵 \(A\) 为:

\[\begin{bmatrix} 0 &0 &0 &··· &c_{k} &1 \\ 1 &0 &0 &··· &c_{k-1} &0 \\ 0 &1 &0 &··· &c_{k-2} &0 \\ ··· &··· &··· &··· &··· &··· \\ 0 &0 &0 &··· &c_{1} &0 \\ 0 &0 &0 &··· &0 &1 \end{bmatrix}\]

那矩乘的柿子就有了:

\[\begin{bmatrix} a_{i}&a_{i+1} &a_{i+2} &··· &a_{i+k} &S_{i-1} \end{bmatrix} \begin{bmatrix} 0 &0 &0 &··· &c_{k} &1 \\ 1 &0 &0 &··· &c_{k-1} &0 \\ 0 &1 &0 &··· &c_{k-2} &0 \\ ··· &··· &··· &··· &··· &··· \\ 0 &0 &0 &··· &c_{1} &0 \\ 0 &0 &0 &··· &0 &1 \end{bmatrix} = \begin{bmatrix} a_{i+1}&a_{i+2} &a_{i+3} &··· &a_{i+k+1} &S_{i} \end{bmatrix} \]

举个栗子,对于 \(k=4\) 时,有:

\[\begin{bmatrix} a_{i}& a_{i+1}& a_{i+2}& a_{i+3}& S_{i-1} \end{bmatrix} \begin{bmatrix} 0& 0& 0& c_{4}& 1\\ 1& 0& 0& c_{3}& 0\\ 0& 1& 0& c_{2}& 0\\ 0& 0& 1& c_{1}& 0\\ 0& 0& 0& 0&1 \end{bmatrix} = \begin{bmatrix} a_{i+1}& a_{i+2}& a_{i+3}& a_{i+4}& S_{i} \end{bmatrix}\]

\(S_{0}\) 计算到 \(S_{n}\) 需要进行矩阵乘法需要使 \(f_{0}\)\(n\)\(A\),

及要计算 \(f_{0}*A^{n}\) 计算 \(A^{n}\) 用矩阵快速幂加速运算即可。

Code

#include<cstdio>
#include<cstring>
#include<algorithm>
#define N 17
#define LL long long 
using namespace std;

LL k,b[N],c[N],n,m,mod;

struct matrix
{
	LL a[N][N];
	void clear()//矩阵清空
	{
		memset(a,0,sizeof(a));
	}
	void build()//建立单位矩阵
	{
		memset(a,0,sizeof(a));
		for(register int i=1;i<=k;i++)
			a[i][i]=1;
	}
	void mat_mod()//矩阵取模
	{
		for(register int i=1;i<=k;i++)
			for(register int j=1;j<=k;j++)
					a[i][j]=a[i][j]%mod;
	}
	void print()//矩阵输出
	{
		for(register int i=1;i<=k;i++)
		{
			for(register int j=1;j<=k;j++)
				printf("%lld ",a[i][j]);
			printf("\n");
		}
	}
};

matrix operator *(const matrix x,const matrix y)//矩阵乘法
{
	matrix ans;
	ans.clear();
	for(register int i=1;i<=k;i++)
		for(register int j=1;j<=k;j++)
			for(register int p=1;p<=k;p++)
				ans.a[i][j]=(ans.a[i][j]+x.a[i][p]*y.a[p][j]%mod)%mod;
	return ans;
}

inline LL qr()
{
	char ch=0;LL w=1,x=0;
	while(ch<'0'||ch>'9'){if(ch=='-')w=-1;ch=getchar();}
	while(ch<='9'&&ch>='0'){x=(x<<3)+(x<<1)+(ch^48);ch=getchar();}
	return x*w;
}

matrix poww(matrix x,LL opl)//矩阵快速幂
{
	matrix ans;
	ans.build();//创建单位矩阵
	while(opl)
	{
		if(opl&1)
			ans=ans*x;
		x=x*x;
		opl>>=1;
	}
	return ans;
}

int main()
{
	k=qr();
	for(register int i=1;i<=k;i++)
		b[i]=qr();
	for(register int i=1;i<=k;i++)
		c[i]=qr();
	n=qr();m=qr();mod=qr();
	k++;

	matrix A;//A为转化矩阵
	A.clear();
	for(register int i=2;i<k;i++)//给A填数
		A.a[i][i-1]=1;
	for(register int i=1;i<k;i++)
		A.a[i][k-1]=c[k-i];
	A.a[1][k]=1;
	A.a[k][k]=1;
	A.mat_mod();

	matrix f1,f2;
	f1.clear();//f1,f2为状态矩阵
	f2.clear();
	for(register int i=1;i<k;i++)//给f1,f2填数
		f1.a[1][i]=f2.a[1][i]=b[i];
	f1.mat_mod();
	f2.mat_mod();

	matrix jc1=poww(A,n-1ll);
	matrix jc2=poww(A,m);//A^m
	f1=f1*jc1;//计算前缀和
	f2=f2*jc2;

	printf("%lld\n",((f2.a[1][k]-f1.a[1][k])%mod+mod)%mod);//前缀和减一下即为最终答案
	return 0;
}
posted @ 2021-02-05 08:18  江北南风  阅读(107)  评论(1编辑  收藏  举报