高精度乘法C++

SnowDreamXUE發表於2024-10-21

1.高精度乘高精度的簡單演算法

思想:倒置相乘,統一處理進位,還原。

複雜度:$o(n^2)$

// By SnowDream
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1e5+10;
string s1,s2;
int n1[N],n2[N],n3[N];//n1儲存被乘數,n2儲存乘數,n3儲存積

void mul()
{
    int l1=(int)s1.size();
    int l2=(int)s2.size();
    //讀取字串轉換為int型別並倒置儲存至陣列
    for(int i=0;i<l1;++i)
    {
        n1[l1-1-i]=s1[i]-'0';
    }
    for(int i=0;i<l2;++i)
    {
        n2[l2-1-i]=s2[i]-'0';
    }

    for(int i=0;i<l1;++i)
    {
        for(int j=0;j<l2;++j)
        {
            n3[i+j]+=n1[i]*n2[j];
        }
    }
    //處理進位
    int l_max=l1+l2;
    for(int i=0;i<l_max;++i)
    {
        n3[i+1]+=n3[i]/10;
        n3[i]%=10;
    }
    
    if(n3[l_max])
        l_max++;
    while(n3[l_max-1]>=10)
    {
        n3[l_max]=n3[l_max-1]/10;
        n3[l_max-1]%=10;
        l_max++;
    }
    while(n3[l_max-1]==0&&l_max>1)
        l_max--;

    for(int i=l_max-1;i>=0;i--)
    {
        cout << n3[i];
    }
}

int main()
{
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cin >> s1 >> s2;
    mul();
    return 0;
}

2.高精度乘高精度FFT最佳化演算法

思想:將兩個大整數看作多項式的係數,然後利用FFT演算法在$O(n log n)$的時間複雜度內計算出它們的乘積,並最終得到乘積的各位數字。

複雜度:$o(nlog(n))$

具體步驟:

  1. 將輸入的兩個大整數轉換為對應的多項式表示,其中每個數字位作為多項式的係數。
  2. 對兩個多項式進行補零操作,使其長度變為2的冪,方便進行FFT運算。
  3. 利用FFT演算法對兩個多項式進行快速傅立葉變換,得到它們在頻域上的表示。
  4. 將兩個多項式在頻域上進行點乘操作,得到它們的乘積在頻域上的表示。
  5. 對乘積進行逆FFT操作,得到乘積多項式在時域上的表示。
  6. 對乘積多項式進行進位處理,並輸出最終的乘積結果。
// By SnowDream
#include<bits/stdc++.h>
using namespace std;

#define L(x) (1<< (x))
typedef long long ll;
const int N=1e5+10;

const double PI = acos(-1.0);
string s1,s2;
double ax[N],ay[N],bx[N],by[N];
int n1[N],n2[N],n3[N];

//將一個整數x的二進位制表示進行位逆序排列,並返回結果。
// 函式接受兩個引數:整數x和表示位數的整數bits
int reverseBits(int x, int bits)
{
    int ret = 0;//儲存位逆序排列後的結果,初始化為0
    for(int i=0;i<bits;++i)//迴圈bits次,對x的二進位制表示進行位逆序排列
    {
        ret <<= 1;//將ret左移一位,為下一位的值騰出位置
        ret |=x&1;//將x的最低位與ret進行或操作,將x的最低位值新增到ret的最低位
        x >>=1;//將x右移一位,準備處理下一位
    }
    return ret;
}

