时间戳

首页 > 技术文章 > 多项式总结（unfinished）

多项式总结（unfinished）

ErkkiErkko 2018-11-28 22:07 原文

试试以二级标题为主的格式。

多项式相关

注：本篇博客不包含\(FFT\)基础姿势。如果您想要阅读本篇博客，请确保自己对\(FFT,NTT\)有基本的认识并且能够独立写出代码。

多项式是什么？

左转数学七年级上册课本。

多项式的两种表示法

系数表示法和点值表示法，过于基础不多解释。

多项式的四则运算及扩展

多项式加减法

同类项直接进行加减即可，过于简单不多解释。

多项式乘法

朴素算法：\(O(n^2)\)。

分治乘法：没写过，不常用，不知道。

快速傅里叶变换（\(FFT\)）：\(O(nlogn)\)。

快速数论变换（\(NTT\)）：\(O(nlogn)\)。

（一般认为\(NTT\)的常数比\(FFT\)大，两者的适用情况略有不同。）

多项式求逆

模板题：[洛谷P4238][模板]多项式求逆

给定一个多项式\(F(x)\)，请求出一个多项式\(G(x)\)，满足\(F(x) \times G(x) \equiv 1\ (mod\ x^n)\)。系数对\(998244353\)取模。

推倒过程

显然\(n=1\)时，直接算乘法逆元即可。

假设我们已经求出了\(H(x)\)满足\(F(x) \times H(x) \equiv 1\ (mod\ x^{\lceil \frac{n}{2} \rceil})\)，我们现在想求出\(G(x)\)满足\(F(x) \times G(x) \equiv 1\ (mod\ x^n)\)。

显然有\(F(x) \times G(x) \equiv 1\ (mod\ x^{\lceil \frac{n}{2} \rceil})\)成立，所以我们可以得到\(G(x)-H(x) \equiv 0\ (mod\ x^{\lceil \frac{n}{2} \rceil})\)。

同余号两边同时平方，并且我们发现模数也可以平方，得：

\[(G(x)-H(x))^2 \equiv 0\ (mod\ x^n) \]

\[G(x)^2-2 \times G(x) \times H(x)+H(x)^2 \equiv 0\ (mod\ x^n) \]

同余号两边同时乘\(F(x)\)，得：

\[G(x)-2 \times H(x)+H(x)^2 \times F(x) \equiv 0\ (mod\ x^n) \]

\[G(x) \equiv 2 \times H(x)-H(x)^2 \times F(x)\ (mod\ x^n) \]

从上面的推倒过程我们可以发现，想要求出\(G(x)\)，要先求出\(H(x)\)。求出\(H(x)\)相当于是这个问题的一个规模缩小的子问题，于是我们可以考虑递归求解。

时间复杂度：\(O(nlogn)\)。

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cctype>
#include <algorithm>
#define rin(i,a,b) for(int i=(a);i<=(b);i++)
#define rec(i,a,b) for(int i=(a);i>=(b);i--)
#define trav(i,a) for(int i=head[(a)];i;i=e[i].nxt)
using std::cin;
using std::cout;
using std::endl;
typedef long long LL;

inline int read(){
	int x=0;char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x;
}

const int MAXN=100005;
const int MAXLEN=270005;
const LL MOD=998244353;
const LL G=3;
const LL INVG=332748118;
int n,m;
int a[MAXLEN];
int len,rev[MAXLEN];
LL invn,A[MAXLEN],B[MAXLEN];

inline LL qpow(LL x,LL y){
	LL ret=1,tt=x%MOD;
	while(y){
		if(y&1) ret=ret*tt%MOD;
		tt=tt*tt%MOD;
		y>>=1;
	}
	return ret;
}

inline void ntt(LL *c,int dft){
	rin(i,0,n-1) if(i<rev[i])
		std::swap(c[i],c[rev[i]]);
	for(int mid=1;mid<n;mid<<=1){
		int r=(mid<<1);
		LL u=qpow(dft>0?G:INVG,(MOD-1)/r);
		for(int l=0;l<n;l+=r){
			LL v=1;
			for(int i=0;i<mid;i++,v=v*u%MOD){
				LL x=c[l+i],y=c[l+mid+i]*v%MOD;
				c[l+i]=x+y;
				if(c[l+i]>=MOD) c[l+i]-=MOD;
				c[l+mid+i]=x-y;
				if(c[l+mid+i]<0) c[l+mid+i]+=MOD;
			}
		}
	}
	if(dft<0) rin(i,0,n-1)
		c[i]=c[i]*invn%MOD;
}		

