概念

多项式乘法时,我们发现暴力乘十分缓慢,但是点值乘十分快速。考虑求 \(A\)\(B\) 的卷积。

一个 \(n\) 次多项式可以被 \(n+1\) 个点确定。

设多项式 \(A(x)\) 的系数为 \((a_0,a_1,\cdots,a_n)\)

对其奇偶分类得 \(A(x)=\sum\limits a_{2i}*x^{2i}+\sum a_{2i+1}*x^{2i+1}\)

提取 \(n\)\(A(x)=\sum\limits a_{2i}*x^{2i}+x*\sum a_{2i+1}*x^{2i}\)

化简得 \(A(x)=\sum\limits a_{2i}*(x^2)^i+x*\sum a_{2i+1}*(x^2)^i\)

\((\frac{n}{2}-1)\) 次多项式 \(A_0=\sum\limits a_{2i}*x^{2i},A_1=\sum\limits a_{2i+1}*x^{2i}\)

\(A(x)=A_0(x^2)+x*A_1(x^2)\)

我们发现可以递归

但是时间复杂度不对

然后是一些高深件

复数

定义 \(i^2=-1\)

则所有形如 \(z=a+b*i\) 的数 \(z\) 组成的集合为复数集,记为 \(C\)。而 \(a\) 称为 \(z\) 的实部, \(b\) 为虚部。

复平面

类似坐标系,纵虚横实。复数 \(z\) 对应从原点指向 \((a,b)\) 的向量。幅角为实轴横方向与向量的夹角,记为 \(θ\)

加法:与向量加法一样

乘法:复数相乘,模长相乘,幅角相加。\((a+bi)*(c+di)=(ac-bd)+(bc+ad)i\)

单位根

在复平面上,以原点为圆心,\(1\) 为半径作圆,所得的圆叫单位圆。以圆点为起点,圆的 \(n\) 等分点为终点,做 \(n\) 个向量,设幅角为正且最小的向量对应的复数为 \(\omega_n\),称为 \(n\) 次单位根。剩余 \(n-1\) 个点则为 \(\omega_n^2,\omega_n^3,\omega_n^4,\cdots\)

单位根的幅角为周角的 \(\dfrac{1}{n}\)

欧拉公式有:\(w_n^k=\cos k\frac{2π}{n}+i\sin k\frac{2π}{n}\)

性质 \(1:\) \(\omega_{n}^{k}=\omega_{2n}^{2k}\)

性质 \(2:\) \(\omega _{n}^{k+\frac{n}{2}}=-\omega _{n}^{k}\)

性质 \(3:\) \(w_n^n=1\)

图上自证不难

FFT

承接高深件

高能推柿子

\(A(x)=A_0(x^2)+x*A_1(x^2)\)

带入 \(x=w_n^k(k<\frac{1}{2}n)\)

\[A(w_n^k)=A_0(w_n^{2k})+w_n^k*A_1(w_n^{2k}) \]

带入 \(x=w_n^{k+\frac{n}{2}}(k<\frac{1}{2}n)\)

\[A(w_n^{k+\frac{n}{2}})=A_0(w_n^{2k+n})+w_n^{k+\frac{n}{2}}*A_1(w_n^{2k+n}) \]

\[A(w_n^k)=A_0(w_{\frac{n}{2}}^{k})+w_n^k*A_1(w_{\frac{n}{2}}^{k}) \]

\[A(w_n^k)=A_0(w_{n}^{2k})-w_n^k*A_1(w_{n}^{2k}) \]

不难发现,两式的值可以互相转换,且两室分别处理了一半区间。递归处理,复杂度 \(O(n\log n)\)

IFFT

函数转点并相乘后,需要转化回系数。此时使用傅里叶逆变换。

拉差

拉差不行,还是要推柿子。

设向量 \(c_1,c_2,\cdots c_{n-1}\) 满足 \(c_k=\sum\limits^{n}_{j=1} y_j(w_n^{-k})^j\)

推柿子懒了,多项式早日消失

最终

\(c_k=n*a_k\)

\(a_k=\frac{c_k}{n}\)

代码

递归实现即可

使用STL:complex实现复数

但是但是,递归常数很大。

改成迭代就可以了。

尻尻板子:

#include<iostream>
#include<cstdio>
#include<complex>
#include<cmath>
#include<algorithm>
using namespace std;
double pi=acos(-1);
int tax[5000001];
void Rev(int n)
{
	int d=n>>1;
	int len=-1;
	tax[++len]=0;
	tax[++len]=d;
	for(int i=2;i<=n;i<<=1)
	{
		d>>=1;
		for(int p=0;p<i;p++)
		{
			tax[++len]=tax[p]|d;
		}
	}
}
void FFT(complex<double> *A,int N)
{
	for(int i=1;i<N;i++)
	{
		if(tax[i]>i)
		{
			swap(A[i],A[tax[i]]);
		} 
	}
	for(int len=2,M=1;len<=N;M=len,len<<=1)
	{
		complex<double> W(cos(pi/M),sin(pi/M));
		complex<double> w(1.0,0.0);
		for(int l=0,r=len-1;r<=N;l+=len,r+=len)
		{
			complex<double> w0=w;
			for(int p=l;p<l+M;p++)
			{
				complex<double> x=A[p]+w0*A[p+M];
				complex<double> y=A[p]-w0*A[p+M];
				A[p]=x;
				A[p+M]=y;
				w0*=W;
			}
		}
	}
}
void IFFT(complex<double> *A,int N)
{
	FFT(A,N);
	reverse(A+1,A+N);
}
int n,m;
complex<double> F[5000001],G[5000001];
int main()
{
	scanf("%d%d",&n,&m);
	for(int i=0,_;i<=n;i++)
	{
		scanf("%d",&_);
		F[i]=_;
	}
	for(int i=0,_;i<=m;i++)
	{
		scanf("%d",&_);
		G[i]=_;
	}
	int p2=1;
	while(p2<=m+n) p2<<=1;
	Rev(p2);
	FFT(F,p2);
	FFT(G,p2);
	for(int i=0;i<p2;i++)
	{
		F[i]*=G[i];
	}
	IFFT(F,p2);
	for(int i=0;i<=n+m;i++)
	{
		printf("%d ",int(F[i].real()/p2+0.5));
	}
}

NTT板子

#include<iostream>
#include<cstdio>
#include<complex>
#include<cmath>
#include<algorithm>
#define int long long
#define mod 998244353
using namespace std;
double pi=acos(-1);
int tax[5000001];
int n,m;
void Rev(int pn)
{
	for(int i=0;i<pn;i++)
	{
		tax[i]=tax[i>>1]>>1|((i&1)<<(max((int)ceil(log2(n+m)),1ll)-1));
	}
}
long long poww(int a,int b)
{
	int re=1;
	while(b)
	{
		if(b&1)
		{
			re*=a;
			re%=mod;
		}
		a*=a;
		a%=mod;
		b>>=1;
	}
	return re;
}
void NTT(long long *a,int n,int type) //type:1正0反 
{
    for(int i=0;i<n;++i)
    {
        if(i<tax[i])
        {
            swap(a[i],a[tax[i]]);
        }
    }
    for(int i=1;i<n;i<<=1)
    {
        long long gn=poww(3,(mod-1)/(i<<1));
        for(int j=0;j<n;j+=(i<<1))
        {
            long long g0=1;
            for(int k=0;k<i;++k,g0=g0*gn%mod)
            {
                long long x=a[j+k],y=g0*a[i+j+k]%mod;
                a[j+k]=(x+y)%mod;
                a[i+j+k]=(x-y+mod)%mod;
            }
        }
    }
    if(type!=1)
    {
		int nn=poww(n,mod-2);
		reverse(a+1,a+n);
		for(int i=0;i<n;i++)
		{
			a[i]=1ll*a[i]*nn%mod;
		}
	}
}
long long F[5000001],G[5000001];
signed main()
{
	scanf("%lld%lld",&n,&m);
	for(int i=0,_;i<=n;i++)
	{
		scanf("%lld",&_);
		F[i]=_;
	}
	for(int i=0,_;i<=m;i++)
	{
		scanf("%lld",&_);
		G[i]=_;
	}
	int p2=1;
	while(p2<=m+n) p2<<=1;
	Rev(p2);
	NTT(F,p2,1);
	NTT(G,p2,1);
	for(int i=0;i<p2;i++)
	{
		F[i]*=G[i];
	}
	NTT(F,p2,0);
	for(int i=0;i<=n+m;i++)
	{
		printf("%lld ",F[i]);
	}
}
posted on 2023-05-04 18:23  lizhous  阅读(10)  评论(0编辑  收藏  举报