模板题：传送门
分析
严格规定 (f_k eq 0)
不难列出转移矩阵：(left( egin{matrix} a_{n} \ a_{n-1} \ a_{n-2} \ vdots \ a_{n-k} end{matrix} ight)= left( egin{matrix} f_1 & f_2 & cdots & f_{k-1} & f_k \ 1 \ & 1 \ && ddots \ &&& 1&0 end{matrix} ight)cdot left( egin{matrix} a_{n-1} \ a_{n-2} \ a_{n-3} \ vdots \ a_{n-k-1} end{matrix} ight))
记该转移为 (vec A_n=F_kcdot vec A_{n-1})
不难得出 (vec A_n=F_k^n cdot vec A_0)
这里形式化定义了 $a_{-1}, a_{-2}, a_{-3},cdots $
现在来进行一个求解：
假设我们能找到一个次数尽可能小的多项式函数 (G) 使得 (G(F_k)=O_k) 即 (k) 阶零矩阵
则我们可以得到： (vec A_n=F_k^ncdot vec A_0=[ Q(F_k)G(F_k)+R(F_k) ]cdot vec A_0=R(F_k)cdot vec A_0)
其中，(R) 为 (x^n mod G(x))
例如 (n=2,G(x)=x^2-x-1) 时，(R(x)=x^2mod (x^2-x-1)=x+1)
那么，显然有 (displaystyle vec A_n=(sum_{i=0}^{ ext{deg }G-1} r_icdot F_k^i)cdot vec A_0=sum_{i=0}^{ ext{deg }G-1} r_icdot (F_k^icdot vec A_0)=sum_{i=0}^{ ext{deg }G-1}r_icdot vec A_i)
( ext{deg }G) 表示多项式 (G(x)) 的次数
由于我们最后所求为 (a_n) ，仅为 (vec A_n) 的第一个元素；故对等式右边向量 (vec A_i) ，均只需取第一个元素 (a_i)
如果觉得不严谨，可以认为是方程两边同时右乘一个向量 (vec v=(1, 0, 0, cdots , 0))
故 (displaystyle a_n=sum_{i=0}^{ ext{deg }G-1}r_icdot a_i)
现在我们考虑如何构造 (G(x)) 使得 (G(F_k)=O_k)
先给出结论，后续提供证明：
( ext{deg }G=k, g_n=egin{cases} -f_{k-n}, n<k \ 1, n=k end{cases})
故原式变换为 (displaystyle a_n=sum_{i=0}^{k-1} r_icdot a_i)
其中，(a_0)~(a_{k-1}) 为题目给定的初值，(r_0)~(r_{k-1}) 为 (x^nmod G(x)) 的多项式取余
对于求解 (x^nmod G(x)) 可以由倍增实现，故复杂度为
无 FFT 优化：(O(k^2log n))；有 FFT 优化：(O(klog kcdot log n))
附 1 ：AC 代码
FFT 版：
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef double db;
typedef pair<int,int> pii;
#define fi first
#define se second
#define de(x) cout<<#x<<" = "<<x<<endl
#define dd(x) cout<<#x<<" = "<<x<<" "
const int LimBit=18, M=1<<LimBit<<1, MOD=998244353;
int a[M], b[M], c[M];
inline int Normal(int n) { return (n%=MOD)>=0?n:n+MOD; }
struct NTT{
int G=3, P=998244353;
int N, na, nb, W[2][M+M], Rev[M+M], *w[2], *rev, invN, InV[M];
inline ll kpow(ll a, ll x) { ll ans=1; for(;x;x>>=1,a=a*a%P) if(x&1) ans=ans*a%P; return ans; }
inline ll inv(ll a) { return kpow(a, P-2); }
inline int add(int a, int b) { return (a+b>=P)?(a+b-P):(a+b); }
inline int dis(int a, int b) { return (a-b<0)?(a-b+P):(a-b); }
inline void display(int *f,int n){
for(int i=0;i<n;++i) cout<<f[i]<<" "; cout<<"
";
}
NTT(){
InV[1] = 1;
for(int i=2;i<P&&i<M;++i)
InV[i]=P-(ll)P/i*InV[P%i]%P;
w[0]=W[0], w[1]=W[1], rev=Rev;
w[0][0]=w[1][0]=1;
for(int i=1, x=kpow(G,P>>LimBit), y=inv(x); !(i>>LimBit); ++i){
rev[i]=(rev[i>>1]>>1)|((i&1)<<LimBit-1);
w[0][i]=(ll)w[0][i-1]*x%P, w[1][i]=(ll)w[1][i-1]*y%P;
}
int *lstw[2]={w[0], w[1]};
w[0]+=1<<LimBit, w[1]+=1<<LimBit, rev+=1<<LimBit;
for(int d=LimBit-1; d>=0; --d){
for(int i=0; !(i>>d); ++i){
rev[i]=(rev[i>>1]>>1)|((i&1)<<d-1);
w[0][i]=lstw[0][i<<1], w[1][i]=lstw[1][i<<1];
}
lstw[0]=w[0], lstw[1]=w[1];
w[0]+=1<<d, w[1]+=1<<d, rev+=1<<d;
}
}
inline void work(){
w[0]=W[0], w[1]=W[1], rev=Rev;
for(int d=LimBit;1<<d>N;--d)
w[0]+=1<<d, w[1]+=1<<d, rev+=1<<d;
invN=inv(N);
}
inline void FFT(int *a,int f){
for(int i=0;i<N;++i) if(i<rev[i]) swap(a[i], a[rev[i]]);
for(int i=1;i<N;i<<=1)
for(int j=0, t=N/(i<<1); j<N; j+=i<<1)
for(int k=0, l=0, x; k<i; ++k, l+=t)
x=(ll)w[f][l]*a[j+k+i]%P, a[j+k+i]=dis(a[j+k], x), a[j+k]=add(a[j+k], x);
if(f) for(int i=0;i<N;++i) a[i]=(ll)a[i]*invN%P;
}
inline void doit(int *a,int *b,int na,int nb,int *c=0){
static int x[M], y[M];
if(!c) c=a;
for(N=1;N<na+nb-1;N<<=1);
memset(x, 0, sizeof(x[0])*N); memset(y, 0, sizeof(y[0])*N);
memcpy(x, a, sizeof(a[0])*na); memcpy(y, b, sizeof(b[0])*nb);
work(); FFT(x, 0); FFT(y, 0);
for(int i=0;i<N;++i) x[i]=(ll)x[i]*y[i]%P;
FFT(x, 1);
memcpy(c, x, sizeof(x[0])*N);
}
inline void get_inv(int *f,int *g,int n){
g[0]=inv(f[0]);
if(n==1) return ;
for(int l=0; (1<<l)<n; ++l){
memcpy(a, f, sizeof(a[0])<<l+1);
memset(a+(1<<l+1), 0, sizeof(a[0])<<l+1);
memset(g+(1<<l), 0, sizeof(g[0])*3<<l);
N=1<<l+2;
work(); FFT(a, 0); FFT(g, 0);
for(int i=0; i<N; ++i) g[i]=(ll)g[i]*dis(2, (ll)g[i]*a[i]%P)%P;
FFT(g, 1);
}
memset(g+n, 0, sizeof(g[0])*(N-n));
}
inline void get_div(int *f,int *g,int *q,int *r, int n,int m){
while(f[n-1]==0&&n) --n;
while(g[m-1]==0&&m) --m;
if(n<m||m==0){
memcpy(r, f, sizeof(f[0])*n);
memset(q, 0, sizeof(q[0])*n);
return ;
}
reverse(g, g+m);
get_inv(g, q, n-m+1);
reverse(f, f+n);
doit(q, f, n-m+1, n-m+1);
reverse(q, q+n-m+1);
reverse(g, g+m);
reverse(f, f+n);
memset(q+n-m+1, 0, sizeof(q[0])*(m-1) );
memcpy(r, g, sizeof(g[0])*m);
doit(r, q, m, n);
for(int i=0;i<m-1;++i) r[i]=dis(f[i], r[i]);
}
inline void get_xpow(int *g,int *m,int x,int nm){
if(x==0){
g[0]=1;
return ;
}
get_xpow(g, m, x>>1, nm);
for(N=1;N<nm+nm-3;N<<=1);
for(int i=nm-1;i<N;++i) g[i]=0;
work(); FFT(g, 0);
for(int i=0;i<N;++i) g[i]=(ll)g[i]*g[i]%MOD;
FFT(g, 1);
if(x&1){
for(int i=N;i>=1;--i) g[i]=g[i-1];
g[0]=0;
get_div(g, m, c, b, N+1, nm);
}
else get_div(g, m, c, b, N, nm);
for(N=1;N<nm+nm-3;N<<=1);
for(int i=nm-1;i<N;++i) g[i]=0;
for(int i=0;i<nm-1;++i) g[i]=b[i];
}
}ntt;
int n, k, f[M], g[M], r[M];
inline int ans(){
int res=0;
for(int i=0;i<k;++i) res=ntt.add(res, (ll)g[i]*r[i]%MOD);
return res;
}
inline void init(){
cin>>n>>k;
for(int i=1;i<=k;++i) cin>>f[k-i], f[k-i]=Normal(-f[k-i]);
for(int i=0;i<k;++i) cin>>g[i], g[i]=Normal(g[i]);
f[k]=1;
}
int main(){
ios::sync_with_stdio(0);
cin.tie(0); cout.tie(0);
init();
ntt.get_xpow(r, f, n, k+1);
cout<<ans();
// cout<<fixed<<setprecision(6);
cout.flush();
return 0;
}
附 2 ：多项式取余的理论
设给定 (F(x), G(x)) ，求解 (F(x)=Q(x)cdot G(x)+R(X))
其中，(Q(x)) 为多项式 (F(x)) 除以 (G(x)) 的商，(R(x)) 是除法的余数
观察可以发现一个很棒的性质：当 ( ext{deg }Fgeq ext{deg }G) 时，有 ( ext{deg }Q= ext{deg }F- ext{deg G}, ext{deg }R= ext{deg }G-1)
这里的 (R(x)) 强制设置为最高可能的位数
故我们考虑将 (x) 换元为 ({1over x}) 且在方程两边同时乘上 (x^{ ext{deg } F}) 得到：
(x^{ ext{deg }F}F({1over x})=x^{ ext{deg }F- ext{deg }G}Q({1over x})cdot x^{ ext{deg }G}G({1over x})+x^{ ext{deg }F-G+1}cdot x^{ ext{deg }G-1}R({1over x}))
让我们再考虑一下 (x^{ ext{deg }F}F({1over x})) 的含义：
原来的 (n) 次项 (f_ncdot x^n) 转变为了 (f_ncdot {1over x^n}cdot x^{ ext{deg }F}=f_ncdot x^{ ext{deg }F-n}) 成为了现在的 (( ext{deg }F-n)) 次项
相当于原来的多项式函数 (F(x)) 的后 ( ext{deg }F) 次项镜像翻转了一下，故我们记这种多项式为 (F^R(x))
代入上面得到 (F^R(x)=Q^R(x)cdot G^R(x)+x^{ ext{deg }F- ext{deg }G+1}cdot R^R(x))
两边同余 (x^{ ext{deg }F- ext{deg }G+1}) 得 (F^R(x)equiv Q^R(x)cdot G^R(x)pmod{x^{ ext{deg }F- ext{deg }G+1}})
移项得到 (Q^R(x)equiv [F^R(x)]^{-1}cdot G^R(x)pmod{x^{ ext{deg }F- ext{deg }G+1}})
通过多项式的取逆，我们可以得到 (Q^R(x)) 的后 (( ext{deg }F- ext{deg }G)) 项；不过因为 ( ext{deg }Q= ext{deg }F- ext{deg }G) ，我们直接得到了 (Q) 的所有项
故两多项式乘积后，直接翻转得到 (Q(x))
再通过多项式乘法和多项式减法求得 (R(x)=F(x)-Q(x)cdot G(x))，复杂度为 (O(nlog n))
附 3 ： (G(x)) 的构造
根据 Cayley-Hamilton 定理，对于 (F_k) 的特征多项式 (p(lambda)) 有 (p(F_k)=O_k) 即 (k) 阶零矩阵
故只要令 (G(x)) 为 (F_k) 的特征多项式即可
现考虑如何求解 (F_k) 的特征多项式：
(p(lambda)=det (F_k-lambdacdot E_k)=det left( egin{matrix} f_1-lambda & f_2 & cdots & f_{k-1} & f_k \ 1 & -lambda \ & 1 & -lambda \ && ddots & ddots \ &&& 1&-lambda end{matrix} ight))
据说有人手动展开后，得到 (p(lambda)=(-1)^ncdot (lambda^n-f_1cdot lambda^{n-1}-f_2cdot lambda^{n-2}-cdots -f_kcdot lambda^0))
由于 (G(F_k)=O_k) ，所以前面那个 ((-1)^n) 也不怎么重要了，直接得到：
(displaystyle G(x)=x^n-sum_{i=1}^kf_ix^{n-i})
于是就有了上面那个式子
拓展
有个神仙算法叫 Berlekamp-Massey 算法
将某个序列的前若干项丢进去，它会贪心地得出满足这些项的最短递推式
设丢进去的项数为 (k) ，则复杂度为 (O(k^2))
结合该算法，可以得出一个非常优（暴）美（力）的算法：
先暴力打表，尽可能多打几项，然后丢进 BM 算法，跑出最短递推式，然后用常系数齐次线性递推算法直接得出第 (n) 项
由于最短递推式最多为 (k) 项，故总复杂度为 (O(k^2+klog klog n)) 或 (O(k^2log n))
相比于拉格朗日插值，优势非常明显：
对于类似卡特兰数 (displaystyle c_n={1over n+1}dbinom{2n}{n}) ，其线性项数不定，拉格朗日插值无从下手
但其最短递推式并不长：(displaystyle c_n={1over n+1}dbinom{2n}{n}={(2n)!over n!(n+1)!}={(2n)cdot (2n-1)over (n+1)cdot n}cdot {(2n-2)!over (n-1)!n!}={4n-2over n+1}cdot {1over n}dbinom{2n-2}{n-1}={4n-2over n+1}c_{n-1})
一项的递推式，显然跑得飞快
板子的话，显然网上有，我就不弄了
题意：给定一条 (n) 元链，链上相邻元素相连，共有 ((n-1)) 条边。链外有 (3) 个点，每个点连向链上的 (n) 个点。整个图共 ((4n-1)) 条边。(nleq 10^9) ，问生成树方案数，答案对 ((10^9+7)) 取模。
当图中点数为 (|V|) 时，根据矩阵树定理：可以建立基尔霍夫矩阵，用高斯消元法 (O(|V|^3)) 求出生成树的方案数。
由于点数为 (n+3) ，故我们可以暴力打出约前 (100) 项的生成树方案数
接着，丢进 BM 算法，跑出最短递推式，发现只有 (6) 项，为：
(displaystyle a_n=sum_{i=1}^6 f_ia_{n-i})
其中，(f={-1, 15, -78, 155, -78, 15})
直接用常系数齐次线性递推求出答案即可