分治FFT:
解决的是形似以下的问题:
给定n次多项式 g ( x ) g(x) g(x),求多项式 f ( x ) f(x) f(x),其中 f f f的第 i i i项系数的表达式为 ∑ f j × g i − j \sum f_{j}\times g_{i-j} ∑fj×gi−j。
解法:
不难发现式子也是卷积的形式,但是与普通多项式乘法不一样的是,每一项的系数依赖前面的项的系数,使得普通的FFT无法起作用。
考虑分治,将区间 [ l , r ] [l,r] [l,r]分为 [ l , m i d ] , [ m i d + 1 , r ] [l,mid],[mid+1,r] [l,mid],[mid+1,r]两个区间计算,计算完 [ l , m i d ] [l,mid] [l,mid]中的多项式的系数之后,可以很方便的将 [ l , m i d ] [l,mid] [l,mid]对 [ m i d + 1 , r ] [mid+1,r] [mid+1,r]的贡献给算出来,贡献的计算同上面的递推式。
由于贡献的计算也需要依赖前面的系数,所以分治时的进程是:
divide(l,mid);
calc([l,mid]->[mid+1,r]);
divide(mid+1,r);
不难发现每个点之前的部分会被拆成若干个区间累加到这个位置上,所以正确性显而易见。
/*=======================================
* Author : ylsoi
* Time : 2019.1.29
* Problem : luogu4721
* E-mail : ylsoi@foxmail.com
* ====================================*/
#include<bits/stdc++.h>
#define REP(i,a,b) for(int i=a,i##_end_=b;i<=i##_end_;++i)
#define DREP(i,a,b) for(int i=a,i##_end_=b;i>=i##_end_;--i)
#define debug(x) cout<<#x<<"="<<x<<" "
#define fi first
#define se second
#define mk make_pair
#define pb push_back
typedef long long ll;
using namespace std;
void File(){
freopen("luogu4721.in","r",stdin);
freopen("luogu4721.out","w",stdout);
}
template<typename T>void read(T &_){
_=0; T fl=1; char ch=getchar();
for(;!isdigit(ch);ch=getchar())if(ch=='-')fl=-1;
for(;isdigit(ch);ch=getchar())_=(_<<1)+(_<<3)+(ch^'0');
_*=fl;
}
const int maxn=1e5+10;
const int mod=998244353;
ll qpow(ll x,ll y){
ll ret=1; x%=mod;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;
y>>=1;
}
return ret;
}
struct ksslbh{
int lim,cnt,dn[maxn<<2];
ll g[maxn<<2],ig[maxn<<2];
void ntt(ll *A,int ty){
REP(i,0,lim-1)if(i<dn[i])swap(A[i],A[dn[i]]);
for(int len=1;len<lim;len<<=1){
ll w= ty==1 ? g[len<<1] : ig[len<<1];
for(int L=0;L<lim;L+=len<<1){
ll wk=1;
REP(i,L,L+len-1){
ll u=A[i],v=A[i+len]*wk%mod;
A[i]=(u+v)%mod;
A[i+len]=(u-v)%mod;
wk=wk*w%mod;
}
}
}
}
void work(int n,int m,ll *a,ll *b){
lim=1,cnt=0;
while(lim<=n+m)lim<<=1,++cnt;
if(!cnt)cnt=1;
REP(i,0,lim){
dn[i]=dn[i>>1]>>1|((i&1)<<(cnt-1));
if(i>n)a[i]=0;
if(i>m)b[i]=0;
}
g[lim]=qpow(3,(mod-1)/lim);
ig[lim]=qpow(g[lim],mod-2);
for(int i=lim>>1;i;i>>=1){
g[i]=g[i<<1]*g[i<<1]%mod;
ig[i]=ig[i<<1]*ig[i<<1]%mod;
}
ntt(a,1),ntt(b,1);
REP(i,0,lim-1)a[i]=a[i]*b[i]%mod;
ntt(a,-1);
ll inv=qpow(lim,mod-2);
REP(i,0,lim-1)a[i]=a[i]*inv%mod;
}
}T;
int n;
ll f[maxn<<2],g[maxn<<2],a[maxn<<2],b[maxn<<2];
void divide(int l,int r){
if(l==r)return;
int mid=(l+r)>>1;
divide(l,mid);
int nn=mid-l,mm=r-l-1;
//nn means the pows of f/a
//mm means the pows of g/b
REP(i,0,nn)a[i]=f[l+i];
REP(i,0,mm)b[i]=g[i+1];
T.work(nn,mm,a,b);
REP(i,mid+1,r)f[i]=(f[i]+a[i-1-l])%mod;
divide(mid+1,r);
}
int main(){
//File();
read(n);
REP(i,1,n-1)read(g[i]);
f[0]=1;
divide(0,n);
REP(i,0,n-1)printf("%lld ",(f[i]+mod)%mod);
return 0;
}
多项式求逆:
对于一个多项式 A ( x ) A(x) A(x),找到一个多项式 A − 1 ( x ) A^{-1}(x) A−1(x),使得 A ( x ) × A − 1 ( x ) A(x)\times A^{-1}(x) A(x)×A−1(x)只有常数项的系数为1,其余项的系数全部为0,那么我们称 A − 1 ( x ) A^{-1}(x) A−1(x)为多项式 A ( x ) A(x) A(x)的逆。
当常数项存在逆元的时候,利用递推的方法去求解多项式的逆,不难发现所得的多项式是一个无限项的级数,所以我们一般只取前面 n − 1 n-1 n−1次,即对 x n x^n xn取模。
暴力求法:
从低次依次向高次递推,每次求解 n n n次项的系数 a n a_n an时,将前 n n n项带入,所得的和需为 0 0 0,直接解出来即可。
优化:
观察上面的递推过程,不难发现后面的项求解需要用到前面的项的信息,利用分治FFT可以优化到 n log 2 n n\log^2n nlog2n。
更加优秀的解法:
直接一项一项地递推显然太浪费时间,考虑利用倍增的思想来利用前面的项推出后面的项。
假设我们现在得到了前 n 2 \frac{n}{2} 2n项的系数,那么我们就要想办法得出前 n n n项的系数。
假设前
n
2
\frac{n}{2}
2n项的多项式为
B
(
x
)
B(x)
B(x),最终的多项式为
A
−
1
(
x
)
A^{-1}(x)
A−1(x),那么我们可以得到
A
−
1
(
x
)
−
B
(
x
)
≡
0
m
o
d
  
