代码

// 1004a 高精度乘法2 - double FFT, base 1e5, size 4096, DIF/DIT no bit-reversal
#include <stdio.h>
#include <stdint.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>

typedef uint32_t u32;
typedef uint64_t u64;

#define N 4096
static double fre[N], fim[N];
static double gre[N], gim[N];
static double wcos[N/2], wsin[N/2];

static void fft_fwd(double* re, double* im){
    for(int len=N; len>1; len>>=1){
        int half=len>>1, step=N/len;
        for(int i=0;i<N;i+=len){
            for(int j=0;j<half;j++){
                double wr = wcos[j*step];
                double wi = wsin[j*step];
                double ur = re[i+j], ui = im[i+j];
                double vr = re[i+j+half], vi = im[i+j+half];
                double sr = ur+vr, si = ui+vi;
                double dr = ur-vr, di = ui-vi;
                re[i+j] = sr; im[i+j] = si;
                re[i+j+half] = dr*wr - di*wi;
                im[i+j+half] = dr*wi + di*wr;
            }
        }
    }
}

static void fft_inv(double* re, double* im){
    for(int len=2; len<=N; len<<=1){
        int half=len>>1, step=N/len;
        for(int i=0;i<N;i+=len){
            for(int j=0;j<half;j++){
                double wr = wcos[j*step];
                double wi = -wsin[j*step];
                double ur = re[i+j], ui = im[i+j];
                double vr = re[i+j+half]*wr - im[i+j+half]*wi;
                double vi = re[i+j+half]*wi + im[i+j+half]*wr;
                re[i+j] = ur+vr; im[i+j] = ui+vi;
                re[i+j+half] = ur-vr; im[i+j+half] = ui-vi;
            }
        }
    }
    double invn = 1.0/N;
    for(int i=0;i<N;i++){ re[i]*=invn; im[i]*=invn; }
}

static char ibuf[30000];
static int ilen;
static char obuf[21000];
static int olen;
static char tab3[1000][4];
static void build_tab3(void){
    for(int i=0;i<1000;i++){ int v=i; tab3[i][2]='0'+v%10; v/=10; tab3[i][1]='0'+v%10; v/=10; tab3[i][0]='0'+v%10; }
}

int main(){
    build_tab3();
    const double PI = acos(-1.0);
    for(int k=0;k<N/2;k++){ wcos[k]=cos(2*PI*k/N); wsin[k]=sin(2*PI*k/N); }

    ilen = fread(ibuf, 1, sizeof(ibuf), stdin);
    int pa=0;
    while(pa<ilen && (ibuf[pa]==' '||ibuf[pa]=='\n'||ibuf[pa]=='\r'||ibuf[pa]=='\t')) pa++;
    int a_start=pa; while(pa<ilen && ibuf[pa]!=' '&&ibuf[pa]!='\n'&&ibuf[pa]!='\r'&&ibuf[pa]!='\t') pa++;
    int a_end=pa; while(pa<ilen && (ibuf[pa]==' '||ibuf[pa]=='\n'||ibuf[pa]=='\r'||ibuf[pa]=='\t')) pa++;
    int b_start=pa; while(pa<ilen && ibuf[pa]!=' '&&ibuf[pa]!='\n'&&ibuf[pa]!='\r'&&ibuf[pa]!='\t') pa++;
    int b_end=pa;

    int sa=a_start; while(sa<a_end-1 && ibuf[sa]=='0') sa++;
    int sb=b_start; while(sb<b_end-1 && ibuf[sb]=='0') sb++;

    int na=0, nb=0;
    { int pos=a_end; while(pos>sa){ int st=pos-5; if(st<sa)st=sa; u32 v=0; for(int i=st;i<pos;i++) v=v*10+(ibuf[i]-'0'); fre[na++]=(double)v; pos=st; } }
    { int pos=b_end; while(pos>sb){ int st=pos-5; if(st<sb)st=sb; u32 v=0; for(int i=st;i<pos;i++) v=v*10+(ibuf[i]-'0'); gre[nb++]=(double)v; pos=st; } }

    fft_fwd(fre, fim);
    fft_fwd(gre, gim);
    for(int i=0;i<N;i++){
        double ar=fre[i], ai=fim[i], br=gre[i], bi=gim[i];
        fre[i] = ar*br - ai*bi;
        fim[i] = ar*bi + ai*br;
    }
    fft_inv(fre, fim);

    int outlen = na + nb - 1;
    u64 carry = 0;
    static u32 outlimbs[4002];
    int nout = 0;
    for(int i=0;i<outlen;i++){
        u64 c = (u64)llround(fre[i]);
        c += carry;
        outlimbs[nout++] = (u32)(c % 100000ULL);
        carry = c / 100000ULL;
    }
    while(carry){ outlimbs[nout++] = (u32)(carry % 100000ULL); carry /= 100000ULL; }
    int hi = nout-1;
    while(hi>0 && outlimbs[hi]==0) hi--;

    {
        u32 v = outlimbs[hi];
        char tmp[6]; int t=0;
        do { tmp[t++]='0'+v%10; v/=10; } while(v);
        while(t>0) obuf[olen++] = tmp[--t];
    }
    for(int i=hi-1;i>=0;i--){
        u32 v = outlimbs[i];
        u32 g1 = v / 1000;
        u32 g0 = v % 1000;
        obuf[olen++] = '0' + g1/10;
        obuf[olen++] = '0' + g1%10;
        const char* p = tab3[g0];
        obuf[olen++]=p[0]; obuf[olen++]=p[1]; obuf[olen++]=p[2];
    }
    obuf[olen++] = '\n';
    fwrite(obuf, 1, olen, stdout);
    return 0;
}

评测结果