概述
题目:bzoj_4161
题目描述:
给定数列 h n {h_n} hn前k项(下标从0开始),其后每一项满足 h n = a 1 ∗ h n − 1 + a 2 ∗ h n − 2 + . . . + a k ∗ h n − k h_n = a_1*h_{n-1} + a_2*h_{n-2} + ... + a_k*h_{n-k} hn=a1∗hn−1+a2∗hn−2+...+ak∗hn−k其中 a1,a2…ak 为给定数列。请计算 h n h_n hn,并将结果对 1000000007 取模输出。
题解:
不妨将h的下标从1开始计,则所求值为
h
n
+
1
h_{n+1}
hn+1首先,最容易想出的办法是暴力,在暴力的基础上,我们可以引入矩阵快速幂,容易想到转移矩阵如下:
A
=
(
a
1
a
2
a
3
⋯
a
k
−
1
a
k
1
0
0
⋯
0
0
0
1
0
⋯
0
0
0
0
1
⋯
0
0
⋮
⋮
⋮
⋱
⋮
⋮
0
0
0
⋯
1
0
)
A= left( begin{matrix} a_1&a_2&a_3&cdots&a_{k-1}&a_k\ 1&0&0&cdots&0&0\ 0&1&0&cdots&0&0\ 0&0&1&cdots&0&0\ vdots&vdots&vdots&ddots&vdots&vdots\ 0&0&0&cdots&1&0 end{matrix} right)
A=⎝⎜⎜⎜⎜⎜⎜⎜⎛a1100⋮0a2010⋮0a3001⋮0⋯⋯⋯⋯⋱⋯ak−1000⋮1ak000⋮0⎠⎟⎟⎟⎟⎟⎟⎟⎞
而初始矩阵如下:
H
1
=
(
h
k
h
k
−
1
⋮
h
2
h
1
)
H_1= left( begin{matrix} h_k\h_{k-1}\ vdots\h_2\h_1 end{matrix} right)
H1=⎝⎜⎜⎜⎜⎜⎛hkhk−1⋮h2h1⎠⎟⎟⎟⎟⎟⎞
则转移过后的矩阵为
H
n
+
1
=
A
n
∗
H
1
=
(
h
k
+
n
h
k
+
n
−
1
⋮
h
n
+
2
h
n
+
1
)
H_{n+1}=A^n*H_1= left( begin{matrix} h_{k+n}\h_{k+n-1}\ vdots\h_{n+2}\h_{n+1} end{matrix} right)
Hn+1=An∗H1=⎝⎜⎜⎜⎜⎜⎛hk+nhk+n−1⋮hn+2hn+1⎠⎟⎟⎟⎟⎟⎞
然而这样是会超时的。计算的瓶颈在于
A
n
A^n
An的求法,因此需要用到矩阵的特征多项式取模来优化矩阵快速幂。将矩阵的乘法转成多项式的乘法。
前置知识:
1.特征值:
设 A A A是n阶矩阵,如果数 λ lambda λ和n维向量 x x x使关系式 (1) A x = λ x Ax = lambda x tag{1} Ax=λx(1)成立,那么,这样的数 λ lambda λ称为矩阵 A A A的 特征值
2.特征多项式:
(1)式也可写成
(
A
−
λ
E
)
x
=
0
,
(A-lambda E)x=0,
(A−λE)x=0,这是n个未知数n个方程的齐次线性方程组,它有非零解的充分必要条件是系数行列式
∣
A
−
λ
E
∣
=
0
,
|A-lambda E|=0,
∣A−λE∣=0,即
(
a
11
−
λ
a
12
⋯
a
1
n
a
21
a
22
−
λ
⋯
a
2
n
⋮
⋮
⋮
a
n
1
a
n
2
⋯
a
n
n
−
λ
)
=
0
left( begin{matrix} a_{11}-lambda&a_{12}&cdots&a_{1n}\ a_{21}&a_{22}-lambda&cdots&a_{2n}\ vdots&vdots&&vdots\ a_{n1}&a_{n2}&cdots&a_{nn}-lambda end{matrix} right)=0
⎝⎜⎜⎜⎛a11−λa21⋮an1a12a22−λ⋮an2⋯⋯⋯a1na2n⋮ann−λ⎠⎟⎟⎟⎞=0
上式是以
λ
lambda
λ为未知数的一元n次方程,成为矩阵
A
A
A的***特征方程***,其左端
∣
A
−
λ
E
∣
|A-lambda E|
∣A−λE∣是
λ
lambda
λ的n次多项式,记作
f
(
λ
)
f(lambda)
f(λ),称为矩阵
A
A
A的***特征多项式***。
3.Cayley-Hamilton定理:
对 A A A的特征多项式 f ( λ ) f(lambda) f(λ),将矩阵 A A A作为自变量带入,可得 f ( A ) = 0 f(A)=0 f(A)=0详细证明可以参考Wikipedia,这里是传送门。有一个很简单的理解方式就是直接代入,得到 f ( A ) = ∣ A − A E ∣ = 0 f(A)=|A-AE|=0 f(A)=∣A−AE∣=0,看起来还是挺显然的。
4.多项式取模:
对于一个多项式
A
(
x
)
A(x)
A(x),称其最高项的次数为这个多项式的度(degree),记作
d
e
g
(
A
)
deg(A)
deg(A)。设多项式
g
(
x
)
g(x)
g(x)为被除式,
d
e
g
(
g
)
=
n
deg(g)=n
deg(g)=n,多项式
f
(
x
)
f(x)
f(x)为除式,
d
e
g
(
f
)
=
m
(
n
>
m
)
deg(f)=m(n>m)
deg(f)=m(n>m)。只要多项式
f
(
x
)
f(x)
f(x)不为常数0,则存在
g
(
x
)
=
f
(
x
)
∗
p
(
x
)
+
r
(
x
)
g(x)=f(x)*p(x)+r(x)
g(x)=f(x)∗p(x)+r(x),其中
d
e
g
(
r
)
<
m
deg(r)<m
deg(r)<m。其中
p
(
x
)
p(x)
p(x)为商式,
r
(
x
)
r(x)
r(x)为余式。已知
f
(
x
)
f(x)
f(x)和
g
(
x
)
g(x)
g(x)来求
r
(
x
)
r(x)
r(x)的过程就叫做多项式取模。暴力取模可以直接模拟竖式除法。复杂度是
o
(
n
2
)
o(n^2)
o(n2)(n与m同阶)。也可以用FFT来加速,达到log的复杂度,但是常数巨大,此题不适用。
有了这四样东西后,我们就可以用来做题了。先看原来的转移矩阵的特征多项式,即
f
(
λ
)
=
∣
A
−
λ
E
∣
=
∣
a
1
−
λ
a
2
a
3
⋯
a
k
−
1
a
k
1
−
λ
0
⋯
0
0
0
1
−
λ
⋯
0
0
0
0
1
⋯
0
0
⋮
⋮
⋮
⋱
⋮
⋮
0
0
0
⋯
1
−
λ
∣
f(lambda)=|A-lambda E|= left| begin{matrix} a_1-lambda&a_2&a_3&cdots&a_{k-1}&a_k\ 1&-lambda&0&cdots&0&0\ 0&1&-lambda&cdots&0&0\ 0&0&1&cdots&0&0\ vdots&vdots&vdots&ddots&vdots&vdots\ 0&0&0&cdots&1&-lambda end{matrix} right|
f(λ)=∣A−λE∣=∣∣∣∣∣∣∣∣∣∣∣∣∣a1−λ100⋮0a2−λ10⋮0a30−λ1⋮0⋯⋯⋯⋯⋱⋯ak−1000⋮1ak000⋮−λ∣∣∣∣∣∣∣∣∣∣∣∣∣
为了展开更方便,我们给等式乘上
(
−
1
)
k
(-1)^k
(−1)k即
(
−
1
)
k
f
(
λ
)
=
(
−
1
)
k
∣
A
−
λ
E
∣
=
(
−
1
)
k
∣
a
1
−
λ
a
2
a
3
⋯
a
k
−
1
a
k
1
−
λ
0
⋯
0
0
0
1
−
λ
⋯
0
0
0
0
1
⋯
0
0
⋮
⋮
⋮
⋱
⋮
⋮
0
0
0
⋯
1
−
λ
∣
=
∣
λ
−
a
1
−
a
2
−
a
3
⋯
−
a
k
−
1
−
a
k
−
1
λ
0
⋯
0
0
0
−
1
λ
⋯
0
0
0
0
−
1
⋯
0
0
⋮
⋮
⋮
⋱
⋮
⋮
0
0
0
⋯
−
1
λ
∣
begin{aligned} &(-1)^kf(lambda)\ =&(-1)^k|A-lambda E|\ =&(-1)^k left| begin{matrix} a_1-lambda&a_2&a_3&cdots&a_{k-1}&a_k\ 1&-lambda&0&cdots&0&0\ 0&1&-lambda&cdots&0&0\ 0&0&1&cdots&0&0\ vdots&vdots&vdots&ddots&vdots&vdots\ 0&0&0&cdots&1&-lambda end{matrix} right| \=& left| begin{matrix} lambda-a_1&-a_2&-a_3&cdots&-a_{k-1}&-a_k\ -1&lambda&0&cdots&0&0\ 0&-1&lambda&cdots&0&0\ 0&0&-1&cdots&0&0\ vdots&vdots&vdots&ddots&vdots&vdots\ 0&0&0&cdots&-1&lambda end{matrix} right| end{aligned}
===(−1)kf(λ)(−1)k∣A−λE∣(−1)k∣∣∣∣∣∣∣∣∣∣∣∣∣a1−λ100⋮0a2−λ10⋮0a30−λ1⋮0⋯⋯⋯⋯⋱⋯ak−1000⋮1ak000⋮−λ∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣∣λ−a1−100⋮0−a2λ−10⋮0−a30λ−1⋮0⋯⋯⋯⋯⋱⋯−ak−1000⋮−1−ak000⋮λ∣∣∣∣∣∣∣∣∣∣∣∣∣
将行列式按照第一行展开,可以得到
(
−
1
)
k
f
(
λ
)
=
λ
k
−
∑
i
=
1
k
a
i
∗
x
k
−
i
(-1)^kf(lambda)=lambda^k-sum_{i=1}^k{a_i*x^{k-i}}
(−1)kf(λ)=λk−i=1∑kai∗xk−i
不妨令该式为
g
(
λ
)
g(lambda)
g(λ),代入
A
A
A,进行多项式取模,有如下等式
A
n
=
g
(
A
)
∗
p
(
A
)
+
r
(
A
)
d
e
g
(
g
)
=
k
A^n=g(A)*p(A)+r(A)\ deg(g)=k
An=g(A)∗p(A)+r(A)deg(g)=k则有
d
e
g
(
r
)
<
k
deg(r)<k
deg(r)<k,不妨看作
d
e
g
(
r
)
=
k
−
1
deg(r)=k-1
deg(r)=k−1,取模后得到r的系数
c
0
,
c
1
,
c
2
⋯
c
k
−
1
c_0,c_1,c_2cdots c_{k-1}
c0,c1,c2⋯ck−1。由于
g
(
A
)
=
0
g(A)=0
g(A)=0,得到
A
n
=
r
(
A
)
=
c
0
+
c
1
∗
A
1
+
⋯
+
c
k
−
1
∗
A
k
−
1
A^n=r(A)=c_0+c_1*A^1+cdots+c_{k-1}*A^{k-1}
An=r(A)=c0+c1∗A1+⋯+ck−1∗Ak−1
则
H
n
+
1
=
A
n
∗
H
1
=
∑
i
=
0
k
−
1
c
i
∗
A
i
∗
H
1
=
∑
i
=
0
k
−
1
c
i
∗
H
i
+
1
H_{n+1}=A^n*H_1=sum_{i=0}^{k-1}{c_i*A_i*H_1}=sum_{i=0}^{k-1}{c_i*H_{i+1}}
Hn+1=An∗H1=i=0∑k−1ci∗Ai∗H1=i=0∑k−1ci∗Hi+1
即
H
n
+
1
=
c
0
∗
(
h
k
⋮
h
2
h
1
)
+
c
1
∗
(
h
k
+
1
⋮
h
3
h
2
)
+
⋯
+
c
k
−
1
∗
(
h
2
k
⋮
h
k
+
1
h
k
)
H_{n+1}=c_0* left( begin{matrix} h_k\ vdots\h_2\h_1 end{matrix} right)+c_1* left( begin{matrix} h_{k+1}\ vdots\h_3\h_2 end{matrix} right)+ cdots+ c_{k-1}* left( begin{matrix} h_{2k}\ vdots\h_{k+1}\h_k end{matrix} right)
Hn+1=c0∗⎝⎜⎜⎜⎛hk⋮h2h1⎠⎟⎟⎟⎞+c1∗⎝⎜⎜⎜⎛hk+1⋮h3h2⎠⎟⎟⎟⎞+⋯+ck−1∗⎝⎜⎜⎜⎛h2k⋮hk+1hk⎠⎟⎟⎟⎞
那么
h
n
+
1
=
∑
i
=
0
k
−
1
c
i
h
i
+
1
h_{n+1}=sum_{i=0}^{k-1}{c_ih_{i+1}}
hn+1=∑i=0k−1cihi+1,即为答案。
求
c
n
c_n
cn数列的过程就是用
x
n
x^n
xn来模上
g
(
x
)
g(x)
g(x)所得的系数,这个过程可以用快速幂来实现。代码如下
#include<cstdio>
#include<algorithm>
using namespace std;
const int N=2005,mod=1000000007;
int n,k,a[N],h[N],c[N],b[N],tmp[N<<1];
void mul(int *a1,int *a2){
for(int i=0;i<=2*k;i++)tmp[i]=0;
for(int i=0;i<=k;i++){
for(int j=0;j<=k;j++){
tmp[i+j]=(tmp[i+j]+1ll*a1[i]*a2[j])%mod;
}
}
for(int i=2*k;i>=k;i--){
for(int j=0;j<k;j++){
tmp[i-k+j]=(tmp[i-k+j]+1LL*tmp[i]*a[k-j])%mod;
}
tmp[i]=0;
}
for(int i=0;i<k;i++)a1[i]=tmp[i];
}
void pow(int *a1,int t,int *a2){
for(;t;t>>=1){
if(t&1)mul(a2,a1);
mul(a1,a1);
}
}
int main(){
#ifndef ONLINE_JUDGE
freopen("1.in","r",stdin);
#endif
scanf("%d%d",&n,&k);
for(int i=1;i<=k;i++)scanf("%d",&a[i]);
for(int i=1;i<=k;i++){
scanf("%d",&h[i]);
h[i]=(1ll*h[i]%mod+mod)%mod;
}
if(n<k)return printf("%dn",h[n+1])*0;
b[1]=1;c[0]=1;
pow(b,n,c);
int ans=0;
for(int i=0;i<k;i++){
ans=(ans+1LL*c[i]*h[i+1]%mod+mod)%mod;
}
printf("%dn",ans);
return 0;
}
最后
以上就是从容小丸子为你收集整理的常系数线性递推——多项式取模优化矩阵快速幂题目:bzoj_4161的全部内容,希望文章能够帮你解决常系数线性递推——多项式取模优化矩阵快速幂题目:bzoj_4161所遇到的程序开发问题。
如果觉得靠谱客网站的内容还不错,欢迎将靠谱客网站推荐给程序员好友。
发表评论 取消回复