x
n
2
A^{-1}(x)-B(x)\equiv 0 \mod x^{\frac{n}{2}}
A−1(x)−B(x)≡0modx2n
为了使模数提升至
x
n
x^{n}
xn,可以将等式两边平方得:
(
A
−
1
(
x
)
−
B
(
x
)
)
2
≡
0
m
o
d
  
x
n
(A^{-1}(x)-B(x))^2\equiv 0 \mod x^{n}
(A−1(x)−B(x))2≡0modxn
考虑证明这一步,任意一个
i
(
i
<
n
)
i (i < n)
i(i<n)次项的系数的拆分必定有一个数
<
n
2
< \frac{n}{2}
<2n,而任意
<
n
2
< \frac{n}{2}
<2n次项的系数都为0,所以
i
(
i
<
n
)
i (i<n)
i(i<n)次项的系数都为
0
0
0。
同时注意这种变化只有在等号右边等于 0 0 0时才成立,这也是为什么要拿两式相减的理由。
然后将平方项拆开:
A
−
2
(
x
)
−
2
×
A
−
1
(
x
)
×
B
(
x
)
+
B
2
(
x
)
≡
0
m
o
d
  
x
n
A^{-2}(x)-2\times A^{-1}(x)\times B(x)+B^{2}(x)\equiv 0 \mod x^n\\
A−2(x)−2×A−1(x)×B(x)+B2(x)≡0modxn
两边同时乘以
A
(
x
)
A(x)
A(x):
A
−
1
(
x
)
−
2
×
B
(
x
)
+
A
(
x
)
×
B
2
(
x
)
≡
0
m
o
d
  
x
n
A
−
1
(
x
)
≡
2
×
B
(
x
)
−
A
(
x
)
×
B
2
(
x
)
m
o
d
  
