多项式求逆元模板【NTT/拆系数FFT】

最新推荐文章于 2023-12-16 21:46:55 发布
最新推荐文章于 2023-12-16 21:46:55 发布
阅读量434 收藏 3
点赞数
分类专栏： FFT 多项式
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/C20181220_xiang_m_y/article/details/85288785
版权

在这里插入图片描述

在这里插入图片描述
在这里插入图片描述
板题：

一. 洛谷4238：

模998244353，NTT：

#include<cstdio>
#include<algorithm>
#define maxn 400005
#define LL long long
using namespace std;
const int mod = 998244353, G = 3;
LL wn,w;
inline LL ksm(LL a,LL b){
	LL s=1;
	for(;b;b>>=1,a=a*a%mod) if(b&1) s=s*a%mod;
	return s;
}
inline void change(LL *a,int len)
{
	for(int i=1,j=len/2,k;i<len-1;i++)
	{
		if(i<j) swap(a[i],a[j]);
		for(k=len/2;j>=k;j-=k,k>>=1);
		j+=k;
	}
}
inline void ntt(LL *a,int len,int flg)
{
	change(a,len);
	for(int i=2;i<=len;i<<=1)
	{
		if(flg==1) wn=ksm(G,(mod-1)/i);
		else wn=ksm(G,mod-1-(mod-1)/i);
		for(int j=0;j<len;j+=i)
		{
			w=1;
			for(int k=j;k<j+i/2;k++)
			{
				LL u=a[k],v=w*a[k+i/2]%mod;
				a[k]=(u+v)%mod,a[k+i/2]=(u-v+mod)%mod;
				w=w*wn%mod;
			}
		}
	}
	if(flg==-1){
		LL ni=ksm(len,mod-2);
		for(int i=0;i<len;i++) a[i]=a[i]*ni%mod;
	}
}
int n,len;
LL a[maxn],b[maxn],tmp[maxn];
inline void poly_inv(int n,LL *a,LL *b)
{
	if(n==1) {b[0]=ksm(a[0],mod-2);return;}
	poly_inv((n+1)>>1,a,b);
	len=1;while(len<n+n-1) len<<=1;
	for(int i=0;i<len;i++) tmp[i]=i<n?a[i]:0;
	ntt(tmp,len,1),ntt(b,len,1);
	for(int i=0;i<len;i++) {b[i]=(2-tmp[i]*b[i]%mod)*b[i]%mod;if(b[i]<0) b[i]+=mod;}
	ntt(b,len,-1);
	for(int i=n;i<len;i++) b[i]=0;
}
int main()
{
	scanf("%d",&n);
	for(int i=0;i<n;i++) scanf("%lld",&a[i]);
	poly_inv(n,a,b);
	for(int i=0;i<n;i++) printf("%lld%c",b[i],i==n-1?10:32);
}
二. 洛谷4239

模109+7，拆系数FFT（7次fft版本）

