高精度アルゴリズム(C/C++)

高精度アルゴリズム(C/C++)
ACMの問題をする時、いつも大きい数の加減乗除に出会って、べき乗、階乗の計算、この時与えられたデータのタイプは往々にして最後の結果を表すことができなくて、この時高精度のアルゴリズムを使う必要があります.高精度アルゴリズムの本質は,大数をいくつかの固定長のブロックに分解し,各ブロックに対して対応する演算を行うことである.ここでは、4桁の数字を1ブロックとして考える例として、入力の大数は正の整数(他のビットも考えられるが、各ブロックで対応する演算を行う場合はデータ型の数値範囲を超えてはいけないことに注意し、負の整数があれば読み込み時に正負の番号を判断して演算を決定する)とする.
1.高精度加算
3479957928375817+897259321542445を例に挙げます.
3479
9579
2837
5817
+897
+2593
+2154
+4245
=
=
=
=
4376
12172
4991
10062
キャリー0
キャリー1
キャリー0
キャリー1
4377
2172
4992
0062
C言語実装コードは以下の通りである.

#include 
#include 
#include 
#define N 200

//        
int Pow(int a, int b)
{
    int i = 0, result = 1;
    for(i = 0; i < b; ++i)
    {
        result *= a;
    }
    return result;
}