x
n
A^{-1}(x)-2\times B(x)+A(x)\times B^2(x)\equiv 0 \mod x^n\\ A^{-1}(x)\equiv 2\times B(x)-A(x)\times B^{2}(x) \mod x^n\\
A−1(x)−2×B(x)+A(x)×B2(x)≡0modxnA−1(x)≡2×B(x)−A(x)×B2(x)modxn
然后我们就成功地得到了在模
x
n
x^n
xn意义下的多项式
A
−
1
(
x
)
A^{-1}(x)
A−1(x)。
具体实现:
考虑递归实现上述过程。
实现过程中,理论上 A ( x ) A(x) A(x)和 B ( x ) B(x) B(x)项数都很多,但为了保证复杂度, A ( x ) A(x) A(x)每次取对计算有贡献的前 n n n项,而 B ( x ) B(x) B(x)原本是在模 x n 2 x^\frac{n}{2} x2n意义下的多项式,只取前 n 2 \frac{n}{2} 2n项即可,所以最后乘出来的多项式为 2 n 2n 2n次的,时间复杂度满足 T ( n ) = n log n + T ( n 2 ) T(n)=n\log n+T(\frac{n}{2}) T(n)=nlogn+T(2n)。
所以最终时间复杂度为 n log n n\log n nlogn。
/*=======================================
* Author : ylsoi
* Time : 2019.2.2
* Problem : luogu4238
* E-mail : ylsoi@foxmail.com
* ====================================*/
#include<bits/stdc++.h>
#define REP(i,a,b) for(int i=a,i##_end_=b;i<=i##_end_;++i)
#define DREP(i,a,b) for(int i=a,i##_end_=b;i>=i##_end_;--i)
#define debug(x) cout<<#x<<"="<<x<<" "
#define fi first
#define se second
#define mk make_pair
#define pb push_back
typedef long long ll;
using namespace std;
void File(){
freopen("luogu4238.in","r",stdin);
freopen("luogu4238.out","w",stdout);
}
template<typename T>void read(T &_){
_=0; T fl=1; char ch=getchar();
for(;!isdigit(ch);ch=getchar())if(ch=='-')fl=-1;
for(;isdigit(ch);ch=getchar())_=(_<<1)+(_<<3)+(ch^'0');
_*=fl;
}
const int maxn=1e5+10;
const int mod=998244353;
int lim,cnt,dn[maxn<<2];
ll g[maxn<<2],ig[maxn<<2];
ll qpow(ll x,ll y){
ll ret=1; x%=mod;
while(y){
if(y&1)ret=ret*x%mod;
x=x*x%mod;
y>>=1;
}
return ret;
}
void ntt(ll *A,int ty){
REP(i,0,lim-1)if(i<dn[i])swap(A[i],A[dn[i]]);
for(int len=1;len<lim;len<<=1){
ll w= ty==1 ? g[len<<1] : ig[len<<1];
for(int L=0;L<lim;L+=len<<1){
ll wk=1;
REP(i,L,L+len-1){
ll u=A[i],v=A[i+len]*wk%mod;
A[i]=(u+v)%mod;
A[i+len]=(u-v)%mod;
wk=wk*w%mod;
}
}
}
}
int n;
ll f[maxn<<2],a[maxn<<2],b[maxn<<2];
void init(){
read(n);
REP(i,0,n-1)read(f[i]);
lim=1;
while(lim<=n+n)lim<<=1;
g[lim]=qpow(3,(mod-1)/lim);
ig[lim]=qpow(g[lim],mod-2);
for(int i=lim>>1;i;i>>=1){
g[i]=g[i<<1]*g[i<<1]%mod;
ig[i]=ig[i<<1]*ig[i<<1]%mod;
}
}
void poly_inv(int k){
if(k==1){
b[0]=qpow(f[0],mod-2);
return;
}
poly_inv((k+1)>>1);
lim=1,cnt=0;
while(lim<=k+k)lim<<=1,++cnt;
if(!cnt)cnt=1;
REP(i,0,lim-1)dn[i]=dn[i>>1]>>1|((i&1)<<(cnt-1));
REP(i,0,k-1)a[i]=f[i];
REP(i,k,lim-1)a[i]=0;
ntt(a,1),ntt(b,1);
REP(i,0,lim-1)b[i]=b[i]*(2ll-a[i]*b[i]%mod)%mod;
ntt(b,-1);
ll inv=qpow(lim,mod-2);
REP(i,0,k-1)b[i]=b[i]*inv%mod;
REP(i,k,lim)b[i]=0;
}
int main(){
File();
init();
poly_inv(n);
REP(i,0,n-1)printf("%lld ",(b[i]%mod+mod)%mod);
return 0;
}