// luogu-judger-enable-o2
#include<cstdio>
#include<cmath>
#include<algorithm>
#define maxn 400005
#define LL long long
using namespace std;
const int mod = 1e9+7, M1 = 31623, M2 = 14122;
const double Pi = acos(-1);
struct complex
{
	double r,i;
	complex(double _r=0,double _i=0):r(_r),i(_i){}
	complex operator + (const complex &t)const{return complex(r+t.r,i+t.i);}
	complex operator - (const complex &t)const{return complex(r-t.r,i-t.i);}
	complex operator * (const complex &t)const{return complex(r*t.r-i*t.i,r*t.i+i*t.r);}
	complex conj(){return complex(r,-i);}
}w[20][maxn/2];
inline void change(complex *a,int len)
{
	for(int i=1,j=len/2,k;i<len-1;i++)
	{
		if(i<j) swap(a[i],a[j]);
		for(k=len/2;j>=k;j-=k,k>>=1);
		j+=k;
	}
}
inline void fft(complex *a,int len,int flg)
{
	change(a,len);
	for(int i=2,o=0;i<=len;i<<=1,o++)
		for(int j=0;j<len;j+=i)
			for(int k=j;k<j+i/2;k++)
			{
				complex u=a[k],v=(flg==1?w[o][k-j]:w[o][k-j].conj())*a[k+i/2];
				a[k]=u+v,a[k+i/2]=u-v;
			}
	if(flg==-1) for(int i=0;i<len;i++) a[i].r/=len;
}
inline void MTT(int *a,int *b,int n,int *ret)
{
	static complex A[maxn],B[maxn],C[maxn],D[maxn],rt[3][maxn];
	int len=1;while(len<n+n-1) len<<=1;
	for(int i=0;i<len;i++)
		if(i<n)
		{
			A[i]=a[i]/M1,B[i]=a[i]%M1;
			C[i]=b[i]/M1,D[i]=b[i]%M1;
		}
		else A[i]=B[i]=C[i]=D[i]=0;
	fft(A,len,1),fft(B,len,1),fft(C,len,1),fft(D,len,1);
	for(int i=0;i<len;i++) rt[2][i]=A[i]*C[i],rt[1][i]=A[i]*D[i]+B[i]*C[i],rt[0][i]=B[i]*D[i];
	fft(rt[0],len,-1),fft(rt[1],len,-1),fft(rt[2],len,-1);
	for(int i=0;i<n;i++) ret[i]=((llround(rt[0][i].r)%mod+llround(rt[1][i].r)*M1%mod)%mod+llround(rt[2][i].r)*M2%mod)%mod;
}
inline int ksm(int a,int b){
	int s=1;
	for(;b;b>>=1,a=1ll*a*a%mod) if(b&1) s=1ll*s*a%mod;
	return s;
}
int n,a[maxn],b[maxn],tmp[maxn];
inline void poly_inv(int n,int *a,int *b)
{
	if(n==1) {b[0]=ksm(a[0],mod-2);return;}
	poly_inv((n+1)>>1,a,b);
	MTT(a,b,n,tmp);
	MTT(b,tmp,n,tmp);
	for(int i=0;i<n;i++) b[i]=(2*b[i]%mod-tmp[i]+mod)%mod;
}
int main()
{
	scanf("%d",&n);int len=1;while(len<n+n-1) len<<=1;
	for(int i=2,o=0;i<=len;i<<=1,o++)
		for(int j=0;j<i/2;j++)
			w[o][j]=complex(cos(2*Pi*j/i),sin(2*Pi*j/i));
	for(int i=0;i<n;i++) scanf("%d",&a[i]);
	poly_inv(n,a,b);
	for(int i=0;i<n;i++) printf("%d%c",b[i],i==n-1?10:32);
}

FFT（4次fft版本）
算法解析可以看这两篇 ↓ \darr ↓（膜拜大佬）
Miskcoo’s Space
Clin233

#include<cstdio>
#include<cmath>
#include<algorithm>
#define maxn 400005
#define LL long long
using namespace std;
const int mod = 1e9+7, M = (1<<15)-1;
const double Pi = acos(-1);
struct complex
{
	double r,i;
	complex(double _r=0,double _i=0):r(_r),i(_i){}
	complex operator + (const complex &t)const{return complex(r+t.r,i+t.i);}
	complex operator - (const complex &t)const{return complex(r-t.r,i-t.i);}
	complex operator * (const complex &t)const{return complex(r*t.r-i*t.i,r*t.i+i*t.r);}
	complex conj(){return complex(r,-i);}
}w[20][maxn/2];
inline void change(complex *a,int len)
{
	for(int i=1,j=len/2,k;i<len-1;i++)
	{
		if(i<j) swap(a[i],a[j]);
		for(k=len/2;j>=k;j-=k,k>>=1);
		j+=k;
	}
}
inline void fft(complex *a,int len,int flg)
{
	change(a,len);
	for(int i=2,o=0;i<=len;i<<=1,o++)
		for(int j=0;j<len;j+=i)
			for(int k=j;k<j+i/2;k++)
			{
				complex u=a[k],v=(flg==1?w[o][k-j]:w[o][k-j].conj())*a[k+i/2];
				a[k]=u+v,a[k+i/2]=u-v;
			}
	if(flg==-1) for(int i=0;i<len;i++) a[i].r/=len,a[i].i/=len;//注意这里的虚部也要除
}
inline void MTT(int *a,int *b,int n,int *ret)
{
	static complex A[maxn],B[maxn],rt[2][maxn],da,db,dc,dd;
	int len=1;while(len<n+n-1) len<<=1;
	for(int i=0;i<len;i++)
		if(i<n) A[i]=complex(a[i]>>15,a[i]&M),B[i]=complex(b[i]>>15,b[i]&M);
		else A[i]=B[i]=0;
	fft(A,len,1),fft(B,len,1);
	for(int i=0,j;i<len;i++)
	{
		j=(i==0?0:len-i);
		da=(A[i]+A[j].conj())*complex(0.5,0);
		db=(A[i]-A[j].conj())*complex(0,-0.5);
		dc=(B[i]+B[j].conj())*complex(0.5,0);
		dd=(B[i]-B[j].conj())*complex(0,-0.5);
		rt[0][i]=da*dc+db*dd*complex(0,1),rt[1][i]=db*dc+da*dd*complex(0,1);
	}
	fft(rt[0],len,-1),fft(rt[1],len,-1);
	for(int i=0;i<n;i++) ret[i]=((llround(rt[0][i].r)%mod<<30)+llround(rt[0][i].i)+((llround(rt[1][i].r)+llround(rt[1][i].i))%mod<<15))%mod;
}
inline int ksm(int a,int b){
	int s=1;
	for(;b;b>>=1,a=1ll*a*a%mod) if(b&1) s=1ll*s*a%mod;
	return s;
}
int n,a[maxn],b[maxn],tmp[maxn];
inline void poly_inv(int n,int *a,int *b)
{
	if(n==1) {b[0]=ksm(a[0],mod-2);return;}
	poly_inv((n+1)>>1,a,b);
	MTT(a,b,n,tmp);
	MTT(b,tmp,n,tmp);
	for(int i=0;i<n;i++) b[i]=(2*b[i]%mod-tmp[i]+mod)%mod;
}
int main()
{
	scanf("%d",&n);int len=1;while(len<n+n-1) len<<=1;
	for(int i=2,o=0;i<=len;i<<=1,o++)
		for(int j=0;j<i/2;j++)
			w[o][j]=complex(cos(2*Pi*j/i),sin(2*Pi*j/i));
	for(int i=0;i<n;i++) scanf("%d",&a[i]);
	poly_inv(n,a,b);
	for(int i=0;i<n;i++) printf("%d%c",b[i],i==n-1?10:32);
}

