1.求逆
设要对多项式A求逆,逆为B:
求在模ceil(n/2)的逆,将原式与之相减后平方,发现可得在模n意义下也是0了,然后两边同乘A,移项即可倍增了(注意此时A也是在模当前长度意义下的)。
最后倍增形式为B'=2*B-B*B*A (b’等于2b减b方a),边界为常数取逆元,也可发现多项式有无逆取决于常数有无逆元
2.开根
与求逆类似,先求ceil(n/2),然后将A移到左边使等式右边变成0,平方,加上4b^2a,再把右边b^2除到左边就发现ok了。
最终倍增形式:B'=(A+B*B)/2B(B'等于A加B方后除以2B),除法就上求逆,边界是常数开根,暴力枚举(并不会二次剩余...)?模板题保证是1...
(开根里面套求逆,虽说复杂度仍然是nlog,但已经比大部分两个log慢了。。。)
3.Ln
两边求导得到所求的导数B’=A'/A,再积分即可。
4.Exp
B=e^(A)
根据牛顿迭代的那套理论,设一个函数C,使C(B)=0,那么C(B)=LnB-A
假设已知在mod x^n的B0,要求在mod x^2n的B1
有C(B)=C(B0)/0!+C'(B0)/1!*(B-B0)+...
发现...可以省略,因为在mod x^n下B是唯一确定的,而mod x^2n下一定也满足mod x^n下,所以在mod x^2n时后面低的n位与mod x^n相等,当求的导数大于等于二阶发现(B-B0)那一项的幂次大于等于2,也就是最小的x项幂次大于等于2n,省略即可。
最终推得
B1=B0-C(B0)/C'(B0)
而C求导时是在把B0这个多项式看成变量求导(不要把每一项的x看成变量),A是已知多项式要看成一个常数
故B1=B0(1-Ln(B0)+A)
倍增即可。(注意能Exp一定要满足0次项为0,不然无法赋初值)
先贴4个的代码:
#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int mod=998244353;
const int N=1e5+100;
void Ad(int &x,int y)
{if((x+=y)>=mod)x-=mod;}
void Mul(int &x,int y)
{x=1LL*x*y%mod;}
int qpow(int x,int y)
{
int res=1;
while(y)
{
if(y&1)res=1LL*res*x%mod;
x=1LL*x*x%mod,y>>=1;
}
return res;
}
namespace Poly{
vector<int>A,B,C;
int wh[N*3],cc,len,n,inv[N*3];
void Pre_Ntt(int l)
{
cc=0,len=1;
while(len<=l)++cc,len<<=1;
}
void Cal_wh()
{for(int i=1;i<len;i++)wh[i]=(wh[i>>1]>>1)|((i&1)<<(cc-1));}
void Ntt(vector<int>&a,bool inv)
{
for(int i=0;i<len;i++)
if(i<wh[i])swap(a[i],a[wh[i]]);
int tp,mo,ha;
for(int l=2,md;l<=len;l<<=1)
{
md=l>>1,tp=qpow(3,(mod-1)/l);
for(int i=0;i<len;i+=l)
{
mo=1;
for(int j=0;j<md;j++,mo=1LL*mo*tp%mod)
{
ha=1LL*mo*a[i+j+md]%mod;
a[i+j+md]=(a[i+j]-ha+mod)%mod;
a[i+j]=(a[i+j]+ha)%mod;
}
}
}
if(inv)
{
tp=qpow(len,mod-2);
for(int i=1;i<len/2;i++)swap(a[i],a[len-i]);
for(int i=0;i<len;i++)Mul(a[i],tp);
}
}
void Get_Inv(vector<int>&b,int n)
{
if(n==1){b[0]=qpow(A[0],mod-2);return;}
int md=(n+1)/2;
vector<int>a,c;
Pre_Ntt(n+md+md),a.resize(len),c.resize(len);
Get_Inv(c,md);
Pre_Ntt(n+md+md),Cal_wh();
for(int i=0;i<n;i++)a[i]=A[i];
Ntt(a,0),Ntt(c,0);
for(int i=0;i<len;i++)a[i]=(2*c[i]-1LL*c[i]*c[i]%mod*a[i]%mod+mod)%mod;
Ntt(a,1);
for(int i=0;i<n;i++)b[i]=a[i];
}
void Inv(vector<int> &a,vector<int> &res,int n)
{
A=a,res.clear(),res.resize(n);
Get_Inv(res,n);
}
void Get_Sqrt(vector<int>&b,int n)
{
if(n==1){b[0]=1;return;}
int md=(n+1)/2;
vector<int>a,c;
Pre_Ntt(n+n),a.resize(len),c.resize(len);
Get_Sqrt(c,md);
for(int i=0;i<n;i++)b[i]=2*c[i]%mod;
Inv(b,b,n);
Pre_Ntt(n+n),Cal_wh();
Ntt(c,0);
for(int i=0;i<len;i++)c[i]=1LL*c[i]*c[i]%mod;
Ntt(c,1);
for(int i=0;i<n;i++)a[i]=(B[i]+c[i])%mod;
Ntt(a,0),Ntt(b,0);
for(int i=0;i<len;i++)b[i]=1LL*b[i]*a[i]%mod;
Ntt(b,1);
}
void Sqrt(vector<int> &a,vector<int> &res,int n)
{
Pre_Ntt(n+n);
B=a,res.clear(),res.resize(len);
Get_Sqrt(res,n);
}
void Dao(vector<int> &a,vector<int> &res,int n)
{
res.clear(),res.resize(n);
for(int i=1;i<n;i++)
res[i-1]=1LL*a[i]*i%mod;
res[n-1]=0;
}
void Jf(vector<int> &a,vector<int> &res,int n)
{
res.clear(),res.resize(n+1),inv[1]=1;
for(int i=2;i<=n;i++)
inv[i]=1LL*(mod-mod/i)*inv[mod%i]%mod;
for(int i=n-1;i>=0;i--)
res[i+1]=1LL*a[i]*inv[i+1]%mod;
res[0]=0;
}
void Ln(vector<int>&a,vector<int> &res,int n)
{
static vector<int>mo,ha;
Dao(a,mo,n),Inv(a,ha,n);
Pre_Ntt(n+n),Cal_wh();
mo.resize(len),Ntt(mo,0);
ha.resize(len),Ntt(ha,0);
for(int i=0;i<len;i++)mo[i]=1LL*mo[i]*ha[i]%mod;
Ntt(mo,1),Jf(mo,res,n);
}
void Get_Exp(vector<int>&b1,int n)
{
if(n==1){b1[0]=1;return;}
int md=(n+1)/2;
vector<int>b0;
b0.resize(n),b1.resize(n);
Get_Exp(b0,md),Ln(b0,b1,n);
for(int i=0;i<n;i++)b1[i]=((i==0)-b1[i]+C[i]+mod)%mod;
Pre_Ntt(n+n),b0.resize(len),b1.resize(len);
Cal_wh(),Ntt(b0,0),Ntt(b1,0);
for(int i=0;i<len;i++)
b1[i]=1LL*b0[i]*b1[i]%mod;
Ntt(b1,1);
}
void Exp(vector<int>&a,vector<int> &res,int n)
{C=a,res.resize(n),Get_Exp(res,n);}
}
int main()
{
int n;
scanf("%d",&n);
vector<int>a,b;
a.resize(n);
for(int i=0;i<n;i++)
scanf("%d",&a[i]);
Poly::Exp(a,b,n);
for(int i=0;i<n;i++)
printf("%d ",b[i]);
puts("");
}
刷到一道题要多项式取模/除法,补上:
https://www.luogu.org/problemnew/show/P4512
A(x)=B(x)D(x)+R(x)
翻转,即最高n次的F(x) -> x^n*F(1/x)
整个等式操作一下:
x^n*A(1/x)=x^n(B(1/x)D(1/x)+R(1/x))
拆开:
x^nA(1/x)=x^mB(1/x)x^(n-m)D(1/x)+x^(m-1)R(1/x)x^(n-m+1)
换个元:
A'=B'D' mod(x n-m+1)
即可求D反转后的东西,求完转回来即可,R也随之解决
代码:
#include<bits/stdc++.h>
#define ll long long
using namespace std;
const ll mod=998244353;
const int N=1e5+100;
ll qpow(ll x,ll y=mod-2)
{
ll res=1;
while(y)
{
if(y&1)res=res*x%mod;
x=x*x%mod,y>>=1;
}
return res;
}
inline ll Add(ll x,ll y)
{return x+y>=mod?x+y-mod:x+y;}
inline ll Sub(ll x,ll y)
{return x-y<0?x-y+mod:x-y;}
inline ll Mul(ll x,ll y)
{return x*y%mod;}
namespace Poly{
int wh[N<<2],len,cc;
void pre_ntt(int nn)
{
len=1,cc=0;
while(len<nn)len<<=1,++cc;
for(int i=1;i<len;i++)
wh[i]=(wh[i>>1]>>1)|((i&1)<<(cc-1));
}
void ntt(vector<ll> &a,bool inv)
{
for(int i=0;i<len;i++)
if(i<wh[i])swap(a[i],a[wh[i]]);
ll tp,mo,ha;
for(int l=2,md;l<=len;l<<=1)
{
md=l>>1,tp=qpow(3,(mod-1)/l);
for(int i=0;i<len;i+=l)
{
mo=1;
for(int j=0;j<md;j++,mo=Mul(mo,tp))
{
ha=Mul(a[i+j+md],mo);
a[i+j+md]=Sub(a[i+j],ha);
a[i+j]=Add(a[i+j],ha);
}
}
}
if(inv)
{
tp=qpow(len);
for(int i=1;i<len/2;i++)swap(a[i],a[len-i]);
for(int i=0;i<len;i++)a[i]=Mul(a[i],tp);
}
}
vector<ll> Poly_Mul(vector<ll> x,vector<ll> y)
{
pre_ntt(x.size()+y.size());
x.resize(len),y.resize(len);
ntt(x,0),ntt(y,0);
for(int i=0;i<len;i++)
x[i]=Mul(x[i],y[i]);
ntt(x,1);
return x;
}
vector<ll> Inv(vector<ll> x)
{
int lim=1;vector<ll>A,B;
while(lim<x.size())lim<<=1;
for(int i=1;i<=lim;i<<=1)
{
if(i==1){B.resize(1),B[0]=qpow(x[0]);continue;}
A.resize(i);
for(int j=0;j<i&&j<x.size();j++)A[j]=x[j];
pre_ntt(i+B.size());
A.resize(len),B.resize(len);
ntt(A,0),ntt(B,0);
for(int i=0;i<len;i++)
B[i]=(2*B[i]-A[i]*B[i]%mod*B[i]%mod+mod)%mod;
ntt(B,1),B.resize(i);
}
B.resize(x.size());
return B;
}
vector<ll> Div(vector<ll> A,vector<ll> B)
{
int n=A.size()-1,m=B.size()-1;
reverse(A.begin(),A.end());
reverse(B.begin(),B.end());
A.resize(n-m+1),B.resize(n-m+1);
B=Inv(B),A=Poly_Mul(A,B),A.resize(n-m+1);
reverse(A.begin(),A.end());
return A;
}
vector<ll> Poly_Sub(vector<ll> x,vector<ll> y)
{
if(x.size()<y.size())x.resize(y.size());
for(int i=0;i<y.size();i++)
x[i]=Sub(x[i],y[i]);
return x;
}
}
int main()
{
vector<ll>x,y,res;
int n,m;
scanf("%d%d",&n,&m),++n,++m;
x.resize(n),y.resize(m);
for(int i=0;i<n;i++)
scanf("%lld",&x[i]);
for(int i=0;i<m;i++)
scanf("%lld",&y[i]);
res=Poly::Div(x,y);
for(int i=0;i<res.size();i++)
printf("%lld ",res[i]);
puts("");
res=Poly::Poly_Sub(x,Poly::Poly_Mul(y,res));
for(int i=0;i<m-1;i++)
printf("%lld ",res[i]);
puts("");
}