void solve(int Floor){
	if(Floor==1){
		B[0]=qpow(a[0],MOD-2);
		return;
	}
	solve((Floor+1)>>1);
	m=Floor-1+((((Floor+1)>>1)-1)<<1);
	len=0;
	for(n=1;n<=m;n<<=1) len++;
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
	rin(i,0,n-1) A[i]=(i<Floor?a[i]:0);
	invn=qpow(n,MOD-2);
	ntt(A,1);
	ntt(B,1);
	rin(i,0,n-1) B[i]=((2*B[i]-A[i]*B[i]%MOD*B[i])%MOD+MOD)%MOD;
	ntt(B,-1);
	rin(i,Floor,n-1) B[i]=0;
}

int main(){
	n=read();
	n--;
	rin(i,0,n){
		a[i]=read();
	}
	int nn=n;
	solve(n+1);
	rin(i,0,nn){
		printf("%lld ",B[i]);
	}
	printf("\n");
	return 0;
}

多项式除法/多项式取模

模板题：[洛谷P4512][模板]多项式除法

给定一个\(n\)次多项式\(F(x)\)和一个\(m\)次多项式\(G(x)\)，请求出多项式\(Q(x),R(x)\)，满足以下条件：

  1. \(Q(x)\)次数为\(n−m\)，\(R(x)\)次数小于\(m\)。

  2. \(F(x)=Q(x) \times G(x)+R(x)\)。

所有的运算在模\(998244353\)意义下进行。

推倒过程

首先我们要知道对于任意一个\(n\)次多项式\(f(x)\)，将其系数翻转后得到的多项式为\(f_R(x)=x^nf(\frac{1}{x})\)。

回到题目，将等式中的所有多项式的自变量改为\(\frac{1}{x}\)，等号左右同时乘\(x\)，得：

\[x^nF(\frac{1}{x})=x^{n-m}Q(\frac{1}{x}) \times x^mG(\frac{1}{x})+x^{n-m+1} \times x^{m-1}R(\frac{1}{x}) \]

\[F_R(x)=Q_R(x) \times G_R(x)+x^{n-m+1} \times R_R(x) \]

\[F_R(x) \equiv Q_R(x) \times G_R(x)\ (mod\ x^{n-m+1}) \]

\[Q_R(x) \equiv F_R(x) \times G_R^{-1}(x)\ (mod\ x^{n-m+1}) \]

利用多项式求逆可以求出\(Q_R(x)\)，然后根据：

\[R(x)=F(x)-Q(x) \times G(x) \]

可以求出\(R(x)\)。

时间复杂度：\(O(nlogn)\)。

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cctype>
#include <algorithm>
#define rin(i,a,b) for(int i=(a);i<=(b);i++)
#define rec(i,a,b) for(int i=(a);i>=(b);i--)
#define trav(i,a) for(int i=head[(a)];i;i=e[i].nxt)
using std::cin;
using std::cout;
using std::endl;
typedef long long LL;

inline int read(){
	int x=0;char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x;
}

const int MAXN=100005;
const int MAXLEN=270005;
const LL MOD=998244353;
const LL ROOT=3;
const LL INVROOT=332748118;
int n,m;
int len,rev[MAXLEN];
int F[MAXN],G[MAXN];
LL Q[MAXN],R[MAXN];
LL invn;
LL A[MAXLEN],B[MAXLEN];

inline LL qpow(LL x,LL y){
	LL ret=1,tt=x%MOD;
	while(y){
		if(y&1) ret=ret*tt%MOD;
		tt=tt*tt%MOD;
		y>>=1;
	}
	return ret;
}

