多项式全家桶——Part.3 多项式求逆、除法、开根号 多项式全家桶正式进入正片。 多项式求逆 多项式除法 多项式开根
有一说一，感觉之前写的其实都是写不重要或特别基础的东东。
现在才是一些比较有难度的东东。
多项式求逆
介绍
多项式求逆即为求一个多项式的乘法逆元。
它长这样：
然后多项式(B(x))就叫做多项式(A(x))的逆元。
这玩意和普通的乘法逆元有一点不同的就是，它的模是一个(x^n)，其意义为把结果中次数大于等于n的项都给忽略掉。
过程
引理一： 考虑一个乘法逆元满足：(A(x)*B(x)equiv 1(mod x^n))，那么我们可以得到另一个柿子：(A(x)*B(x)equiv 1(mod x^m))其中(1<=m<=n)
证明？把卷积形式写出来就知道了。
考虑设(B(x))表示(A(x))在模(x^n)意义下的逆元。
设(G(x))表示(A(x))在模(x^{lceil frac n 2
ceil})意义下的逆元。
由引理得：
开方（这一步我也不能严谨证明，但是打表发现这是对的）
乘个(A(x))进去，然后由于是模(x^n)意义下，所以可以消去一些(B(x)):
这时，就可以用(G(x))来求出(B(x))。具体过程就递归求解即可（当然不觉得麻烦的话打分治FFT也是可以），中间需要用到NTT。
学习资料：
https://blog.****.net/a_forever_dream/article/details/102483602
https://www.cnblogs.com/LanrTabe/p/11359952.html
https://www.cnblogs.com/guoshaoyang/p/11296025.html
代码
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
using namespace std;
const long long mo=998244353;
const int maxn=300010;
long long a[maxn],b[maxn],rev,A[maxn],G[maxn];
int jl,n,rec[maxn],w[maxn],up,dep;
long long qsm(long long a,long long b)
{
long long t=1;
long long y=a;
while (b>0)
{
if ((b&1)==1) t=t*y%mo;
y=y*y%mo;
b/=2;
}
return t;
}
void FFT(long long a[],int n,int inv)
{
if (n==1) return;
int mid=n/2;
for (int i=0;i<n;i++)
{
rec[i]=(rec[i>>1]>>1)|((i&1)<<(dep-1));
if (i<rec[i]) swap(a[i],a[rec[i]]);
}
for (int len=2;len<=n;len<<=1)
{
int mid=len/2;
for (int j=0;j<mid;j++)
{
for (int k=j;k<n;k+=len)
{
int p,q;
q=w[j*(n/len)]*a[k+mid]%mo,p=a[k];
a[k+mid]=(p-q+mo)%mo;
a[k]=(p+q+mo)%mo;
}
}
}
}
void solve(int n)
{
if (n==1)
{
b[0]=qsm(a[0],mo-2);
}
else
{
double op=n;
int oq=ceil(op/2);
solve(oq);
up=1;dep=0;
while (up<=n+jl) up=up*2,dep++;
w[0]=1;
rev=qsm(3,(mo-1)/up);
for (int i=1;i<=up;i++) w[i]=w[i-1]*rev%mo;
for (int i=0;i<up;i++)
{
G[i]=0;
A[i]=a[i];
if (i<oq) G[i]=b[i];
}
FFT(A,up,1);
FFT(G,up,1);
for (int i=0;i<up;i++)
{
A[i]=G[i]*(2-A[i]*G[i]%mo+mo)%mo;
}
for (int i=0;i<=up/2;i++)
{
swap(w[i],w[up-i]);
}
FFT(A,up,-1);
for (int i=0;i<n;i++) b[i]=A[i]*qsm(up,mo-2)%mo;
int kk=0;
}
}
int main()
{
scanf("%d",&n);
jl=n;
for (int i=0;i<n;i++)
{
scanf("%lld",&a[i]);
}
solve(n);
for (int i=0;i<n;i++)
{
printf("%lld ",b[i]);
}
printf("
");
}
多项式除法
介绍
这玩意是求一个长成这样的柿子的：
其中(A(x))和(B(x))已经给出，且(A(x))次数为n，(B(x))次数为m，求(D(x))和(R(x))
也就是除法运算中求商和余数，话说如果上过初中就应该明白长除法。
过程
引理二：
考虑一个反转操作为将一个多项式的系数调换，用数学形式表达即为：
证明显然。
那么现在就推柿子：
然后再引入引理：
由于我们的(Dr(x))的次数是小于等于(n-m)的，所以如果对原式取模只要大于之即可，而(A(x))、(B(x))都已知，所以不必理会。
然后我们就可以直接求(Br(x))的逆元，再乘(Ar(x))即可得到(Dr(x))，再把系数反过来即可求出(D(x))。
当然，求出(D(x))之后，再求(R(x))就很简单了。
学习资料
https://www.cnblogs.com/knife-rose/p/12056617.html
https://www.cnblogs.com/owenyu/p/6724611.html
代码
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
using namespace std;
const long long mo=998244353;
const int maxn=300010;
long long a[maxn],b[maxn],ar[maxn],br[maxn],d[maxn],D[maxn],r[maxn],rev,A[maxn],G[maxn];
int jl,n,m,rec[maxn],w[maxn],up,dep;
long long qsm(long long a,long long b)
{
long long t=1;
long long y=a;
while (b>0)
{
if ((b&1)==1) t=t*y%mo;
y=y*y%mo;
b/=2;
}
return t;
}
void FFT(long long a[],int n,int inv)
{
if (n==1) return;
int mid=n/2;
for (int i=0;i<n;i++)
{
rec[i]=(rec[i>>1]>>1)|((i&1)<<(dep-1));
if (i<rec[i]) swap(a[i],a[rec[i]]);
}
for (int len=2;len<=n;len<<=1)
{
int mid=len/2;
for (int j=0;j<mid;j++)
{
for (int k=j;k<n;k+=len)
{
int p,q;
q=w[j*(n/len)]*a[k+mid]%mo,p=a[k];
a[k+mid]=(p-q+mo)%mo;
a[k]=(p+q+mo)%mo;
}
}
}
}
void solve(int n)
{
if (n==1)
{
d[0]=qsm(br[0],mo-2);
}
else
{
double op=n;
int oq=ceil(op/2);
solve(oq);
up=1;dep=0;
while (up<=n+jl) up=up*2,dep++;
w[0]=1;
rev=qsm(3,(mo-1)/up);
for (int i=1;i<=up;i++) w[i]=w[i-1]*rev%mo;
for (int i=0;i<up;i++)
{
G[i]=0;
A[i]=br[i];
if (i<oq) G[i]=d[i];
}
FFT(A,up,1);
FFT(G,up,1);
for (int i=0;i<up;i++)
{
A[i]=G[i]*(2-A[i]*G[i]%mo+mo)%mo;
}
for (int i=0;i<=up/2;i++)
{
swap(w[i],w[up-i]);
}
FFT(A,up,-1);
for (int i=0;i<n;i++) d[i]=A[i]*qsm(up,mo-2)%mo;
int kk=0;
}
}
int main()
{
scanf("%d%d",&n,&m);
for (int i=0;i<=n;i++)
{
scanf("%d",&a[i]);
ar[n-i]=a[i];
}
for (int i=0;i<=m;i++)
{
scanf("%d",&b[i]);
br[m-i]=b[i];
}
jl=m;
solve(n-m+1);
up=1;dep=0;
while (up<=2*n-m) up=up*2,dep++;
w[0]=1;
rev=qsm(3,(mo-1)/up);
for (int i=1;i<=up;i++) w[i]=w[i-1]*rev%mo;
FFT(ar,up,1);
FFT(d,up,1);
for (int i=0;i<up;i++)
{
d[i]=d[i]*ar[i]%mo;
}
for (int i=0;i<=up/2;i++)
{
swap(w[i],w[up-i]);
}
FFT(d,up,-1);
for (int i=0;i<=n-m;i++) d[i]=d[i]*qsm(up,mo-2)%mo;
int j=0;
for (int i=n-m;i>=0;i--)
{
printf("%d ",d[i]);
D[j++]=d[i];
}
printf("
");
up=1;dep=0;
while (up<=n) up=up*2,dep++;
w[0]=1;
rev=qsm(3,(mo-1)/up);
for (int i=1;i<=up;i++) w[i]=w[i-1]*rev%mo;
FFT(b,up,1);
FFT(D,up,1);
for (int i=0;i<up;i++)
{
b[i]=b[i]*D[i]%mo;
}
for (int i=0;i<=up/2;i++)
{
swap(w[i],w[up-i]);
}
FFT(b,up,-1);
for (int i=0;i<n;i++) r[i]=(a[i]-b[i]*qsm(up,mo-2)%mo+mo)%mo;
for (int i=0;i<m;i++)
{
printf("%d ",r[i]);
}
printf("
");
}
多项式开根
介绍
这玩意是求一个这样的东东：
其中(A(x))次数为(n-1)，要求一个(B(x))。注意！由于是模意义下的多项式，所以我们的(B(x))的次数不一定为(frac{n-1}2)。（废话）
同样考虑用什么递归分治的思想来搞。
前置芝士
这里就需要一点额外的芝士了。
就这题而言，有极大的概率会用到二次剩余。
那么来康康怎么搞。
首先，二次剩余的定义为：
再定义一个东东：
这个东东叫做勒让德符号（听说还有个叫做高斯记号的东东）。当p为奇质数的时候（(p=2)的时候每个整数就都是其二次剩余了），还有个性质：
那么我们再引入两个东东：（现在所有的p都为奇质数）
定理一： 对于(x^2equiv n(mod p))中有(frac {p-1}2)个n是p的二次剩余。
证明？我不费！
定理二： 首先我们根据定理一，可以考虑random一个数(ain[0,p-1])。设(w=a^2-n)。若找到一个(w)满足：((frac w p)=-1)或者写成这样：(w^{frac{p-1} 2}equiv -1)，则((a+sqrt w)^{frac{p+1} 2})为一组x值。
证明？这个我费！
在模意义下这玩意只有当(i=0)或(i=p)时有值。
那么得到：
根据费马小定理：(a^{p-1}equiv 1)
根据条件：(w^{frac{p+1} 2}equiv 1)
所以可以得到：
那么把原式乘上这玩意即可得到：
证毕。
那么得到定理后，我们可以去计算x值了。
一点注意的地方就是在跑快速幂时，我们需要分实数部分和整数部分来计算，因为(sqrt w)我们并不能很好表示，所以要带进去一起算。
代码：待填。QWQ
过程
n=1
既然我们知道了二次剩余，那么我们可以快速得到(n=1)时的情况。
其实有时候并不需要什么二次剩余，因为可能良心出题人给的(A0)是个比较好看的数，比如说1。
n>1
现在才是正轨。
考虑原式：
用同样的套路设(G(x))表示(A(x))在模(x^{lceil frac n 2
ceil})意义下的开方。
因此得到：
再看到原式，可以得到：
愉快递归求解吧~
学习资料：
百度百科
https://blog.****.net/stevensonson/article/details/85845334
https://www.cnblogs.com/zhangleo/p/11010302.html
https://www.cnblogs.com/xiefengze1/p/9373272.html
代码
话说递归版的常数是真滴大，洛谷反正我是TLE了QWQ。
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
using namespace std;
const long long mo=998244353;
const long long inv2=499122177;
const int maxn=300010;
int a[maxn],b[maxn],d[maxn],rev,A[maxn],G[maxn],D[maxn];
int jl,n,rec[maxn],w[maxn],up,dep;
int jia(int a,int b)
{
return (a+b)%mo;
}
int cheng(int a,int b)
{
return 1ll*a*b%mo;
}
inline long long qsm(long long a,long long b)
{
long long t=1;
long long y=a;
while (b>0)
{
if ((b&1)==1) t=t*y%mo;
y=y*y%mo;
b/=2;
}
return t;
}
inline void FFT(int a[],int n,int inv)
{
if (n==1) return;
int mid=n/2;
for (int i=0;i<n;i++)
{
rec[i]=(rec[i>>1]>>1)|((i&1)<<(dep-1));
if (i<rec[i]) swap(a[i],a[rec[i]]);
}
for (int len=2;len<=n;len<<=1)
{
int mid=len/2;
for (int j=0;j<mid;j++)
{
for (int k=j;k<n;k+=len)
{
int p,q;
q=cheng(w[j*(n/len)],a[k+mid]),p=a[k];
a[k+mid]=jia(p+mo,-q)%mo;
a[k]=jia(p,q)%mo;
}
}
}
}
inline void solve(int n)
{
if (n==1)
{
d[0]=qsm(b[0],mo-2);
}
else
{
double op=n;
int oq=ceil(op/2);
solve(oq);
up=1;dep=0;
while (up<=n+jl) up=up*2,dep++;
w[0]=1;
rev=qsm(3,(mo-1)/up);
for (register int i=1;i<=up;i++) w[i]=cheng(w[i-1],rev)%mo;
for (register int i=0;i<up;i++)
{
G[i]=0;A[i]=0;
if (i<n) A[i]=b[i];
if (i<oq) G[i]=d[i];
}
FFT(A,up,1);
FFT(G,up,1);
for (register int i=0;i<up;i++)
{
A[i]=cheng(G[i],(2-cheng(A[i],G[i])+mo)%mo);
}
for (register int i=0;i<=up/2;i++)
{
swap(w[i],w[up-i]);
}
FFT(A,up,-1);
for (register int i=0;i<n;i++) d[i]=cheng(A[i],qsm(up,mo-2));
for (register int i=n;i<=up;i++) d[i]=0;
int kk=0;
}
}
inline void solve1(int n)
{
if (n==1)
{
b[0]=1;
}
else
{
double op=n;
int oq=ceil(op/2);
solve1(oq);
solve(n);
up=1;dep=0;
while (up<=n+jl) up=up*2,dep++;
w[0]=1;
rev=qsm(3,(mo-1)/up);
for (register int i=1;i<=up;i++) w[i]=cheng(w[i-1],rev);
for (register int i=0;i<up;i++)
{
A[i]=0;
if (i<n) A[i]=a[i];
}
FFT(A,up,1);FFT(b,up,1);FFT(d,up,1);
for (register int i=0;i<up;i++) b[i]=cheng(jia(b[i],cheng(A[i],d[i])),inv2);
for (register int i=0;i<=up/2;i++) swap(w[i],w[up-i]);
FFT(b,up,-1);
for (register int i=0;i<n;i++) b[i]=cheng(b[i],qsm(up,mo-2));
for (register int i=n;i<=up;i++) b[i]=0;
int kk=0;
}
}
int main()
{
scanf("%d",&n);
jl=n;
for (int i=0;i<n;i++)
{
scanf("%d",&a[i]);
}
solve1(n);
for (int i=0;i<n;i++)
{
printf("%d ",b[i]);
}
printf("
");
}