或者还有个优化版，change改成了用r[]，预处理单位根改成O(n)，空间小了很多

wlen=1;while(wlen<n+n-1) wlen<<=1;
for(int i=0;i<wlen;i++) w[i]=complex(cos(2*Pi*i/wlen),sin(2*Pi*i/wlen));

for(int i=0;i<len;i++) r[i]=(r[i>>1]>>1)|(i&1?len>>1:0);

int r[maxn],wlen;
inline void fft(complex *a,int len,int flg)
{
	for(int i=0;i<len;i++) if(i<r[i]) swap(a[i],a[r[i]]);
	for(int i=2;i<=len;i<<=1)
		for(int j=0,t=wlen/i;j<len;j+=i)//t=wlen/i
			for(int k=j,o=0;k<j+i/2;k++,o+=t)
			{
				complex u=a[k],v=(flg==1?w[o]:w[o].conj())*a[k+i/2];
				a[k]=u+v,a[k+i/2]=u-v;
			}
	if(flg==-1) for(int i=0;i<len;i++) a[i].r/=len,a[i].i/=len;//注意这里的虚部也要除
}

多项式求逆NTT最新版：

#include<bits/stdc++.h>
#define maxn 2200005
using namespace std;
const int mod = 998244353;
typedef vector<int> Poly;
inline int Pow(int a,int b){
	int s=1; for(;b;b>>=1,a=1ll*a*a%mod) if(b&1) s=1ll*s*a%mod;
	return s;
}
int w[maxn]={1},r[maxn];
inline int add(int a,int b){return (a+=b)>mod?a-mod:a;}
inline int dec(int a,int b){return (a-=b)<0?a+mod:a;}
void NTT(Poly &a,int len,int flg){
	a.resize(len);
	for(int i=0;i<len;i++) if(i<r[i]) swap(a[i],a[r[i]]);
	for(int i=2,l=1,G=flg==1?3:(mod+1)/3;i<=len;i<<=1,l<<=1){
		for(int k=l-2,wn=Pow(G,(mod-1)/i);k>=0;k-=2) w[k+1]=1ll*(w[k]=w[k>>1])*wn%mod;
		for(int j=0;j<len;j+=i) for(int k=j;k<j+l;k++){
			int u=a[k],v=1ll*w[k-j]*a[k+l]%mod;
			a[k]=add(u,v),a[k+l]=dec(u,v);
		}
	}
	if(flg==-1) for(int i=0,Inv=Pow(len,mod-2);i<len;i++) a[i]=1ll*a[i]*Inv%mod;
}
Poly Inv(Poly A,int n){
	Poly B(1,Pow(A[0],mod-2)),C;
	for(int k=2,len=4;k<n<<1;k=len,len<<=1){
		C=Poly(A.begin(),A.begin()+min(n,k));
		for(int i=0;i<len;i++) r[i]=r[i>>1]>>1|(i&1?len>>1:0);
		NTT(C,len,1),NTT(B,len,1);
		for(int i=0;i<len;i++) B[i]=1ll*B[i]*(2-1ll*C[i]*B[i]%mod+mod)%mod;
		NTT(B,len,-1),B.resize(min(n,k));
	}
	return B;
}
int n; Poly A;
int main()
{
	scanf("%d",&n),A.resize(n);
	for(int i=0;i<n;i++) scanf("%d",&A[i]);
	Poly B(Inv(A,n));
	for(int i=0;i<n;i++) printf("%d%c",B[i],i==n-1?10:32);
}
Master.Yi
关注
专栏目录
android多项式求乘法逆元小软件源码
09-16
这个“android多项式求乘法逆元小软件源码”提供了一个实用的解决方案，尤其对于需要处理数学问题的开发者来说，具有很高的价值。下面我们将详细探讨相关知识点。 1. **多项式**： 多项式是数学中的基本概念，由...
第三次训练赛题解
fdxgcw的博客
08-22 832
题目背景 998244353最优美的性质莫过于它是个完美恶臭数： 998244353=( 114514 * ( 54-1+114 * (1+145+1+4) ) ) + ( 4+11451 * (4-1-15+14) ) + ( 11+4154 + (141+541) ) + (4-1-15+14) 之所以称为完美，如果您仔细观察，您会发现每一个括号里114514都出现了正整数次！ 今天，你需要用这个恶臭数进行一系列的复杂运算(或许?) 题目描述 米4达的替身名为"性感手枪"，可以操纵子弹的飞行轨迹。一天在
斐波那契数列
罗马尼亚硬拉
08-26 1172
一、递归 已知斐波那契数列Fn​=Fn−1​+Fn−2​(n>=3),F1​=1,F2​=1用递归的方法求解该数列的第n项。 输入格式: 输入一个正整数n(1<=n<=40)。 输出格式: 输出一个数，数列的第n项 输入样例1: 1 结尾无空行 输出样例1: 1 结尾无空行 输入样例2: 3 结尾无空行 输出样例2: 2 结尾无空行 答案： #include <bits/stdc++.h> using namespac...
多项式求逆元
zmnbehappy的博客
11-07 7318
Contents [hide] 1 概述2 基本概念3 多项式的逆元 3.1 求解过程3.2 代码实现 4  应用 4.1 预处理 Bernoulli 数4.2 计算有标号简单连通无向图个数 概述 多项式求逆元是多项式除法和多项式取模的必要过程，在竞赛中，有时候多项式求逆元的应用比多项式的除法还要多，用快速傅里叶变换以及倍增算法可以做到在 O(nlogn
拆系数FFT（任意模数NTT）
a_forever_dream的博客
05-26 579
拆系数FFT 任意模数的话，就不能直接用 NTTNTTNTT 了，有一种做法是取三个不同模数做 NTTNTTNTT 然后用 crtcrtcrt 合在一起，但是这种做法要用 999 次 NTTNTTNTT，会比较慢，而另外一种拆系数 FFTFFTFFT 的做法就会快很多。 设 M=m0M=\sqrt {m_0}M=m0​​，m0m_0m0​ 是系数的最大值，对于每个系数，将他拆成 aM+b (b<M)aM+b~(b<M)aM+b (b<M) 的形式，然后以 aaa 和 b
多项式算法5：多项式求逆
ha_ing的博客
06-03 2611
多项式算法5：多项式求逆
jisuanqi.rar_求多项式系数
09-23
求两个一元多项式A（X）=a0+a1x+a2x^2+....+anx^n和B（X）=b0+b1x+b2x^2+....+bnx^n 的和C（X），差D（X），积（X）。 算法的输入：两个多项式（分别输入各项的系数和指数） 算法的输出：多项式的和,差，积。
扩展欧几里得算法的matlab程序（求多项式的乘法逆元）
06-23
函数的输入可能是多项式的系数向量或字符串形式，而输出可能是表示乘法逆元的多项式或其系数向量。 在实际应用中，理解并正确使用这个MATLAB函数能够帮助你解决涉及多项式逆元的问题，例如在计算模逆元的线性同余...
FFT.rar_FFT多项式乘积_fft java_java fft
09-21
标题中的"FFT.rar_FFT多项式乘积_fft java_java fft"揭示了我们即将探讨的核心主题：快速傅里叶变换（FFT）在计算多项式乘法中的应用，以及它与Java编程语言的结合。FFT是一种高效的算法，用于计算离散傅里叶变换...
求逆元的C程序（算法）
01-28
用c语言实现逆元的计算，通过自己设计算法代码实现。
校OJ：树上倍增
魔法猫咪的博客
03-14 811
树上倍增 问题描述   给定一棵树，树的大小为nnn且根节点为节点1，求有 多少组(i,j,k)(i,j,k)(i,j,k)使得节点kkk是节点jjj的祖先且节点jjj是节点iii的祖先。由于这个数可能很大，所以请输出对998244353取模后的结果。 输入 第一行包含一个正整数n(1≤n≤106)n(1\le n \le 10^6)n(1≤n≤106)，表示树的节点个数 接下来n−1n-1n−1行每行包括两个整数u,vu,vu,v，表示节点uuu和节点vvv相连 输出 一行一个整数表示答案对9982443
多项式求逆
a_forever_dream的博客
10-10 8598
介绍 准确来说，是求多项式的逆元。 对于多项式 A(x)A(x)A(x)，如果存在 A(x)B(x)≡1(modxn)A(x)B(x)\equiv 1 \pmod {x^n}A(x)B(x)≡1(modxn)，那么称 B(x)B(x)B(x) 为 A(x)A(x)A(x) 在模 xnx^nxn 意义下的逆元。 注意，这里的模 xnx^nxn，是指舍弃含有 xnx^nxn 及更高次的项。 比如说有一...
逆元模板
weixin_30547797的博客
02-28 115
逆元是求ax≡1(mod p)中的x的最小正整数，常应用于 除法取模。 存在这样的x的条件是a与p互素，即只有当a与p互素时，逆元存在 求逆元有三种常用的方式 扩展欧几里德（要求a与p互素） int extgcd(int a, int b, int &x, int &y) { int d = a; if(b != 0) ...
逆元 笔记 （2）
最新发布
cout_x的博客
12-16 864
逆元 笔记 2
任意模数NTT（拆系数FFT）
tanjunming2020的博客
07-10 287
### 介绍 一般的NTT的模数$P$都要满足$P=r\times 2^k+1$。但如果不满足这个条件，那就不能直接用NTT了，需要用到任意模数NTT。 下面就来介绍一种任意模数NTT——拆系数FFT。
拆系数FFT
dieyi9889的博客
03-10 230
拆系数FFT 表示才发现自己没有掌握这个似乎烂大街了的科技了…… 概念： 应对那种模数比较恶心人的多项式乘法，大概就是吧一个多项式拆成两个，然后让乘法不会爆掉，最后再进行取模。既然拆成了两个多项式，\(DFT\)和\(IDFT\)次数自然就会变多，一共有\(7\)次的和\(4\)次的两种写法，自然是后面的快一些啦，但是后边的精度要求比较高，并且一般也不会卡的这么严重的，这里就只介绍第一种...
C++题解：Darko 的小球
Keven_11的博客
08-12 918
目录 题目 题解 题目 Darko 是一个喜欢数数的女孩子，她喜欢唱，跳，rap 和玩球。 Darko 有n个黑球和m个白球，她想把这些球放进n个不同的盒子。 求满足每一个盒子内黑白球不同时出现的方案数，两个方案不同当且仅当存在对应盒子内黑白球个数不同。 n,m≤106。答案对 998244353取模。 输入格式 输入共11行22个整数n,m。 输出格式 输出共11行一个整数表示方案数对 998244353取模的结果。 数据规模与约定...
PTA 2 斐波那契数列（II）
我的博客
09-20 1000
已知斐波那契数列Fn​=Fn−1​+Fn−2​(n>=3),F1​=1,F2​=1 求解该数列的第n项，结果对998244353取模。 输入格式: 输入一个正整数n(1<=n<=10000000)。 输出格式: 输出一个数，数列的第n项 输入样例1: 1 结尾无空行 输出样例1: 1 结尾无空行 输入样例2: 3 结尾无空行 输出样例2: 2 结尾无空行 #include<stdio.h> #define N_max 99824..
PTA 斐波那契数列 求模
阿菠萝叻的博客
03-06 2748
已知斐波那契数列 F​n​​=F​n−1​​+F​n−2​​(n>=3),F​1​​=1,F​2​​=1 求解该数列的第n项，结果对998244353取模。 输入格式: 输入一个正整数n (1<=n<=10000000)。 输出格式: 输出一个数，数列的第n项 我的原答案如下，最后一个测试点过不了： #include <iostream> using namespace std; int main() { int n; cin>>n; int F(int );
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值