inline void ntt(LL *c,int dft){
	rin(i,0,n-1) if(i<rev[i])
		std::swap(c[i],c[rev[i]]);
	for(int mid=1;mid<n;mid<<=1){
		int r=(mid<<1);
		LL u=qpow(dft>0?ROOT:INVROOT,(MOD-1)/r);
		for(int l=0;l<n;l+=r){
			LL v=1;
			for(int i=0;i<mid;i++,v=v*u%MOD){
				LL x=c[l+i],y=c[l+mid+i]*v%MOD;
				c[l+i]=x+y;
				if(c[l+i]>=MOD) c[l+i]-=MOD;
				c[l+mid+i]=x-y;
				if(c[l+mid+i]<0) c[l+mid+i]+=MOD;
			}
		}
	}
	if(dft<0) rin(i,0,n-1)
		c[i]=c[i]*invn%MOD;
}

void getinv(int Floor){
	if(Floor==1){
		B[0]=qpow(G[0],MOD-2);
		return;
	}
	getinv((Floor+1)>>1);
	m=Floor-1+((((Floor+1)>>1)-1)<<1);
	len=0;
	for(n=1;n<=m;n<<=1) len++;
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
	rin(i,0,n-1) A[i]=(i<Floor?G[i]:0);
	invn=qpow(n,MOD-2);
	ntt(A,1);
	ntt(B,1);
	rin(i,0,n-1) B[i]=((2*B[i]-A[i]*B[i]%MOD*B[i])%MOD+MOD)%MOD;
	ntt(B,-1);
	rin(i,Floor,n-1) B[i]=0;
}

int main(){
	n=read(),m=read();
	rin(i,0,n) F[i]=read();
	rin(i,0,m) G[i]=read();
	std::reverse(F,F+n+1);
	std::reverse(G,G+m+1);
	int nn=n,mm=m;
	getinv(n-m+1);
//求Q
	n=nn,m=mm;
	memset(A,0,sizeof A);
	rin(i,0,n) A[i]=F[i];
	m=n+n-m;
	len=0;
	for(n=1;n<=m;n<<=1) len++;
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
	invn=qpow(n,MOD-2);
	ntt(A,1);
	ntt(B,1);
	rin(i,0,n-1) A[i]=A[i]*B[i]%MOD;
	ntt(A,-1);
	n=nn,m=mm;
	rin(i,0,n-m) Q[i]=A[i];
//求R
	std::reverse(F,F+n+1);
	std::reverse(G,G+m+1);
	std::reverse(Q,Q+n-m+1);
	memset(A,0,sizeof A);
	memset(B,0,sizeof B);
	rin(i,0,n-m) A[i]=Q[i];
	rin(i,0,m) B[i]=G[i];
	m=n;
	len=0;
	for(n=1;n<=m;n<<=1) len++;
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
	invn=qpow(n,MOD-2);
	ntt(A,1);
	ntt(B,1);
	rin(i,0,n-1) A[i]=A[i]*B[i]%MOD;
	ntt(A,-1);
	rin(i,0,m-1) R[i]=(F[i]-A[i]+MOD)%MOD;
	n=nn,m=mm;
	rin(i,0,n-m) printf("%lld ",Q[i]);
	printf("\n");
	rin(i,0,m-1) printf("%lld ",R[i]);
	printf("\n");
	return 0;
}

多项式开方

关于这一部分博主还不能很透彻地理解，所以先放出简要推倒过程。在给出多项式的常数项不为\(0\)或\(1\)时，求开方后多项式的常数项时需要用到二次剩余。但是由于博主目前仍未学习二次剩余（太菜了），所以只能给出给出多项式的常数项为\(0\)或\(1\)时的代码，不久之后会完善这一部分。

对于一个多项式\(A(x)\)，如果存在一个多项式\(B(x)\)，满足\(B(x)^2 \equiv A(x)\ (mod\ x^n)\)，则称\(B(x)\)为\(A(x)\)在\(mod\ x^n\)下的平方根。

所有的运算在模\(998244353\)意义下进行。

推倒过程

\(B(x)\)的常数项可以由\(A(x)\)的常数项二次剩余求出。特别的，当\(A(x)\)的常数项为\(0\)或\(1\)时，\(B(x)\)的常数项与\(A(x)\)相等。

类似多项式求逆，假设我们已经求出了\(C(x)^2 \equiv A(x)\ (mod\ x^{\lceil \frac{n}{2} \rceil})\)，现在我们想要求\(B(x)^2 \equiv A(x)\ (mod\ x^n)\)。

