BZOJ 2179 FFT快速傅立叶 ——FFT

【题目分析】

快速傅里叶变换用于高精度乘法。

其实本质就是循环卷积的计算，也就是多项式的乘法。

两次蝴蝶变换。

二进制取反化递归为迭代。

单位根的巧妙取值，是的复杂度成为了nlogn

范德蒙矩阵计算逆矩阵又减轻了拉格朗日插值法的复杂度。

十分神奇。

【代码】

#include <cstdio>

#include <cstring>

#include <cstdlib>

#include <cmath>

#include <set>

#include <map>

#include <string>

#include <algorithm>

#include <vector>

#include <iostream>

#include <queue>

using namespace std;

#define maxn 200005

#define Complex cp

#define F(i,j,k) for (int i=j;i<=k;++i)

#define D(i,j,k) for (int i=j;i>=k;--i)

#define mlog 16

#define inf (0x3f3f3f3f)

int read()

{

    int x=0,f=1; char ch=getchar();

    while (ch<'0'||ch>'9') {if (ch=='-') f=-1; ch=getchar();}

    while (ch>='0'&&ch<='9') {x=x*10+ch-'0'; ch=getchar();}

    return x*f;

}

struct cp{

    double x,y;//虚数可表示为 x+y*i  i^2=-1

    cp operator + (cp a) {return (cp){x+a.x,y+a.y};}

    cp operator - (cp a) {return (cp){x-a.x,y-a.y};}

    cp operator * (cp a) {return (cp){x*a.x-y*a.y,x*a.y+y*a.x};}

}a[maxn],b[maxn],c[maxn];

double pi=acos(-1.0); // π

int n,m,rev[maxn],len,ans[maxn];

char s[maxn];

void FFT(Complex * x,int n,int f)

{

	F(i,0,n-1) if (rev[i]>i) swap(x[rev[i]],x[i]); //构造迭代的形式

	for (int m=2;m<=n;m<<=1)

	{

		Complex wn=(Complex){cos(2.0*pi/m*f),sin(2.0*pi/m*f)}; //当前的主单位根

		for (int i=0;i<n;i+=m)

		{

			Complex w=(Complex){1.0,0};

			for (int j=0;j<(m>>1);++j)

			{

				Complex u=x[i+j],v=x[i+j+(m>>1)]*w;//计算对应的多项式的值

				x[i+j]=u+v; x[i+j+(m>>1)]=u-v;

				w=w*wn;//在复数域中旋转一个角度

			}

		}

	}

}

int main()

{

    n=read();

    scanf("%s",s); F(i,0,n-1) a[i].x=s[n-1-i]-'0';

    scanf("%s",s); F(i,0,n-1) b[i].x=s[n-1-i]-'0';

    m=1; n=2*n-1;

    while (m<=n) m<<=1,len++; n=m;

    F(i,0,n-1)

    {

        int t=i,ret=0;

        F(j,1,len) ret<<=1,ret|=t&1,t>>=1;

        rev[i]=ret;

    }//二进制翻转，便于化分治为循环迭代

    FFT(a,n,1); FFT(b,n,1); //FFT

    F(i,0,n-1) c[i]=a[i]*b[i];

    FFT(c,n,-1); //IFFT

    F(i,0,n-1) ans[i]=(c[i].x/n)+0.5;//精度QAQ

    F(i,0,n-1) ans[i+1]+=ans[i]/10,ans[i]%=10;//进位QwQ

    n++;

    while (!ans[n]&&n) n--;

    D(i,n,0) printf("%d",ans[i]);

}