多项式求逆

多项式求逆

定义

\(\displaystyle f(x) =\sum^{n-1}_{k=0}a_kx^k\)\(g(x) =\sum^{n-1}_{k=0}b_kx^k\),使得
\(\displaystyle f(x)g(x)\equiv 1 (\mod x^n)\)

\(\displaystyle f(x)g(x)\) 的前\(n\)项中只有常数项为\(1\),其余项均为\(0\),称为

\(f(x)\)\(\displaystyle x^n\)的逆元,或简称求\(f(x)\)的逆元。

多项式可逆当且仅当其常数项可逆,证明见“方法一”\(\displaystyle b_0=\frac{1}{a_0}\)

方法1:分治FFT

从常数项出发,比较系数得\(\displaystyle b_0=\frac{1}{a_0}\)\(\displaystyle \sum^k_{i=0}b_ia_{k-i}=0\)那么有

\[b_k=\sum^{k-1}_{i=0}b_i(-\frac{a_{k-i}}{a_0}) \]

从而可以通过分治FFT 计算,时间复杂度\(\displaystyle \Theta(n \log^2 n)\)

方法2:倍增法

考虑倍增,从局部的逆元(低次)出发,设\(\displaystyle g_0(x) = g(x) \mod x^k\),那么有

\[f(x)g_0(x)-1\equiv 0 (\mod x^k) \]

从而有

\[f(x)g_0(x)-1)^2\equiv 0 (mod x^{2k}) \]

可以发现模的次数增加了,此时展开左边可得

\[g(x)\equiv g_0(x)(2-f(x)g_0(x))(\mod x^{2k}) \]

从而可以算出\(g(x)\)的前\(2k\)项。

通过摊还分析得时间复杂度为\(\displaystyle \Theta(n \log n)\)

实现

多项式求逆的实现不难,两种方法都可以,但第二种更快,也更容易实现,故此处以方法二为例。
我们每次递归算出模\(\displaystyle x^{\frac{n}{2}}\)的逆元,记为\(g_0\),再套用\(\displaystyle g(x)\equiv g_0(x)(2-f(x)g_0(x))(\mod x^{2k})\),用FFT计算即可

void solve(LL *a,LL *b,int p){/*a:seq,b:inv*/
	if(p==1){
		b[0]=fpm(a[0],MOD-2);
		return;
	}
	solve(a,b,(p+1)>>1);
	lim=1; L=0;
	while(lim<(p<<1)){
		lim<<=1;
		L++;
	}
	for(int i=0;i<lim;i++)
		rev[i]=rev[i>>1]>>1|(1&i)<<(L-1);
	for(int i=0;i<p;i++)
		tmp[i]=a[i];
	for(int i=p;i<lim;i++)
		tmp[i]=0;
	NTT(tmp,1); NTT(b,1);
	for(int i=0;i<lim;i++)
		b[i]=(2-b[i]*tmp[i]%MOD+MOD)%MOD*b[i]%MOD;
	NTT(b,-1);
	for(int i=p;i<lim;i++)
		b[i]=0;
}

例题

[P4238 【模板】多项式求逆](%3Ca href="https://www.luogu.org/problem/P4238"%3Ehttps://www.luogu.org/problem/P4238%3C/a%3E)

这是纯粹的模板题,直接给代码

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const int INF=1e9+7,MAXN=3e6+10/*Min:2^20+10*/;
const LL G=3,MOD=998244353,INV=332748118;
inline LL fpm(LL base,LL p){
	LL ret=1;
	while(p){
		if(p&1)
			ret=ret*base%MOD;
		base=base*base%MOD;
		p>>=1;
	}
	return ret%MOD;
}
int N,M,lim=1,L,rev[MAXN];
inline void NTT(LL *a,int type){
	for(int i=0;i<lim;i++)
		if(i<rev[i])
			swap(a[i],a[rev[i]]);
	for(int mid=1;mid<lim;mid<<=1){
		int len=mid<<1/*n*/;
		LL Wn=fpm(G,(MOD-1)/len);
		for(int j=0;j<lim;j+=len){
			LL w=1,t1,t2;
			for(int k=0;k<mid;k++){
				t1=a[j+k],t2=w*a[j+k+mid]%MOD;
				a[j+k]=(t1+t2)%MOD;
				a[j+k+mid]=(t1-t2+MOD)%MOD;
				w=w*Wn%MOD;
			}
		}
	}
	if(type==-1){
		LL lim_inv=fpm(lim,MOD-2);
		reverse(a+1,a+lim);
		for(int i=0;i<lim;i++)
			a[i]=a[i]*lim_inv%MOD;
	}
}
LL tmp[MAXN];
void solve(LL *a,LL *b,int p){/*a:seq,b:inv*/
	if(p==1){
		b[0]=fpm(a[0],MOD-2);
		return;
	}
	solve(a,b,(p+1)>>1);
	lim=1; L=0;
	while(lim<(p<<1)){
		lim<<=1;
		L++;
	}
	for(int i=0;i<lim;i++)
		rev[i]=rev[i>>1]>>1|(1&i)<<(L-1);
	for(int i=0;i<p;i++)
		tmp[i]=a[i];
	for(int i=p;i<lim;i++)
		tmp[i]=0;
	NTT(tmp,1); NTT(b,1);
	for(int i=0;i<lim;i++)
		b[i]=(2-b[i]*tmp[i]%MOD+MOD)%MOD*b[i]%MOD;
	NTT(b,-1);
	for(int i=p;i<lim;i++)
		b[i]=0;
}
LL a[MAXN],b[MAXN];
int main(){
	scanf("%d",&N);
	for(int i=0;i<N;i++){
		scanf("%lld",a+i);
		a[i]=(a[i]+MOD)%MOD;
	}
	solve(a,b,N);
	for(int i=0;i<N;i++)
		printf("%lld ",b[i]);
	return 0;
}
[P4239 【模板】多项式求逆(加强版)](%3Ca href="https://www.luogu.org/problem/P4239"%3Ehttps://www.luogu.org/problem/P4239%3C/a%3E)

这道题是在任意模数下的多项式求逆,需要将上一个板子的NTT改成MTT

posted @ 2019-08-03 19:43  guoshaoyang  阅读(4150)  评论(0编辑  收藏  举报