显然\(B(x)^2 \equiv A(x)\ (mod\ x^{\lceil \frac{n}{2} \rceil})\)。

所以有\(B(x)^2-C(x)^2 \equiv 0\ (mod\ x^{\lceil \frac{n}{2} \rceil})\)。

两边同时平方，得：

\[(B(x)^2-C(x)^2)^2 \equiv 0\ (mod\ x^n) \]

\[(B(x)^2+C(x)^2)^2 \equiv 4B(x)^2C(x)^2\ (mod\ x^n) \]

\[B(x)^2+C(x)^2 \equiv 2B(x)C(x)\ (mod\ x^n) \]

\[A(x)+C(x)^2 \equiv 2B(x)C(x)\ (mod\ x^n) \]

\[B(x) \equiv (A(x)+C(x)^2) \times (2C(x))^{-1}\ (mod\ x^n) \]

然后多项式求逆即可。

时间复杂度：\(O(nlogn)\)。

Update on 2018/12/13：一段更靠谱的证明：

\[(B(x)^2-C(x)^2)^2 \equiv 0\ (mod\ x^n) \]

\[(B(x)^2+C(x)^2)^2 \equiv 4B(x)^2C(x)^2\ (mod\ x^n) \]

\[\frac{(B(x)^2+C(x)^2)^2}{4 \times C(x)^2} \equiv B(x)^2\ (mod\ x^n) \]

\[B(x) \equiv \frac{B(x)^2+C(x)^2}{2 \times C(x)}\ (mod\ x^n) \]

\[B(x) \equiv \frac{A(x)+C(x)^2}{2 \times C(x)}\ (mod\ x^n) \]

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cctype>
#include <algorithm>
#define rin(i,a,b) for(int i=(a);i<=(b);i++)
#define rec(i,a,b) for(int i=(a);i>=(b);i--)
#define trav(i,a) for(int i=head[(a)];i;i=e[i].nxt)
using std::cin;
using std::cout;
using std::endl;
typedef long long LL;

inline int read(){
	int x=0;char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x;
} 

const int MAXN=100005;
const int MAXLEN=270005;
const LL MOD=998244353;
const LL ROOT=3;
const LL INVROOT=332748118;
const LL INV2=499122177;
int n,m;
int len,rev[MAXLEN];
int a[MAXN];
LL invn,A[MAXLEN],B[MAXLEN],C[MAXLEN];

inline LL qpow(LL x,LL y){
	LL ret=1,tt=x%MOD;
	while(y){
		if(y&1) ret=ret*tt%MOD;
		tt=tt*tt%MOD;
		y>>=1;
	}
	return ret;
}

inline void ntt(LL *c,int dft){
	rin(i,0,n-1) if(i<rev[i])
		std::swap(c[i],c[rev[i]]);
	for(int mid=1;mid<n;mid<<=1){
		int r=(mid<<1);
		LL u=qpow(dft>0?ROOT:INVROOT,(MOD-1)/r);
		for(int l=0;l<n;l+=r){
			LL v=1;
			for(int i=0;i<mid;i++,v=v*u%MOD){
				LL x=c[l+i],y=c[l+mid+i]*v%MOD;
				c[l+i]=x+y;
				if(c[l+i]>=MOD) c[l+i]-=MOD;
				c[l+mid+i]=x-y;
				if(c[l+mid+i]<0) c[l+mid+i]+=MOD;
			}
		}
	}
	if(dft<0) rin(i,0,n-1)
		c[i]=c[i]*invn%MOD;
}

void getinv(int Floor){
	if(Floor==1){
		C[0]=qpow(B[0],MOD-2);
		return;
	}
	getinv((Floor+1)>>1);
	m=Floor-1+((((Floor+1)>>1)-1)<<1);
	len=0;
	for(n=1;n<=m;n<<=1) len++;
	invn=qpow(n,MOD-2);
	rin(i,0,n-1) A[i]=(i<Floor?B[i]:0);
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
	ntt(A,1);
	ntt(C,1);
	rin(i,0,n-1) C[i]=(((2*C[i]-A[i]*C[i]%MOD*C[i])%MOD)+MOD)%MOD;
	ntt(C,-1);
	rin(i,Floor,n-1) C[i]=0;
}