void fft(double * a, double * b, int n, bool rev)
{
    int bits = 0;// 計算n的二進位制表示中位數
    while (1 << bits < n) ++bits;// 找到大於n的最小的2的冪
    for (int i = 0; i < n; i++)
    {
        int j = reverseBits(i, bits);// 將i按照bits位反轉後的值賦給j
        if (i < j)// 交換a和b陣列中i和j位置的值
        {
            swap(a[i], a[j]);
            swap(b[i], b[j]);
        }
    }
    for (int len = 2; len <= n; len <<= 1)// 迭代每個長度為len的子序列
    {
        int half = len >> 1;// 子序列長度的一半
        double wmx = cos(2 * PI / len), wmy = sin(2 * PI / len);// 計算旋轉因子的實部和虛部
        if (rev) wmy = -wmy;// 如果是逆變換,則虛部取相反數
        for (int i = 0; i < n; i += len)// 遍歷每個子序列的起始位置
        {
            double wx = 1, wy = 0;// 初始化旋轉因子
            for (int j = 0; j < half; j++)// 遍歷子序列的前半部分
            {
                double cx = a[i + j], cy = b[i + j];// 獲取當前位置的實部和虛部
                double dx = a[i + j + half], dy = b[i + j + half];// 獲取對應的另一半的實部和虛部
                double ex = dx * wx - dy * wy, ey = dx * wy + dy * wx;// 計算旋轉後的值
                a[i + j] = cx + ex, b[i + j] = cy + ey;// 更新前半部分的值
                a[i + j + half] = cx - ex, b[i + j + half] = cy - ey;// 更新後半部分的值
                double wnx = wx * wmx - wy * wmy, wny = wx * wmy + wy * wmx;// 計算下一個旋轉因子
                wx = wnx, wy = wny;// 更新旋轉因子
            }
        }
    }
    if (rev)
    {
        for (int i = 0; i < n; i++)
            a[i] /= n, b[i] /= n;// 對結果進行歸一化
    }
}

int solve(int l1,int l2,int ans[])
{
    int len = max(l1, l2), ln;
    for(ln=0; L(ln)<len; ++ln);// 找到大於len的最小的2的冪
    len=L(++ln);// 計算2的冪
    for (int i = 0; i < len ; ++i)
    {
        if (i >= l1) ax[i] = 0, ay[i] =0;// 如果i超過l1,則ax[i]和ay[i]賦值為0
        else ax[i] = n1[i], ay[i] = 0;// 否則將n1[i]賦值給ax[i],ay[i]賦值為0
    }
    fft(ax, ay, len, false);// 進行快速傅立葉變換
    for (int i = 0; i < len; ++i)
    {
        if (i >= l2) bx[i] = 0, by[i] = 0;// 如果i超過l2,則bx[i]和by[i]賦值為0
        else bx[i] = n2[i], by[i] = 0;// 否則將n2[i]賦值給bx[i],by[i]賦值為0
    }
    fft(bx, by, len, false);// 進行快速傅立葉變換
    for (int i = 0; i < len; ++i)
    {
        double cx = ax[i] * bx[i] - ay[i] * by[i];// 計算乘積的實部
        double cy = ax[i] * by[i] + ay[i] * bx[i];// 計算乘積的虛部
        ax[i] = cx, ay[i] = cy;// 更新ax和ay陣列的值
    }
    fft(ax, ay, len, true);// 進行逆快速傅立葉變換
    for (int i = 0; i < len; ++i)
        ans[i] = (int)(ax[i] + 0.5);// 四捨五入取整
    return len;
}

void mul()
{
    int l;
    int i;
    string ans;
    memset(n3, 0, sizeof(n3));
    int l1 = (int)s1.size();
    int l2 = (int)s2.size();
    for(i = 0; i < l1; i++)
        n1[i] = s1[l1 - i - 1]-'0';
    for(i = 0; i < l2; i++)
        n2[i] = s2[l2-i-1]-'0';
    l = solve(l1,l2, n3);
    for(i = 0; i<l || n3[i] >= 10; i++) // 進位
    {
        n3[i + 1] += n3[i] / 10;
        n3[i] %= 10;
    }
    l = i;
    while(n3[l] <= 0 && l>0) l--; // 檢索最高位
    for(i = l; i >= 0; i--)
        cout << n3[i]; // 倒序輸出
}

int main()
{
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cin >> s1 >> s2;
    mul();
    cout << "\n";
    return 0;
}

3.高精度乘單精度

思想:倒置相乘,統一處理進位,再還原。

複雜度:$o(n)$

// By SnowDream
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1e5+10;
string s1;
int s2;
int n1[N];
void mul()
{
    string ans;
    int l1 = (int)s1.size();
    for(int i=0;i<l1;++i)
    {
        n1[l1-1-i]=s1[i]-'0';
    }
    int w=0;
    for(int i=0;i<l1;++i)
    {
        n1[i]=n1[i]*s2+w;
        w=n1[i]/10;
        n1[i]=n1[i]%10;
    }
    while(w)
    {
        n1[l1++]=w%10;
        w/=10;
    }
    for(int i=l1-1;i>=0;i--)
        cout << n1[i];
}

int main()
{
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cin >> s1 >> s2;
    mul();
    cout << "\n";
    return 0;
}

相關文章