多项式运算
前置知识NTT
所有操作均在对(P= ext{998244353})取模下进行
代码在最下面，由于板子实在有一点长，所以。。。
下文中(pmod {x^n})表示求出了多项式的前(n)项
([x^i]F(x))表示(F(x))第(i)项的系数
每个小问题的模板题都可以在洛谷上找到
[
]
1.多项式求乘法逆
(为什么叫做乘法逆?因为还有求(G(x)=frac{1}{F(x)}pmod {M(x)}的))
求 (G(x)equiv frac{1}{F(x)} pmod {x^n})
形象化的理解就是(F(x)cdot G(x) pmod {x^n})只有第一项是(1)，其他项都是(0)
这个由于是第一个操作，很多人还并不是很能理解多项式操作到底是什么东西，所以讲多一点
Part1 O((n^2))
为了便于理解这个问题，先考虑一个最简单的模拟
([x^i]Fcdot G(x)=sum [x^j]F(x)[x^{i-j}]G(x))
第一项([x^0]G(x)=frac{1}{[x^0]F(0)} pmod P)，因此求逆的前提条件是([x^0]F(x) e 0)
考虑从(1)到(n-1)依次求出每一项，先从前面的项中得到所有(j>0)的和(Sum)，然后带入(j=0)时知道
[[x^i]G(x)=-frac{Sum=sum_{j=1}^{jleq i}[x^j]F(x)[x^{i-j}]G(x)}{[x^0]F(0)}
]
[
]
Part2 O((nlog^2n))
上面这个式子是一个类似(dp)转移的东西，可以直接分治NTT优化掉
[
]
Part3 (O(nlog n))
考虑递归求解，设已经求出了
[H(x)equiv frac{1}{F(x)},pmod {x^{frac{n}{2}}}
]
其中递归边界是(n=1)时,([x^0]G(x)=frac{1}{[x^0]F(0)} pmod P)，因此求逆的前提条件是([x^0]F(x) e 0)
则
[H(x)equiv G(x)pmod {x^{frac{n}{2}}}
]
[H(x)-G(x)equiv 0pmod {x^{frac{n}{2}}}
]
我们对于(H(x)-G(x))平方，结果的前(n)项不可能由两个(ge frac{n}{2})的项相乘得到，而前(frac{n}{2})项都是(0)，所以
[(H(x)-G(x))^2equiv 0pmod {x^n}
]
所以通过平方可以扩大模数，这很常用
展开平方的式子
[H(x)^2-2G(x)H(x)+G(x)^2equiv 0pmod {x^n}
]
两边乘上(F(x))
[H(x)^2F(x)-2H(x)+G(x)equiv 0pmod {x^n}
]
[G(x)equiv 2H(x)-H(x)^2F(x)pmod {x^n}
]
带入这个式子倍增求解即可
分析复杂度，每次有一个(H(x)^2F(x))，可以通过(NTT)求出，倍增过程中访问的长度是(O(n+frac{n}{2}+frac{n}{4}...)=O(n))
所以总复杂度就是(O(nlog n))
[
]
2.多项式开根号
求(G(x)^2equiv F(x) pmod {x^n})
同样的，递归求解，设已经求出了，递归边界是(n=1)时,([x^0]G(x)=sqrt{[x^0]F(x)}pmod P)
可以发现我们需要求二次剩余。。。但是一般题目保证了([x^0]F(x)in{0,1})
设已经求出(H(x)^2equiv F(x) pmod{ x^{lceil frac{n}{2} ceil}})
[H(x)equiv G(x) pmod {x^{lceil frac{n}{2}
ceil}}
]
[H(x)^2-2G(x)H(x)+G(x)^2equiv 0pmod {x^n}
]
[H(x)^2-2G(x)H(x)+F(x)equiv 0 pmod {x^n}
]
[G(x)equiv frac{H(x)^2+F(x)}{2H(x)} pmod {x^n}
]
带入这个式子倍增求解即可
复杂度为(O(nlog n))
[
]
3.多项式求(ln)
对$egin{aligned} G(x)equiv ln F(x) pmod {x^n} end{aligned} $ 两边求导，注意这里是复合函数求导!!!
(egin{aligned} G'(x)equiv F'(x)frac{1}{F(x)} pmod {x^n}end{aligned})
求出(G'(x))，然后求原函数即可
复杂度为(O(nlog n))
[
]
4.多项式求exp
多项式求( ext{exp})即求(G(x)=e^{F(x)} mod x^n)
多项式求( ext{exp})常见的解法有两种
CDQ分治+( ext{NTT})
要求(G(x)=e^{F(x)})
式子两边求导(右边要复合函数求导)，(G'(x)=F'(x) e^{F(x)})
也就是说，(G'(x)=F'(x)G(x))
两边同时积分得到(egin{aligned} G(x)=int{F'(x)G(x)}end{aligned})
我们知道，$ [x^i] egin{aligned}int H(x) =frac{ [x^{i-1}]H(x)}{i}end{aligned} $
带入上面的式子得到([x^i]G(x)=frac{sum_{j=0}^{i-1}[x^j]F'(x)cdot [x^{i-1-j}]G(x)}{i})
那么对于这个式子，直接使用分治NTT求解，其复杂为(O(nlog^2 n))
[
]
牛顿迭代
这是一种渐进意义上更优的做法，但实际在(10^6)以下几乎不可能更快，而且代码难写
但是不管平时用不用，牛顿迭代的知识学习一下肯定是最好的
把题目转化为，对于函数(f(G)=ln G-F)
求出在(mod x^n)意义下的零点
其中(f(x)=ln x-c)
考虑迭代求解，设已经求出(H(x)=e^{F(x)}pmod {x^{frac{n}{2}}})
边界条件是([x^0]H(x)=e^{[x^0]F(x)})(由于没有办法求(e^x)在模意义下的值，所以通常必须要满足([x^0]F(x)=0))
带入牛顿迭代的结果
[G=H-frac{f(H)}{f'(H)}=H(F-ln H+1)
]
每次求(ln) 复杂度和( ext{NTT})相同，所以总复杂度为(O(nlog n))
事实上这个还有优化的余地，就是在求(ln)的时候，多项式逆的部分可以同步倍增求出，不需要每次都倍增一下(但是好像效果并不是特别明显)
[
]
[
]
5.多项式(k)次幂
(G(x)equiv F(x)^kpmod {x^n})
(ln G(x)=k ln F(x) pmod {x^n})
求出(ln G(x))之后，(exp)回来即可
由于要求(ln)，所以这样求的条件是([x^0]F(x)=1)
很显然这个方法对于开根号也是适用的
复杂度(O(nlog n))
[
]
[
]
[
]
6.多项式带余除法
可以参考神仙miskcoo的博客
[
]
以上是基本运算，如果不想继续吸多项式请直接跳到最下面的代码
[
]
[
]
[
]
[
]
[
]
所有的操作均用$ ext{vector} $来实现，主要是为了理清思路，并且清零问题上会比较容易解决，同时对于每次计算完多项式的长度的要求会显得更加严格
稍微整理了一下，没怎么卡过常，所以应该还是比较可读的
代码总览(请使用C++11,O2编译运行)
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef double db;
typedef unsigned long long ull;
typedef pair <int,int> Pii;
#define reg register
#define pb push_back
#define mp make_pair
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))
#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)
template <class T> inline void cmin(T &a,T b){ ((a>b)&&(a=b)); }
template <class T> inline void cmax(T &a,T b){ ((a<b)&&(a=b)); }
char IO;
template <class T=int> T rd(){
T s=0; int f=0;
while(!isdigit(IO=getchar())) if(IO=='-') f=1;
do s=(s<<1)+(s<<3)+(IO^'0');
while(isdigit(IO=getchar()));
return f?-s:s;
}
const int N=1<<17,P=998244353;
int n,k;
ll qpow(ll x,ll k=P-2) {
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}
/*
void NTT(int n,int *a,int f){
rep(i,1,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]);
for(int i=1;i<n;i<<=1) {
int w=qpow(3,(P-1)/i/2);
for(int l=0;l<n;l+=i*2) {
int e=1;
for(int j=l;j<l+i;++j,e=1ll*e*w%P) {
int t=1ll*a[j+i]*e%P;
a[j+i]=a[j]-t,((a[j+i]<0)&&(a[j+i]+=P));
a[j]+=t,((a[j]>=P)&&(a[j]-=P));
}
}
}
if(f==-1) {
reverse(a+1,a+n);
int Inv=qpow(n);
rep(i,0,n-1) a[i]=1ll*a[i]*Inv%P;
}
}
int e[N];
void NTT(int n,int *a,int f){
rep(i,1,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]);
e[0]=1;
for(int i=1;i<n;i<<=1) {
int w=qpow(3,(P-1)/i/2);
for(int j=1;j<i;++j) e[j]=1ll*e[j-1]*w%P;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,((a[j+i]<0)&&(a[j+i]+=P));
a[j]+=t,((a[j]>=P)&&(a[j]-=P));
}
}
}
if(f==-1) {
reverse(a+1,a+n);
int Inv=qpow(n);
rep(i,0,n-1) a[i]=1ll*a[i]*Inv%P;
}
}
int e[N];
void NTT(int n,int *a,int f){
rep(i,1,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]);
e[0]=1;
for(int i=1;i<n;i<<=1) {
int w=qpow(3,(P-1)/i/2);
for(int j=i-2;j>=0;j-=2) e[j+1]=1ll*w*(e[j]=e[j>>1])%P;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,((a[j+i]<0)&&(a[j+i]+=P));
a[j]+=t,((a[j]>=P)&&(a[j]-=P));
}
}
}
if(f==-1) {
reverse(a+1,a+n);
int Inv=qpow(n);
rep(i,0,n-1) a[i]=1ll*a[i]*Inv%P;
}
}
int w[N];
void Init(int N){
w[N>>1]=1;
int t=qpow(3,(P-1)/N);
rep(i,(N>>1)+1,N-1) w[i]=1ll*w[i-1]*t%P;
drep(i,(N>>1)-1,1) w[i]=w[i<<1];
}
void NTT(int n,int *a,int f){
rep(i,1,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]);
for(int i=1;i<n;i<<=1) {
int *e=w+i;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,((a[j+i]<0)&&(a[j+i]+=P));
a[j]+=t,((a[j]>=P)&&(a[j]-=P));
}
}
}
if(f==-1) {
reverse(a+1,a+n);
int Inv=qpow(n);
rep(i,0,n-1) a[i]=1ll*a[i]*Inv%P;
}
}
*/
/*
namespace MTT{
const double PI=acos((double)-1);
int rev[N];
struct Cp{
double x,y;
Cp(){ ; }
Cp(double _x,double _y): x(_x),y(_y){ }
inline Cp operator + (const Cp &t) const { return (Cp){x+t.x,y+t.y}; }
inline Cp operator - (const Cp &t) const { return (Cp){x-t.x,y-t.y}; }
inline Cp operator * (const Cp &t) const { return (Cp){x*t.x-y*t.y,x*t.y+y*t.x}; }
}A[N],B[N],C[N],w[N/2];
#define E(x) ll(x+0.5)%P
void FFT(int n,Cp *a,int f){
rep(i,0,n-1) if(rev[i]<i) swap(a[i],a[rev[i]]);
w[0]=Cp(1,0);
for(reg int i=1;i<n;i<<=1) {
Cp t=Cp(cos(PI/i),f*sin(PI/i));
for(reg int j=i-2;j>=0;j-=2) w[j+1]=t*(w[j]=w[j>>1]);
// 上面提到的最优板子
for(reg int l=0;l<n;l+=2*i) {
for(reg int j=l;j<l+i;j++) {
Cp t=a[j+i]*w[j-l];
a[j+i]=a[j]-t;
a[j]=a[j]+t;
}
}
}
if(f==-1) rep(i,0,n-1) a[i].x/=n,a[i].y/=n;
}
void Multiply(int n,int m,int *a,int *b,int *res,int P){
// [0,n-1]*[0,m-1]->[0,n+m-2]
int S=(1<<15)-1;
int R=1,cc=-1;
while(R<=n+m-1) R<<=1,cc++;
rep(i,1,R) rev[i]=(rev[i>>1]>>1)|((i&1)<<cc);
rep(i,0,n-1) A[i]=Cp((a[i]&S),(a[i]>>15));
rep(i,0,m-1) B[i]=Cp((b[i]&S),(b[i]>>15));
rep(i,n,R-1) A[i]=Cp(0,0);
rep(i,m,R-1) B[i]=Cp(0,0);
FFT(R,A,1),FFT(R,B,1);
rep(i,0,R-1) {
int j=(R-i)%R;
C[i]=Cp((A[i].x+A[j].x)/2,(A[i].y-A[j].y)/2)*B[i];
B[i]=Cp((A[i].y+A[j].y)/2,(A[j].x-A[i].x)/2)*B[i];
}
FFT(R,C,-1),FFT(R,B,-1);
rep(i,0,n+m-2) {
ll a=E(C[i].x),b=E(C[i].y),c=E(B[i].x),d=E(B[i].y);
res[i]=(a+((b+c)<<15)+(d<<30))%P;
}
}
#undef E
}
*/
namespace Polynomial{
typedef vector <int> Poly;
void Show(Poly a,int k=0){
if(!k){ for(int i:a) printf("%d ",i); puts(""); }
else for(int i:a) printf("%d
",i);
}
int rev[N],w[N];
int Inv[N+1],Fac[N+1],FInv[N+1];
void Init_w() {
int t=qpow(3,(P-1)/N);
w[N>>1]=1;
rep(i,(N>>1)+1,N-1) w[i]=1ll*w[i-1]*t%P;
drep(i,(N>>1)-1,1) w[i]=w[i<<1];
Inv[0]=Inv[1]=Fac[0]=Fac[1]=FInv[0]=FInv[1]=1;
rep(i,2,N) {
Inv[i]=1ll*(P-P/i)*Inv[P%i]%P;
FInv[i]=1ll*FInv[i-1]*Inv[i]%P;
Fac[i]=1ll*Fac[i-1]*i%P;
}
}
int Init(int n){
int R=1,c=-1;
while(R<n) R<<=1,c++;
rep(i,1,R-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<c);
return R;
}
#define NTTVersion1
#ifdef NTTVersion1
void NTT(int n,Poly &a,int f){
rep(i,0,n-1) if(rev[i]<i) swap(a[i],a[rev[i]]);
for(int i=1;i<n;i<<=1) {
int *e=w+i;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,Mod2(a[j+i]);
a[j]+=t,Mod1(a[j]);
}
}
}
if(f==-1) {
reverse(a.begin()+1,a.begin()+n);
ll base=Inv[n];
rep(i,0,n-1) a[i]=a[i]*base%P;
}
}
void NTT(int n,int *a,int f){
rep(i,0,n-1) if(rev[i]<i) swap(a[i],a[rev[i]]);
for(int i=1;i<n;i<<=1) {
int *e=w+i;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,Mod2(a[j+i]);
a[j]+=t,Mod1(a[j]);
}
}
}
if(f==-1) {
reverse(a+1,a+n);
ll base=Inv[n];
rep(i,0,n-1) a[i]=a[i]*base%P;
}
}
#else
void NTT(int n,Poly &a,int f){
static int e[N>>1];
rep(i,0,n-1) if(rev[i]<i) swap(a[i],a[rev[i]]);
e[0]=1;
for(int i=1;i<n;i<<=1) {
int t=qpow(f==1?3:(P+1)/3,(P-1)/i/2);
for(int j=i-2;j>=0;j-=2) e[j+1]=1ll*t*(e[j]=e[j>>1])%P;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,Mod2(a[j+i]);
a[j]+=t,Mod1(a[j]);
}
}
}
if(f==-1) {
ll base=Inv[n];
rep(i,0,n-1) a[i]=a[i]*base%P;
}
}
void NTT(int n,int *a,int f){
static int e[N>>1];
rep(i,0,n-1) if(rev[i]<i) swap(a[i],a[rev[i]]);
e[0]=1;
for(int i=1;i<n;i<<=1) {
int t=qpow(f==1?3:(P+1)/3,(P-1)/i/2);
for(int j=i-2;j>=0;j-=2) e[j+1]=1ll*t*(e[j]=e[j>>1])%P;
for(int l=0;l<n;l+=i*2) {
for(int j=l;j<l+i;++j) {
int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,Mod2(a[j+i]);
a[j]+=t,Mod1(a[j]);
}
}
}
if(f==-1) {
ll base=Inv[n];
rep(i,0,n-1) a[i]=a[i]*base%P;
}
}
#endif
Poly operator * (Poly a,Poly b){
int n=a.size()+b.size()-1;
int R=Init(n);
a.resize(R),b.resize(R);
NTT(R,a,1),NTT(R,b,1);
rep(i,0,R-1) a[i]=1ll*a[i]*b[i]%P;
NTT(R,a,-1);
a.resize(n);
return a;
}
Poly operator + (Poly a,Poly b) {
int n=max(a.size(),b.size());
a.resize(n),b.resize(n);
rep(i,0,n-1) a[i]+=b[i],Mod1(a[i]);
return a;
}
Poly operator - (Poly a,Poly b) {
int n=max(a.size(),b.size());
a.resize(n),b.resize(n);
rep(i,0,n-1) a[i]-=b[i],Mod2(a[i]);
return a;
}
Poly Poly_Inv(Poly a) { // 多项式乘法逆，注意这里求出的是前a.size()项
int n=a.size();
if(n==1) return Poly{(int)qpow(a[0],P-2)};
Poly b=a; b.resize((n+1)/2); b=Poly_Inv(b);
int R=Init(n<<1);
a.resize(R),b.resize(R);
NTT(R,a,1),NTT(R,b,1);
rep(i,0,R-1) a[i]=(2-1ll*a[i]*b[i]%P+P)*b[i]%P;
NTT(R,a,-1);
a.resize(n);
return a;
}
Poly operator / (Poly a,Poly b){ // 多项式带余除法
reverse(a.begin(),a.end()),reverse(b.begin(),b.end());
int n=a.size(),m=b.size();
a.resize(n-m+1),b.resize(n-m+1),b=Poly_Inv(b);
a=a*b,a.resize(n-m+1);
reverse(a.begin(),a.end());
return a;
}
Poly operator % (Poly a,Poly b) { // 多项式取模
int n=b.size()-1;
if((int)a.size()<=n) return a;
Poly t=a/b;
if((int)t.size()>n) t.resize(n);
t=t*b; t.resize(n); a.resize(n);
return a-t;
}
int Quad(int a,int k=0) { // 二次剩余(不是原根法),用于求Sqrt
if(a<=1) return a;
ll x;
while(1) {
x=1ll*rand()*rand()%P;
if(qpow((x*x-a+P)%P,(P-1)/2)!=1) break;
}
ll w=(x*x-a+P)%P;
Pii res=mp(1,0),t=mp(x,1);
auto Mul=[&](Pii a,Pii b){
int x=(1ll*a.first*b.first+1ll*a.second*b.second%P*w)%P,y=(1ll*a.first*b.second+1ll*a.second*b.first)%P;
return mp(x,y);
};
int d=(P+1)/2;
while(d) {
if(d&1) res=Mul(res,t);
t=Mul(t,t);
d>>=1;
}
ll r=(res.first%P+P)%P;
if(k) r=min(r,(P-r)%P);
return r;
}
Poly Sqrt(Poly a){ // 多项式开根号
int n=a.size();
if(n==1) return Poly{Quad(a[0],1)};
Poly b=a; b.resize((n+1)/2),b=Sqrt(b),b.resize(n);
Poly c=Poly_Inv(b);
int R=Init(n*2);
a.resize(R),c.resize(R);
NTT(R,a,1),NTT(R,c,1);
rep(i,0,R-1) a[i]=1ll*a[i]*c[i]%P;
NTT(R,a,-1);
a.resize(n);
rep(i,0,n-1) a[i]=1ll*(P+1)/2*(a[i]+b[i])%P;
return a;
}
Poly Deri(Poly a){ //求导
rep(i,1,a.size()-1) a[i-1]=1ll*i*a[i]%P;
a.pop_back();
return a;
}
Poly IDeri(Poly a) { //原函数
a.pb(0);
drep(i,a.size()-1,1) a[i]=1ll*a[i-1]*Inv[i]%P;
a[0]=0;
return a;
}
Poly Ln(Poly a){ // 多项式求Ln
int n=a.size();
a=Poly_Inv(a)*Deri(a),a.resize(n-1);
return IDeri(a);
}
Poly Exp(Poly a){ // 多项式Exp
int n=a.size();
if(n==1) return Poly{1};
Poly b=a; b.resize((n+1)/2),b=Exp(b); b.resize(n);
Poly c=Ln(b);
rep(i,0,n-1) c[i]=a[i]-c[i],Mod2(c[i]);
c[0]++,b=b*c;
b.resize(n);
return b;
}
void Exp_Solve(Poly &A,Poly &B,int l,int r){
static int X[N],Y[N];
if(l==r) {
B[l]=1ll*B[l]*Inv[l]%P;
return;
}
int mid=(l+r)>>1;
Exp_Solve(A,B,l,mid);
int R=Init(r-l+2);
rep(i,0,R) X[i]=Y[i]=0;
rep(i,l,mid) X[i-l]=B[i];
rep(i,0,r-l-1) Y[i+1]=A[i];
NTT(R,X,1),NTT(R,Y,1);
rep(i,0,R-1) X[i]=1ll*X[i]*Y[i]%P;
NTT(R,X,-1);
rep(i,mid+1,r) B[i]+=X[i-l],Mod1(B[i]);
Exp_Solve(A,B,mid+1,r);
}
Poly CDQ_Exp(Poly F){
int n=F.size(); F=Deri(F);
Poly A(n);
A[0]=1;
Exp_Solve(F,A,0,n-1);
return A;
}
Poly Pow(Poly x,int k) { // 多项式k次幂
x=Ln(x);
rep(i,0,x.size()-1) x[i]=1ll*x[i]*k%P;
return Exp(x);
}
Poly EvaluateTemp[N<<1];
void EvaluateSolve1(Poly &a,int l,int r,int p=1){
if(l==r) { EvaluateTemp[p]=Poly{P-a[l],1}; return; }
int mid=(l+r)>>1;
EvaluateSolve1(a,l,mid,p<<1),EvaluateSolve1(a,mid+1,r,p<<1|1);
EvaluateTemp[p]=EvaluateTemp[p<<1]*EvaluateTemp[p<<1|1];
}
void EvaluateSolve2(Poly &res,Poly F,int l,int r,int p=1){
if(l==r){ res[l]=F[0]; return; }
int mid=(l+r)>>1;
EvaluateSolve2(res,F%EvaluateTemp[p<<1],l,mid,p<<1);
EvaluateSolve2(res,F%EvaluateTemp[p<<1|1],mid+1,r,p<<1|1);
}
Poly Evaluate(Poly a,Poly b,int flag=1){ // 多项式多点求值
Poly res(b.size());
if(flag) EvaluateSolve1(b,0,b.size()-1);
EvaluateSolve2(res,a,0,b.size()-1);
return res;
}
Poly InterpolationSolve(Poly &T,int l,int r,int p=1){
if(l==r) return Poly{T[l]};
int mid=(l+r)>>1;
return InterpolationSolve(T,l,mid,p<<1)*EvaluateTemp[p<<1|1]+InterpolationSolve(T,mid+1,r,p<<1|1)*EvaluateTemp[p<<1];
}
Poly Interpolation(Poly X,Poly Y){ // 多项式快速插值
int n=X.size();
EvaluateSolve1(X,0,n-1);
Poly T=Evaluate(Deri(EvaluateTemp[1]),X,0);
rep(i,0,n-1) T[i]=Y[i]*qpow(T[i])%P;
return InterpolationSolve(T,0,n-1);
}
void FFPTrans(Poly &a,int f){ // FFP<->EGF
int n=a.size();
Poly b(n);
if(f==1) rep(i,0,n-1) b[i]=FInv[i];
else rep(i,0,n-1) b[i]=(i&1)?P-FInv[i]:FInv[i];
a=a*b; a.resize(n);
}
Poly FFPMul(Poly a,Poly b){ // FFP卷积
int n=a.size()+b.size()-1;
a.resize(n),b.resize(n);
FFPTrans(a,1),FFPTrans(b,1);
rep(i,0,n-1) a[i]=1ll*a[i]*b[i]%P*Fac[i]%P;
FFPTrans(a,-1);
return a;
}
Poly PolyToFFP(Poly F){ // 多项式转FFP
int n=F.size();
Poly G(n);
rep(i,0,n-1) G[i]=i;
G=Evaluate(F,G);
rep(i,0,n-1) F[i]=1ll*G[i]*FInv[i]%P;
FFPTrans(F,-1);
return F;
}
Poly FFPToPoly(Poly F){ // FFP转多项式
FFPTrans(F,1);
int n=F.size(); Poly X(n);
rep(i,0,n-1) X[i]=i,F[i]=1ll*F[i]*Fac[i]%P;
EvaluateSolve1(X,0,n-1);
rep(i,0,n-1) {
F[i]=1ll*F[i]*FInv[i]%P*FInv[n-i-1]%P;
if((n-i-1)&1) F[i]=(P-F[i])%P;
}
return InterpolationSolve(F,0,n-1);
}
}
using namespace Polynomial;
Poly Lag(int n,Poly X,Poly Y){
Poly T(n+1),R(n+1),A(n+1);
T[0]=1;
rep(i,0,n) drep(j,i+1,0) T[j]=(1ll*T[j]*(P-X[i])+(j?T[j-1]:0))%P;
rep(i,0,n) {
ll t=1;
rep(j,0,n) if(i!=j) t=t*(X[i]-X[j]+P)%P;
t=qpow(t)*Y[i]%P,R[n+1]=T[n+1];
drep(j,n,0) A[j]=(A[j]+t*R[j+1])%P,R[j]=(T[j]+1ll*R[j+1]*X[i]%P+P)%P;
}
return A;
}
int main(){
int n=rd();
Init_w();
Poly F(n);
rep(i,0,n-1) F[i]=rd();
Show(CDQ_Exp(F));
}
[
]
[
]
[
]