void solve(int Floor){
	if(Floor==1){
		B[0]=a[0];
		return;
	}
	solve((Floor+1)>>1);
	getinv(Floor);
	m=Floor-1+((((Floor+1)>>1)-1)<<1);
	len=0;
	for(n=1;n<=m;n<<=1) len++;
	invn=qpow(n,MOD-2);
	rin(i,0,n-1) A[i]=(i<Floor?a[i]:0);
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
//	ntt(A,1);
//	ntt(B,1);
//	ntt(C,1);
//	rin(i,0,n-1) B[i]=(A[i]+B[i]*B[i])%MOD*C[i]%MOD*INV2%MOD;
//	ntt(B,-1);
	ntt(A,1);
	ntt(C,1);
	rin(i,0,n-1) A[i]=A[i]*C[i]%MOD;
	ntt(A,-1);
	rin(i,0,Floor-1) B[i]=(A[i]+B[i])*INV2%MOD;
	rin(i,Floor,n-1) B[i]=0;
	rin(i,0,n-1) C[i]=0;
}

int main(){
	n=read();
	n--;
	rin(i,0,n) a[i]=read();
	int nn=n;
	solve(n+1);
	rin(i,0,nn) printf("%lld ",B[i]);
	printf("\n");
	return 0;
}

拉格朗日插值

模板题：[洛谷P4781][模板]拉格朗日插值

由初中知识可知，\(n\)个点\((x_i,y_i)\)可以唯一地确定一个多项式。

现在，给定\(n\)个点，请你确定这个多项式，并将\(k\)代入求值，无需输出多项式每项的系数。

求出的值对\(998244353\)取模。

推倒过程

令\(f(x)=\sum_{i=1}^{n}y_il_i(x)\)，\(l_i(x)=\prod_{j=1,j \neq i}^{n}\frac{x-x_j}{x_i-x_j}\)。

容易发现\(l_i(x_i)=1\)，\(l_i(x_j)=0\ (j \neq i)\)。

显然\(f(x)\)过所有\(n\)个点，\(f(x)\)为所求多项式，求出\(f(k)\)即可。

时间复杂度：\(O(n^2)\)。

Update on 2019/2/14：如果点值连续或等差，时间复杂度可以优化到\(O(n)\)。

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cctype>
#include <algorithm>
#define rin(i,a,b) for(int i=(a);i<=(b);i++)
#define rec(i,a,b) for(int i=(a);i>=(b);i--)
#define trav(i,a) for(int i=head[(a)];i;i=e[i].nxt)
using std::cin;
using std::cout;
using std::endl;
typedef long long LL;

inline int read(){
	int x=0;char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x;
}

const int MAXN=100005;
const LL MOD=998244353;
int n;
LL k,X[MAXN],Y[MAXN];

inline LL qpow(LL x,LL y){
	LL ret=1,tt=x%MOD;
	while(y){
		if(y&1) ret=ret*tt%MOD;
		tt=tt*tt%MOD;
		y>>=1;
	}
	return ret;
}

int main(){
	n=read(),k=read();
	rin(i,1,n){
		X[i]=read();
		Y[i]=read();
	}
	LL ans=0;
	rin(i,1,n){
		LL son=1,mot=1;
		rin(j,1,n){
			if(i==j) continue;
			son=son*(k-X[j]+MOD)%MOD;
			mot=mot*(X[i]-X[j]+MOD)%MOD;
		}
		ans=(ans+Y[i]*son%MOD*qpow(mot,MOD-2))%MOD;
	}
	printf("%lld\n",ans);
	return 0;
}

所有给出点的横坐标是连续的话，还可以做到\(O(n)\)求解。但是博主太菜了，你们懂的。

任意模数\(NTT\)

三模数\(NTT\)

由于两个多项式卷积后每项系数的大小是\(O(na_i^2)\)级别的，所以我们只需选取两到三个大小合适的可以进行\(NTT\)的素数分别\(NTT\)后\(CRT\)合并答案即可，可能需要用到快速乘防止爆\(long\ long\)。

需要做\(9\)遍\(NTT\)，常数巨大。

拆系数\(FFT\)

我们可以把给出的多项式系数\(a_i\)写成\(a_i=p_i \times 32768+q_i\)。同理\(b_i=r_i \times 32768+s_i\)。

