FFT 高精度乘法

最新推荐文章于 2023-04-05 16:06:55 发布

RicardoWWW

最新推荐文章于 2023-04-05 16:06:55 发布

阅读量781

点赞数

分类专栏： FFT 模板

本文链接：https://blog.csdn.net/qq_34637390/article/details/51743006

版权

模板同时被 2 个专栏收录

8 篇文章 0 订阅

订阅专栏

FFT

2 篇文章 0 订阅

订阅专栏

#include<cstdlib>
#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<queue>
#include<vector>
using namespace std;
#define MAXN 150005
char s[MAXN];
struct C
{
	double a,b;
	C(double _=0,double __=0) : a(_) , b(__)  { }
	friend C operator + (C x,C y){return C(x.a+y.a,x.b+y.b);}
	friend C operator - (C x, C y){return C(x.a-y.a,x.b-y.b);}
	friend C operator * (C x,C y){return C(x.a*y.a-x.b*y.b,x.a*y.b+x.b*y.a);}
}A[MAXN],B[MAXN];
#define pi 3.141592653589793238462
void FFT(C *x , int l,int on)
{
	int c=0,tmp=1; while(tmp<l){c++; tmp=tmp<<1;}
	for(int i=0;i<l;i++)
	{
		tmp=0;
		for(int j=0;j<c;j++)if(i&(1<<j))tmp|=(1<<(c-1-j));
		if(i<tmp)swap(x[i],x[tmp]);
	}
	C wn,w,u,v;
	for(int h=2;h<=l;h=h<<1)
	{
		wn=C(cos(2.0*pi*on/h),sin(2.0*pi*on/h));
		for(int i=0;i<l;i+=h)
		{
			w=C(1.0,0.0);
			for(int j=0;(j<<1)<h;j++)
			{
				u=x[i+j]; v=w*x[i+j+h/2];
				x[i+j]=u+v; x[i+j+h/2]=u-v;
				w=w*wn;
			}
		}
	}
	if(on==-1)
	{
		for(int i=0;i<l;i++)x[i].a=x[i].a/l;
	}
	return ;
}
int n1,n2,l,Ans[MAXN];
int main()
{
	freopen("in.txt","r",stdin);
	while(scanf("%s",s)==1)
	{
		n1=strlen(s);
		for(int i=0;i<n1;i++)A[n1-1-i]=C(s[i]-'0',0);
		scanf("%s",s); n2=strlen(s);
		for(int i=0;i<n2;i++)B[n2-1-i]=C(s[i]-'0',0);
		l=1;
		while(l<2*n1 || l<2*n2)l=l<<1;
		for(int i=n1;i<l;i++)A[i]=C(0,0);
		for(int i=n2;i<l;i++)B[i]=C(0,0);
		FFT(A,l,1); FFT(B,l,1);
		for(int i=0;i<l;i++)A[i]=A[i]*B[i];
		FFT(A,l,-1);
		for(int i=0;i<l;i++)Ans[i]=(int)(A[i].a+0.4);
		for(int i=0;i<l;i++)if(Ans[i]>=0){Ans[i+1]+=Ans[i]/10; Ans[i]=Ans[i]%10;}
		for(int i=l-1;i>=0;i--)
		if(Ans[i])
		{
			for(int j=i;j>=0;j--)putchar('0'+Ans[j]);
			putchar('\n');
			break;
		}
		else if(i==0)
		{
			putchar('0'); putchar('\n');
		}
	}
	return 0;
}