//High Precision Of Addition
int main()
{
    char stra[N], strb[N];      //     ，           ；
    int i = 0, step = 4, carry = 0; //step    ，carry    ；
    int lengtha, lengthb, maxlength, resultsize;    //maxlength  stra strb         ；
    int numa[N], numb[N],numc[N];   //       ，  ， ；
    memset(numa, 0, sizeof(numa));
    memset(numb, 0, sizeof(numb));
    memset(numc, 0, sizeof(numc));         //     ；
    scanf("%s%s", stra, strb);
    lengtha = strlen(stra);
    lengthb = strlen(strb);     //         
    //               
    for(i = lengtha-1; i >= 0; --i)
    {
        numa[(lengtha-1-i)/step] += (stra[i]-'0')*Pow(10,(lengtha-1-i)%step);
    }
    for(i = lengthb-1; i >= 0; --i)
    {
        numb[(lengthb-1-i)/step] += (strb[i]-'0')*Pow(10,(lengthb-1-i)%step);
    }
    maxlength = lengtha > lengthb ? lengtha : lengthb;

    //    ，   
    for(i = 0; i <= maxlength/step; ++i)
    {
        numc[i] = (numa[i] + numb[i])%Pow(10, step) + carry;    //   
        carry = (numa[i] + numb[i])/Pow(10, step);  //    
    }

    //          
    resultsize = numc[maxlength/step] > 0 ? maxlength/step : maxlength/step - 1;
    printf("%d", numc[resultsize]);
    for(i = resultsize-1; i >= 0; --i)
    {
        printf("%04d", numc[i]);    //   ，    ；
    }
    printf("
");
    return 0;
}

2.高精度減算
加算と同様に、プラスマイナス記号や表示桁数の変化に注意します.9999037289799-100044215000を例にとると、まず被減数と減数のどちらが大きいかを判断し、明らかにここで減数が大きいため、出力結果は負である.小数を大数で減算し、(あるブロックが負に減算された場合、上位ブロックにビットを借りる)次の表
100
0046
4201
5000
-99
-9990
-3728
-9799
1
56
473
5201
借方0
借方1
借方0
借方1
0
56
472
5201
C言語実装コードは以下の通りである.

#include 
#include 
#include 
#define N 200

//        
int Pow(int a, int b)
{
    int i = 0, result = 1;
    for(i = 0; i < b; ++i)
    {
        result *= a;
    }
    return result;
}

//High Precision Of Subtraction
int main()
{
    char stra[N], strb[N];     //     ，           ；
    int i = 0, step = 4, borrow = 0, mark = 0; //step    ，borrow    , mark      ；
    int lengtha, lengthb, maxlength, resultsize;    //maxlength  stra strb         ；
    int numa[N], numb[N],numc[N],  *maxnum, *minnum;   //       ，  ， ；
    memset(stra, 0, sizeof(stra));
    memset(strb, 0, sizeof(strb));
    memset(numa, 0, sizeof(numa));
    memset(numb, 0, sizeof(numb));
    memset(numc, 0, sizeof(numc));         //     ；
    scanf("%s%s", stra, strb);
    lengtha = strlen(stra);
    lengthb = strlen(strb);     //         
    maxlength = lengtha >= lengthb ? lengtha : lengthb;

    //               
    for(i = lengtha-1; i >= 0; --i)
    {
        numa[(lengtha-1-i)/step] += (stra[i]-'0')*Pow(10,(lengtha-1-i)%step);
    }
    for(i = lengthb-1; i >= 0; --i)
    {
        numb[(lengthb-1-i)/step] += (strb[i]-'0')*Pow(10,(lengthb-1-i)%step);
    }

    //      
    maxnum = numa;
    minnum = numb;
    mark = 1;
    for(i = (maxlength-1)/step; i >= 0; --i)
    {
        if(numa[i] > numb[i])
        {
            maxnum = numa;
            minnum = numb;
            mark = 1;
            break;
        }
        else if(numa[i] < numb[i])
        {
            maxnum = numb;
            minnum = numa;
            mark = -1;
            break;
        }
    }

    //    ，   
    for(i = 0; i <= maxlength/step; ++i)
    {
        numc[i] = (maxnum[i] - minnum[i] + Pow(10, step) + borrow)%Pow(10,step);    //   
        borrow = (maxnum[i] - minnum[i] + Pow(10, step) + borrow)/Pow(10, step) - 1;  //    
    }

    //          
    resultsize = maxlength/step;
    while(!numc[resultsize])    --resultsize;
    printf("%d", mark*numc[resultsize]);
    for(i = resultsize-1; i >= 0; --i)
    {
        printf("%04d", numc[i]);    //   ，    ；
    }
    printf("
");
    return 0;
}

3.高精度乗算
乗算は、4296556241 x 56241を例に、乗数の各ビットが被乗数に乗算された後に加算されると見なすことができる.
ふとう乗数
42
9655
6241
乗数
5
6
2
4
1
被乗数x乗数
42
9655
6241
1
42
9655
6241
4
168*10
38620*10
24964*10
2
84*100
19310*100
12482*100
6
252*1000
57930*1000
37446*1000
5
210*10000
48275*10000
31205*10000
累加和
2362122
543006855
351000081
けたあげ
241
54304
35100
積
241
6426
1955
0081
C言語実装コードは以下の通りである.

#include 
#include 
#include 
#define N 200

//        
int Pow(int a, int b)
{
    int i = 0, result = 1;
    for(i = 0; i < b; ++i)
    {
        result *= a;
    }
    return result;
}


//High Precision Of Multiplication
int main()
{
    char stra[N], strb[N];      //     ，           ；
    int i = 0, j = 0, k = 0, step = 4, carry = 0; //step    ，carry    ；
    int lengtha, lengthb, resultsize, tmpsize, eachnum;  //resultsize      ，eachnum          
    int numa[N], numb[N], numc[N], tmp[N];   //        & ，  ；
    memset(numa, 0, sizeof(numa));
    memset(numb, 0, sizeof(numb));
    memset(numc, 0, sizeof(numc));  //     ；
    scanf("%s%s", stra, strb);
    lengtha = strlen(stra);
    lengthb = strlen(strb);     //         
    //                   
    for(i = lengtha-1; i >= 0; --i)
    {
        numa[(lengtha-1-i)/step] += (stra[i]-'0')*Pow(10,(lengtha-1-i)%step);
    }
    //                    
    for(i = lengthb-1; i >= 0; --i)
    {
        numb[lengthb-1-i] = strb[i]-'0';
    }

    resultsize = tmpsize = (lengtha-1)/step;
    //                ，   ；
    for(i = 0; i < lengthb; ++i)
    {
        memcpy(tmp, numa, sizeof(numa));    // numa     tmp  ；

        k = i/step;     //k                ；
        if(k)
        {
            for(j = tmpsize; j >= 0; --j)
            {
                tmp[j+k] = tmp[j];
                tmp[j] = 0;
            }
            tmpsize += k;
        }

        //            ；
        eachnum = numb[i]*Pow(10, i%step);
        for(j = 0; j <= tmpsize; ++j)
        {
            tmp[j] *= eachnum;
        }

        //    
        carry = 0;  //    ；
        for(j = 0; j <= resultsize; ++j)
        {
            numc[j] += tmp[j] + carry;
            carry = numc[j]/Pow(10,step);
            numc[j] %= Pow(10, step);
        }
        if(carry)
        {
            ++resultsize;
            numc[j] += carry;
        }
    }

    //  
    printf("%d", numc[resultsize]);
    for(i = resultsize-1; i >= 0; --i)
    {
        printf("%04d", numc[i]);    //   ，    ；
    }
    printf("
");
    return 0;
}

4.高精度除算
高精度除算には2種類あり、1つは高精度を低精度で割ったもので、もう1つは高精度を高精度で割ったものである.前者は各ブロックを低精度で除算するだけでよい.後者は,高精度減算で実現することを考慮し,すなわち,高精度除数を毎回減算し,除数より小さくなるまで減算すると,減算回数は商であり,残りは残数である.

高精度を低精度で割った例:9876342876/343:

除算される
98
7634
2876
除算
343
ローブロックへのキャリー
98
137
190
商
0
2879
4002
余剰
190
C言語コードは以下のように実現される.

#include 
#include 
#include 
#define N 200

//        
int Pow(int a, int b)
{
    int i = 0, result = 1;
    for(i = 0; i < b; ++i)
    {
        result *= a;
    }
    return result;
}

//High Precision Of division
//(1)        
int main()
{
    char stra[N];      //     ，             ；
    int i = 0, step = 4, carry = 0; //step    ，carry         ；
    int lengtha, resultsize;
    int numa[N], numb, numc[N], numd;   //       ，  ， ,   ；
    memset(numa, 0, sizeof(numa));
    memset(numc, 0, sizeof(numc));         //     ；
    scanf("%s%d", stra, &numb);
    lengtha = strlen(stra);   //        

    //               
    for(i = lengtha-1; i >= 0; --i)
    {
        numa[(lengtha-1-i)/step] += (stra[i]-'0')*Pow(10,(lengtha-1-i)%step);
    }

    carry = 0;  //          
    resultsize = (lengtha-1)/step;
    //    ，       
    for(i = resultsize; i >= 0; --i)
    {
        numc[i] = (numa[i] + carry*Pow(10,step))/numb;    //   
        carry = (numa[i] + carry*Pow(10,step))%numb;  //    
    }
    numd = carry;   //                

    //          
    while(!numc[resultsize])    --resultsize;
    //   
    printf("%d", numc[resultsize]);
    for(i = resultsize-1; i >= 0; --i)
    {
        printf("%04d", numc[i]);    //   ，    ；
    }
    //    
    printf("
%d
", numd);
    return 0;
}

高精度を高精度で割った例:176342876/3453452:

除算される
176342876
-(51 x除数)
51 x 3453452
余剰
216824
商
51
C言語コードは以下のように実現される.

#include 
#include 
#include 
#define N 200

//        
int Pow(int a, int b)
{
    int i = 0, result = 1;
    for(i = 0; i < b; ++i)
    {
        result *= a;
    }
    return result;
}

//High Precision Of division
//(2)        
int main()
{
    char stra[N], strb[N];     //     ，           ；
    int i = 0, step = 4, borrow = 0; //step    ，borrow    ；
    int lengtha, lengthb, tmpnum, numbsize, numcsize, numdsize, maxsize, mark;    //maxlength  stra strb         ；
    int numa[N], numb[N], numc[N], numd[N];   //        ，   ， ，  ；
    memset(stra, 0, sizeof(stra));
    memset(strb, 0, sizeof(strb));
    memset(numa, 0, sizeof(numa));
    memset(numb, 0, sizeof(numb));
    memset(numc, 0, sizeof(numc));
    memset(numd, 0, sizeof(numd));      //     ；
    scanf("%s%s", stra, strb);
    lengtha = strlen(stra);
    lengthb = strlen(strb);     //         

    //               
    for(i = lengtha-1; i >= 0; --i)
    {
        numa[(lengtha-1-i)/step] += (stra[i]-'0')*Pow(10,(lengtha-1-i)%step);
    }
    for(i = lengthb-1; i >= 0; --i)
    {
        numb[(lengthb-1-i)/step] += (strb[i]-'0')*Pow(10,(lengthb-1-i)%step);
    }
    memcpy(numd, numa, sizeof(numa));
    numbsize = (lengthb-1)/step;
    numcsize = 0;
    numdsize = (lengtha-1)/step;

    do
    {
        maxsize = numdsize > numbsize ? numdsize : numbsize;
        //           
        mark = 1;
        for(i = maxsize; i >= 0; --i)
        {
            if(numd[i] > numb[i])
            {
                mark = 1;
                break;
            }
            else if(numd[i] < numb[i])
            {
                mark = -1;
                break;
            }
        }

        //            
        if(!(mark+1))   break;

        borrow = 0; //    ；
        //    ，   
        for(i = 0; i <= maxsize; ++i)
        {
            tmpnum = (numd[i] - numb[i] + Pow(10, step) + borrow)%Pow(10,step);    //   
            borrow = (numd[i] - numb[i] + Pow(10, step) + borrow)/Pow(10,step) - 1;  //    
            numd[i] = tmpnum;
        }
        while(!numd[numdsize])  --numdsize;

        //      ，   ；
        borrow = 1;
        for(i = 0; i <= numcsize; ++i)
        {
            numc[i] += borrow;
            borrow = numc[i]/Pow(10,step);
            numc[i] %= Pow(10,step);
        }
        if(borrow)
        {
            ++numcsize;
            numc[i] += borrow;
        }
    }while(1);

    printf("%d", numc[numcsize]);
    for(i = numcsize-1; i >= 0; --i)
    {
        printf("%04d", numc[i]);    //   ，    ；
    }
    printf("
");
    printf("%d", numd[numdsize]);
    for(i = numdsize-1; i >= 0; --i)
    {
        printf("%04d", numd[i]);    //   ，    ；
    }
    printf("
");
    return 0;
}