快速傅里叶变换（FFT）

原理部分

代码实现

多项式乘法

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<vector>
#include<cmath>
#define N 4000009
#define reg register
using namespace std;
typedef long long ll;
const double pi=3.141592653589793;
struct complex
{
	double x,y;
	complex(double x=0,double y=0):x(x),y(y) {}

	inline complex operator + (const complex b) const
	{
		return complex(x+b.x,y+b.y);
	}

	inline complex operator - (const complex b) const
	{
		return complex(x-b.x,y-b.y);
	}

	inline complex operator * (const complex b) const
	{
		complex res;
		res.x = x*b.x-y*b.y;
		res.y = x*b.y+y*b.x;
		return res;
	}
};

int n,m;
complex F[N];
int rev[N];
void print(int x);
inline int read()
{
	int x=0,w=0;
	char ch=0;
	while(!isdigit(ch)) w|=ch=='-',ch=getchar();
	//while(isdigit(ch)) x=(x<<1)+(x<<3)+(ch^48),ch=getchar();以前以为位运算能快点，但问了大佬才发现，其实没啥差别
	while(isdigit(ch)) x=x*10+ch-'0',ch=getchar();
	return w?-x:x;
}
inline void FFT(complex *a,int type,int lim)  //没什么好说的,FFT的模板
{
	for(reg int i=1; i<=lim; ++i)
	{
		if(i>=rev[i]) continue;
		swap(a[i],a[rev[i]]);
	}
	reg complex rt,w,x,y;
	for(reg int mid=1; mid<lim; mid<<=1)
	{
		reg int r = mid<<1;
		rt = complex(cos(pi/mid),type*sin(pi/mid));
		for(reg int j=0; j<lim; j+=r)
		{
			w = complex(1,0);
			for(reg int k=0; k<mid; ++k)
			{
				x = a[j|k];
				y = w*a[j|k|mid];
				a[j|k] = x+y;
				a[j|k|mid] = x-y;
				w = w*rt;
			}
		}
	}
	if(type==1) return;
	for(reg int i=0; i<=lim; ++i)
	{
		a[i].y = a[i].y/lim;
		a[i].x = a[i].x/lim;
	}
}

int main()
{
	int qwq,t,lim = 1,l = -1;
	n=read(),m=read();
	for(reg int i=0; i<=n; ++i) F[i].x=read();
	for(reg int i=0; i<=m; ++i) F[i].y=read(); //把G(x)放到F(x)的虚部上
	t = n+m;
	n = max(n,m);
	while(lim<=(n<<1))
	{
		lim <<= 1;
		++l;
	}
	for(reg int i=1; i<=lim; ++i)
		rev[i] = (rev[i>>1]>>1)|((i&1)<<l);
	FFT(F,1,lim);//FFT
	for(reg int i=0; i<=lim; ++i)
		F[i] = F[i]*F[i]; //求出F(x)^2
	FFT(F,-1,lim);//IFFT
	for(reg int i=0; i<=t; ++i)
	{
		qwq = F[i].y/2+0.5; //虚部取出来除2，注意要+0.5，否则精度会有问题
		print(qwq);
		putchar(' ');
	}
	return 0;
}

void print(int x)
{
	if(x>9) print(x/10);
	putchar(x%10+'0');
}

超大数乘法

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<vector>
#include<cmath>
#define reg register
using namespace std;
typedef long long ll;
const int maxn=4000009;
const double pi=3.141592653589793;
struct complex
{
	double x,y;
	complex(double x=0,double y=0):x(x),y(y) {}

	inline complex operator + (const complex b) const
	{
		return complex(x+b.x,y+b.y);
	}

	inline complex operator - (const complex b) const
	{
		return complex(x-b.x,y-b.y);
	}

	inline complex operator * (const complex b) const
	{
		complex res;
		res.x = x*b.x-y*b.y;
		res.y = x*b.y+y*b.x;
		return res;
	}
};

long long n,m;
complex F[maxn];
long long rev[maxn];
void print(long long x);
inline long long read()
{
	long long x=0,w=0;
	char ch=0;
	while(!isdigit(ch)) w|=ch=='-',ch=getchar();
	//while(isdigit(ch)) x=(x<<1)+(x<<3)+(ch^48),ch=getchar();以前以为位运算能快点，但问了大佬才发现，其实没啥差别
	while(isdigit(ch)) x=x*10+ch-'0',ch=getchar();
	return w?-x:x;
}
inline void FFT(complex *a,long long type,long long lim)  //没什么好说的,FFT的模板
{
	for(reg long long i=1; i<=lim; ++i)
	{
		if(i>=rev[i]) continue;
		swap(a[i],a[rev[i]]);
	}
	reg complex rt,w,x,y;
	for(reg long long mid=1; mid<lim; mid<<=1)
	{
		reg long long r = mid<<1;
		rt = complex(cos(pi/mid),type*sin(pi/mid));
		for(reg long long j=0; j<lim; j+=r)
		{
			w = complex(1,0);
			for(reg long long k=0; k<mid; ++k)
			{
				x = a[j|k];
				y = w*a[j|k|mid];
				a[j|k] = x+y;
				a[j|k|mid] = x-y;
				w = w*rt;
			}
		}
	}
	if(type==1) return;
	for(reg long long i=0; i<=lim; ++i)
	{
		a[i].y = a[i].y/lim;
		a[i].x = a[i].x/lim;
	}
}
char s1[maxn],s2[maxn];
int output[maxn],cnt;
int main()
{
	long long t,lim = 1,l = -1;
	cin>>s1>>s2;
	n=strlen(s1),m=strlen(s2);
	for(reg long long i=0; i<n; ++i) F[i].x=s1[n-i-1]-'0';
	for(reg long long i=0; i<m; ++i) F[i].y=s2[m-i-1]-'0'; //把G(x)放到F(x)的虚部上
	t = n+m;
	n = max(n,m);
	while(lim<=(n<<1))
	{
		lim <<= 1;
		++l;
	}
	for(reg long long i=1; i<=lim; ++i)
		rev[i] = (rev[i>>1]>>1)|((i&1)<<l);
	FFT(F,1,lim);//FFT
	for(reg long long i=0; i<=lim; ++i)
		F[i] = F[i]*F[i]; //求出F(x)^2
	FFT(F,-1,lim);//IFFT
	for(reg long long i=0; i<=t; ++i)
	{
		output[i]+= F[i].y/2+0.5; //虚部取出来除2，注意要+0.5，否则精度会有问题
		output[i+1]+=output[i]/10;
		output[i]%=10;

	}
	int i;
	for(i=n+m; !output[i]&&i>=0; i--); //去掉前导零
	if(i==-1)printf("0");//特判长度为0的情况
	for(; i>=0; i--) //输出这个十进制数
	{
		printf("%d",output[i]);
	}
	putchar('\n');
	return 0;
}