大數乘法

mul_fft.cpp

#include <bits/stdc++.h>
using namespace std;

using BigNum = vector<int>;
BigNum& trim(BigNum &b) {
    while(b.size() > 1 && b.back()==0)
        b.pop_back();
    return b;
}

BigNum stob(const string &s) {
    BigNum b;
    for(char c:s | views::reverse) b.emplace_back(c - '0');
    return trim(b);
}

void show(const BigNum &b) {
    for(auto i : b | views::reverse) cout << i;
    cout << '\n';
}

BigNum& advance(BigNum& b) {
    for(int i=0;i+1<b.size();++i){
        b[i+1] += b[i]/10;
        b[i]%=10;
    }
    return b;
}

BigNum mul(const BigNum &a, const BigNum &b)
{
    BigNum c(a.size() + b.size());
    for(int i=0;i<a.size();++i)
        for(int j=0;j<b.size();++j)
            c[i+j] += a[i] * b[j];
    
    return trim(advance(c));
}

using CBigNum = vector<complex<double>>;
CBigNum btoc(const BigNum &b, int N) {
    CBigNum c;
    for(auto i:b) c.emplace_back(i);
    while(c.size() < N) c.emplace_back(0);
    return c;
}

BigNum ctob(const CBigNum &c) {
    BigNum b;
    for(auto i:c) b.emplace_back(round(i.real()));
    return b;
}

BigNum mul2(const BigNum &a_, const BigNum &b_)
{
    int N = a_.size() + b_.size();

    CBigNum a = btoc(a_, N);
    CBigNum b = btoc(b_, N);

    auto fft = [=](const CBigNum &a, bool inv)-> CBigNum {
        constexpr double PI = std::numbers::pi; //C++20, or acos(-1)
        constexpr auto I = complex<double>(0, 1);
        int flag = inv ? -1 : 1;
        
        CBigNum res(a.size());
        for(int i=0;i<N;++i)
        {
            for(int j=0;j<N;++j)
                res[i] += a[j] * exp( flag* i*j*(2*PI/N) * I);
        }
        if (inv) for(auto &c:res) c/=N;

        return res;
    };

    a = fft(a, false);
    b = fft(b, false);

    for(int i=0;i<N;++i)
        a[i] *= b[i];

    a = fft(a, true);
    auto res = ctob(a);
    return trim(advance(res));
}

CBigNum fft_r(const CBigNum &a, bool inv, int N, bool first = true) {
    constexpr double PI = std::numbers::pi; //C++20, or acos(-1)
    constexpr auto I = complex<double>(0, 1);
    int flag = inv ? -1 : 1;

    CBigNum res(a);
    if (N==1) 
        return a;

    CBigNum even, odd;
    for(int i=0;i<N;++i) {
        if (i % 2 == 0) even.emplace_back(a[i]);
        else odd.emplace_back(a[i]);
    }

    even = fft_r(even, inv, N/2, false);
    odd  = fft_r(odd , inv, N/2, false);

    for(int i=0;i<N/2;++i)
    {
        res[i]     = even[i] + odd[i] * exp( flag * i * (2*PI/N) * I);
        res[i+N/2] = even[i] - odd[i] * exp( flag * i * (2*PI/N) * I);
    }

    if (inv&&first) for(auto &c:res) c/=N;
    return res;
}

BigNum mul3(const BigNum &a_, const BigNum &b_)
{
    int N = a_.size() + b_.size();
    int p2 = 1;
    while (p2 < N) p2*=2;
    N = p2;

    CBigNum a = btoc(a_, N);
    CBigNum b = btoc(b_, N);

    a = fft_r(a, false, N);
    b = fft_r(b, false, N);

    for(int i=0;i<N;++i)
        a[i] *= b[i];

    a = fft_r(a, true, N);
    auto res = ctob(a);
    return trim(advance(res));
}

inline unsigned int bit_rev(unsigned int a,int N) {
    a=((a&0x55555555U)<<1)|((a&0xAAAAAAAAU)>>1);
    a=((a&0x33333333U)<<2)|((a&0xCCCCCCCCU)>>2);
    a=((a&0x0F0F0F0FU)<<4)|((a&0xF0F0F0F0U)>>4);
    a=((a&0x00FF00FFU)<<8)|((a&0xFF00FF00U)>>8);
    a=((a&0x0000FFFFU)<<16)|((a&0xFFFF0000U)>>16);
    return a>>(32-N);
}

CBigNum fft(const CBigNum &a, bool inv, int N) {
    constexpr double PI = std::numbers::pi; //C++20, or acos(-1)
    constexpr auto I = complex<double>(0, 1);
    int flag = inv ? -1 : 1;

    int bit_len = __lg(N); // log2
    CBigNum in(N);
    CBigNum out(N);

    for(int i=0;i<N;++i)
        in[bit_rev(i, bit_len)] = a[i];

    for(int n=2 ; n<=N ; n*=2) {
        for(int base=0 ; base<N ; base+=n)
        {
            for(int i=0;i<n/2;++i)
            {
                out[base+i]     = in[base+i] + in[base+i+n/2] * exp( flag * i * (2*PI/n) * I);
                out[base+i+n/2] = in[base+i] - in[base+i+n/2] * exp( flag * i * (2*PI/n) * I);
            }
        }
        in.swap(out);
    }

    if (inv) for(auto &c:in) c/=N;
    return in;
}

BigNum mul4(const BigNum &a_, const BigNum &b_)
{
    int N = a_.size() + b_.size();
    int p2 = 1;
    while (p2 < N) p2*=2;
    N = p2;

    CBigNum a = btoc(a_, N);
    CBigNum b = btoc(b_, N);

    a = fft(a, false, N);
    b = fft(b, false, N);

    for(int i=0;i<N;++i)
        a[i] *= b[i];

    a = fft(a, true, N);
    auto res = ctob(a);
    return trim(advance(res));
}

int main()
{
    string a, b;
    while( cin >> a >> b )
    {
        auto A = stob(a);
        auto B = stob(b);

        show(mul(A,B));
        show(mul2(A,B));
        show(mul3(A,B));
        show(mul4(A,B));
        cout << "long long =" << stoll(a) * stoll(b) << endl;
    }
    return 0;
}