那么就有：

\[a(x) \times b(x)=(p(x) \times 32768+q(x)) \times (r(x) \times 32768+s(x)) \]

整理，得：

\[a(x) \times b(x)=p(x) \times r(x) \times 1073741824+(p(x) \times s(x)+q(x) \times r(x)) \times 32768+q(x) \times s(x) \]

这样只需要做\(4\)遍\(DFT\)和\(3\)遍\(IDFT\)了。

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cctype>
#include <algorithm>
#define rin(i,a,b) for(int i=(a);i<=(b);i++)
#define rec(i,a,b) for(int i=(a);i>=(b);i--)
#define trav(i,a) for(int i=head[(a)];i;i=e[i].nxt)
using std::cin;
using std::cout;
using std::endl;
typedef long long LL;

inline int read(){
	int x=0;char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x;
}

const int MAXN=100005;
const int MAXLEN=270005;
const long double pi=std::acos(-1);
int n,m;LL MOD;
int len,rev[MAXLEN];
int a[MAXN],b[MAXN];
struct Complex{
	long double real,imag;
	inline friend Complex operator + (Complex x,Complex y){
		return (Complex){x.real+y.real,x.imag+y.imag};
	}
	inline friend Complex operator - (Complex x,Complex y){
		return (Complex){x.real-y.real,x.imag-y.imag};
	}
	inline friend Complex operator * (Complex x,Complex y){
		return (Complex){x.real*y.real-x.imag*y.imag,x.real*y.imag+x.imag*y.real};
	}
};
Complex A[MAXLEN],B[MAXLEN],C[MAXLEN],D[MAXLEN],E[MAXLEN];

inline void fft(Complex *c,int dft){
	rin(i,0,n-1) if(i<rev[i])
		std::swap(c[i],c[rev[i]]);
	for(int mid=1;mid<n;mid<<=1){
		int r=(mid<<1);
		Complex u=(Complex){std::cos(pi/mid),dft*std::sin(pi/mid)};
		for(int l=0;l<n;l+=r){
			Complex v=(Complex){1,0};
			for(int i=0;i<mid;i++,v=v*u){
				Complex x=c[l+i],y=c[l+mid+i]*v;
				c[l+i]=x+y;
				c[l+mid+i]=x-y;
			}
		}
	}
	if(dft<0) rin(i,0,n-1)
		c[i].real/=n;
}

int main(){
	n=read(),m=read(),MOD=read();
	rin(i,0,n){
		a[i]=read();
		A[i].real=(a[i]>>15);
		B[i].real=(a[i]&32767);
	}
	rin(i,0,m){
		b[i]=read();
		C[i].real=(b[i]>>15);
		D[i].real=(b[i]&32767);
	}
	int nn=n,mm=m;
	m+=n;
	for(n=1;n<=m;n<<=1) len++;
	rin(i,1,n-1) rev[i]=((rev[i>>1]>>1)|((i&1)<<(len-1)));
	fft(A,1);
	fft(B,1);
	fft(C,1);
	fft(D,1);
	rin(i,0,n-1){
		E[i]=A[i]*D[i]+B[i]*C[i];
		A[i]=A[i]*C[i];
		B[i]=B[i]*D[i];
	}
	fft(A,-1);
	fft(B,-1);
	fft(E,-1);
	n=nn,m=mm;
	rin(i,0,n+m) printf("%lld ",((((LL)(A[i].real+0.5)%MOD)<<30)%MOD+(((LL)(E[i].real+0.5)%MOD)<<15)%MOD+(LL)(B[i].real+0.5))%MOD);
	printf("\n");
	return 0;
}

拆系数\(FFT\)（\(4\)遍\(DFT\)）

利用\(FFT\)的虚部空间可以实现一遍\(DFT\)对两个多项式进行求值，一遍\(IDFT\)对两组点值进行插值。

只需做\(2\)遍\(DFT\)和\(2\)遍\(IDFT\)即可完成，代码另外放出。

之后要补充的东西

牛顿迭代法

多项式多点求值

多项式指数，对数函数

*多项式快速插值

update on 2019.3.23：填了一些坑，看这里。

参考博客

http://blog.miskcoo.com/

推荐阅读