A*B problem（FFT）

A*B problem（FFT）

设两个多项式(A(x))和(B(x))，它们的系数镜像反转一下，得到的多项式是(A'(x))和(B'(x))。那么(C(x)=A(x)*B(x))和(C'(x)=A'(x)*B'(x))的系数也是镜像反转的。这个，，感性理解一下吧。

于是倒过来搞会很方便。由于是大整数乘法，算完从个位到最高位进一下位就行了。注意数组要开四倍！

#include <cmath>
#include <cctype>
#include <cstdio>
using namespace std;

const int maxn=6e4+5;
const double pi=3.1415926535898;

struct Cpx{
    double x, y;
    Cpx (double t1=0, double t2=0){ x=t1, y=t2; }
}A[maxn*4], B[maxn*4], C[maxn*4];
Cpx operator +(Cpx &a, Cpx &b){ return Cpx(a.x+b.x, a.y+b.y); }
Cpx operator -(Cpx &a, Cpx &b){ return Cpx(a.x-b.x, a.y-b.y); }
Cpx operator *(Cpx &a, Cpx &b){ return Cpx(a.x*b.x-a.y*b.y, a.x*b.y+a.y*b.x); }
void swap(Cpx &x, Cpx &y){ Cpx t=x; x=y; y=t; }

int n, r[maxn*4], limit=1, l, ans[maxn*4];  //这里也要乘4！

void fdft(Cpx *a, int n, int flag){
    for (int i=0; i<n; ++i) if (i<r[i]) swap(a[i], a[r[i]]);
    for (int mlen=1; mlen<n; mlen<<=1){
        Cpx w1(cos(pi/mlen), flag*sin(pi/mlen)), x, y;
        for (int i=0; i<n; i+=(mlen<<1)){  //起点
            Cpx w(1, 0);
            for (int j=i; j<i+mlen; ++j, w=w*w1){  //遍历区间系数转点值
                x=a[j], y=w*a[j+mlen];
                a[j]=x+y; a[j+mlen]=x-y; }
        }
    }
}

int main(){
    scanf("%d", &n);  //倒着存方便补零 两个n-1次的多项式
    char c; while (!isdigit(c=getchar()));
    for (int i=0; i<n; ++i, c=getchar()) A[n-i-1].x=c-48;
    while (!isdigit(c=getchar()));
    for (int i=0; i<n; ++i, c=getchar()) B[n-i-1].x=c-48;
    while (limit<=2*n-2) limit<<=1, ++l;  //需要2×n-2个点，但是有2×n-1个系数
    for (int i=1; i<limit; ++i)
        r[i]=(r[i>>1]>>1)+((i&1)<<(l-1));  //1<<r就是倒数第r+1位
    fdft(A, limit, 1); fdft(B, limit, 1);
    for (int i=0; i<limit; ++i) C[i]=A[i]*B[i];
    fdft(C, limit, -1);
    for (int i=0; i<limit; ++i) ans[i]+=+C[i].x/limit+0.5;
    for (int i=0; i<limit; ++i)
        ans[i+1]+=ans[i]/10, ans[i]%=10;
    bool flag=false;
    for (int i=limit; i>=0; --i){
        if (ans[i]>0) flag=true;
        if (flag) printf("%d", ans[i]);
    }
    return 0;
}