代码

// 1004a 高精度乘法2 - single 64-bit NTT (Goldilocks prime 2^64-2^32+1), base 1e8
#include <stdio.h>
#include <stdint.h>
#include <string.h>
#include <stdlib.h>

typedef uint64_t u64;
typedef __uint128_t u128;
typedef __int128_t i128;

static const u64 P = 18446744069414584321ULL; // 2^64 - 2^32 + 1
static const i128 P128 = (i128)P;
static const u64 W4096 = 0xF2C35199959DFCB6ULL; // primitive 4096-th root of unity
static const u64 WINV   = 0x9AF01E431FBD6EA0ULL; // W4096^{-1}

// Goldilocks reduction: (a*b) mod P, a,b in [0,P)
static inline u64 reduce(u64 a, u64 b){
    u128 x = (u128)a * b;
    u64 x0 = (u64)x;
    u64 x1 = (u64)(x >> 64);
    u64 h = x1 >> 32;
    u64 l = x1 & 0xFFFFFFFFULL;
    i128 r = (i128)x0 - (i128)h - (i128)l + ((i128)l << 32);
    if (r < 0) r += P128;
    if (r >= P128) r -= P128;
    return (u64)r;
}

// (u+v) mod P, u,v in [0,P)  (u+v may overflow 2^64)
static inline u64 addmod(u64 u, u64 v){
    u64 s = u + v;
    u64 t = s - P;
    u64 m = (s < u) | (s >= P);
    return m ? t : s;
}
// (u-v) mod P, u,v in [0,P)
static inline u64 submod(u64 u, u64 v){
    u64 d = u - v;
    return (u < v) ? (d + P) : d;
}

#define N 4096
static u64 A[N], B[N];
static u64 roots[N], iroots[N]; // roots[k] = w^k, iroots[k] = w^{-k}

static void build_roots(void){
    roots[0] = 1;
    for(int i=1;i<N;i++) roots[i] = reduce(roots[i-1], W4096);
    iroots[0] = 1;
    for(int i=1;i<N;i++) iroots[i] = reduce(iroots[i-1], WINV);
}

// forward DIF NTT, in-place, output bit-reversed
static void ntt_fwd(u64* a){
    for(int len=N; len>1; len>>=1){
        int half = len>>1;
        int step = N / len;
        for(int i=0;i<N;i+=len){
            for(int j=0;j<half;j++){
                u64 u = a[i+j];
                u64 v = a[i+j+half];
                a[i+j] = addmod(u, v);
                a[i+j+half] = reduce(submod(u, v), roots[j*step]);
            }
        }
    }
}

// inverse DIT NTT, in-place, input bit-reversed, output natural
static void ntt_inv(u64* a){
    for(int len=2; len<=N; len<<=1){
        int half = len>>1;
        int step = N / len;
        for(int i=0;i<N;i+=len){
            for(int j=0;j<half;j++){
                u64 u = a[i+j];
                u64 v = reduce(a[i+j+half], iroots[j*step]);
                a[i+j] = addmod(u, v);
                a[i+j+half] = submod(u, v);
            }
        }
    }
    u64 ninv = 0; // n^{-1} mod P, n=4096. compute pow.
    // 4096 = 2^12; compute inverse via reduce
    u64 e = 4096;
    // compute e^{-1} mod P using Fermat (P is prime)
    // but simpler: precompute. 4096^{-1} mod P
    // We'll compute via extended approach at runtime once.
    static u64 ninv_cache = 0;
    if(!ninv_cache){
        u128 r=1, b=4096, ee=P-2;
        while(ee){ if(ee&1) r=(u128)r*b%P; b=(u128)b*b%P; ee>>=1; }
        ninv_cache = (u64)r;
    }
    ninv = ninv_cache;
    for(int i=0;i<N;i++) a[i] = reduce(a[i], ninv);
}

static char ibuf[30000];
static int ilen = 0;

static char obuf[21000];
static int olen = 0;

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_roots();
    build_tab3();

    // read all input
    ilen = fread(ibuf, 1, sizeof(ibuf), stdin);

    // find two tokens
    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;

    // strip leading zeros
    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++;

    // parse to base 1e8 limbs little-endian
    int na=0, nb=0;
    {
        int pos = a_end;
        while(pos > sa){
            int start = pos-8; if(start < sa) start = sa;
            u64 v=0;
            for(int i=start;i<pos;i++) v = v*10 + (ibuf[i]-'0');
            A[na++] = v; pos = start;
        }
    }
    {
        int pos = b_end;
        while(pos > sb){
            int start = pos-8; if(start < sb) start = sb;
            u64 v=0;
            for(int i=start;i<pos;i++) v = v*10 + (ibuf[i]-'0');
            B[nb++] = v; pos = start;
        }
    }

    ntt_fwd(A);
    ntt_fwd(B);
    for(int i=0;i<N;i++) A[i] = reduce(A[i], B[i]);
    ntt_inv(A);

    // carry in base 1e8
    u64 carry = 0;
    int outlen = na + nb - 1;
    // coefficients A[0..outlen-1] are the result limbs (base 1e8) after carry
    for(int i=0;i<outlen;i++){
        u64 c = A[i] + carry;
        A[i] = c % 100000000ULL;
        carry = c / 100000000ULL;
    }
    while(carry){
        A[outlen++] = carry % 100000000ULL;
        carry /= 100000000ULL;
    }
    int hi = outlen-1;
    while(hi>0 && A[hi]==0) hi--;

    // output: most significant limb, then 8-digit groups
    // build output string into obuf
    // most significant limb (no padding)
    {
        u64 v = A[hi];
        // print up to 8 digits, no leading zero
        char tmp[9];
        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--){
        u64 v = A[i];
        u64 g2 = v / 1000000;        // 2 digits
        u64 g1 = (v/1000) % 1000;    // 3 digits
        u64 g0 = v % 1000;           // 3 digits
        obuf[olen++] = '0' + g2/10;
        obuf[olen++] = '0' + g2%10;
        const char* p = tab3[g1];
        obuf[olen++]=p[0]; obuf[olen++]=p[1]; obuf[olen++]=p[2];
        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;
}

评测结果