聲明,這裡講的效能其實都很低,如果某個操作是 O(n^2),這裡只是盡可能將 c1 * O(n^2) + c2 之常數項盡可能壓低,對於一些真正高效的算法有空再聊。另外,大數除以大數 for beginner ,筆者會花較長篇幅做說明,所以不在這篇提起。

 高效能部份 (其實很多 "最高效" 的筆者做得有點差) 整理起來很花時間,這部份只能等筆者有空再補上。

初學者若沒概念、讀起來吃力,前兩篇可先閱過。

[大數] C 語言大數演算法 for beginner

[大數] C 語言大數演算法 for genera

  

 

[0] 一次性進位問題 - 不用

 

 

這問題筆者細思過,我們以一般加法為例,一種是逐次進位,一種是一次性調整,兩種寫法如下所示。

 

Code Snippet
  1. // 逐次進位
  2. void big_add(int * A , int * B , int * C)
  3. {
  4.     int i, j, c;
  5.     for(c=i=0; i<SIZE; ++i) {
  6.         C[i] = A[i] + B[i] + c;
  7.         if(C[i]>=10) c=1, C[i]-=10;
  8.         else c=0;
  9.     }
  10. }
  11.  
  12. // 一次進位
  13. void big_add(int * A , int * B , int * C)
  14. {
  15.     int i;
  16.     for(i=0; i<SIZE; ++i)
  17.         C[i] = A[i] + B[i];
  18.     for(i=0; i<SIZE-1; ++i) {
  19.         C[i+1] += C[i]/10;
  20.         C[i] %= 10;
  21.     }
  22. }

 

在加、減法可能比較看不出來,但在乘法的時候,一次性進位缺點會不少。在萬進位情況下,使用一次性進位容易溢位;65536進位下更容易溢位,所以筆者打算全都用逐次進位之方式。但對乘法而言,似乎不多人會寫逐次進位,所以寫起來有三個 loop,其實一樣只要兩個 loop 就行了。

 

Code Snippet
  1. void big_mul(int * A , int * B , int * C)
  2. {
  3.     int i, j, ij, carry=0;
  4.     memset( (void*)C, 0, BSIZE);
  5.     for(i=0; i<SIZE; ++i){
  6.         for(j=0; (ij=i+j) < SIZE; ++j){
  7.             C[ij] += carry + A[i] * B[j];
  8.             carry = C[ij] / 10;
  9.             C[ij] %= 10;
  10.         }
  11.     }
  12. }

 

上述逐次進位方式並不會造成任何 ov 可能性,接下來我們全都用逐次進位方式,不再一口氣進位。

 

 

[1] 改善資料結構

 

 

一般我們從字串「讀入」開始,是直接轉成大數,這其實也沒問題,但一般大數都是直接以  int arr[SIZE] 方式表達,這也合理。而一份簡單降低一點 BigO 常數的技巧是 : 紀錄用了多少位數。

ok,我不確定這是十分有效的,只是直覺它是可以有效的。為方便說明,我還是給靜態陣列大小,不做動態配置。但事實上這篇講完之後,其實要實作動態自動調整陣列大小,已不是問題,只是為避免內文冗長,所以用靜態說明。

在定義大數時候,我給了一份 struct

 

Code Snippet
  1. #define BIG_SIZE 20
  2. #define BIG_CAP  10
  3. #define BIG_LEN  1
  4.  
  5. typedef struct tagBig{
  6.     int arr[BIG_SIZE]; // 大數陣列
  7.     int dig;           // 使用位數
  8. }Big;

 

dig 代表這個陣列裡面實際使用的位數,試想一下,如果是 00001 + 000002 ,實際上是只要加一個位數就好 (他們兩個的 dig 都是 1),但一般傳統給的是把整個 BIG_SIZE 都算完。換而言之,傳統之方式, BIG_SIZE 給的愈大,計算愈久,但表示範圍廣;BIG_SIZE 給的愈小,計算時間不會拉長,表示範圍長。

另外這份程式其實可以再加一個 int cur_size,做動態方式調整,這裡就不再敘述。而傳遞方式以 pass by pointer ( C++ 可用 pass by reference) 方式傳遞,否則會造成 struct 裡之 array 初始化時間過長問題。

接下來所有的作法都只要稍想一下就可以了。然後有四份 macro 先定義下來,方便使用。

 

Code Snippet
  1. #define MAX(a,b) ( ((a)>=(b)) ? (a) : (b) )
  2. #define MIN(a,b) ( ((a)<=(b)) ? (a) : (b) )
  3. #define BSIZE(n) ( sizeof(int) * n)
  4. #define SWAP(TYPE, a, b) { TYPE t; t=a; a=b; }

 

 

[2] 從字串讀入大數

 

 

和之前程式碼雷同,不過為了省時間,並沒有調用 memset、填零之動作。但多了「去字串前導零」動作,另外也紀錄了使用了多少位數。

Code Snippet
  1. void big_from_string(Big * big , const char * str)
  2. {
  3.     int i, j, len, pwr;
  4.     int * arr = big->arr;
  5.     // 去前導0
  6.     while(*str=='0' && *str!='\0') ++str;
  7.     // 轉入大數
  8.     len = strlen(str);
  9.     arr[0] = 0;
  10.     for(j=0, i=len-1; i>=0; ++j, arr[j]=0)
  11.         for(pwr=1; pwr!=BIG_CAP && i>=0; pwr*=10, --i)
  12.             arr[j]+=(str[i]-'0')*pwr;
  13.     // 計算位數
  14.     big->dig = j;
  15. }

 

換而言之,即使 BIG_SIZE = 100,實際使用 dig = 10,但其實不會把其他 90 個大小填入 0 ,以期節省時間。另注意件事,假設 capacity = 10 (10進位),當轉出來之數值為 0 時, dig = 0;1~9 時,dig=1;10~19 時,dig=2.. 依此類推。這麼設計理由,是認為在做以下的運算會比較方便。

 

 

[3] 輸出大數結果

 

 

若 dig = 0 時直接輸出 0 ,否則從 dig 開始輸出即可。 

Code Snippet
  1. // 輸出
  2. void big_print(Big * big)
  3. {
  4.     int i = big->dig;
  5.  
  6.     //printf("(dig = %d)", i);
  7.     if(!i) printf("0");
  8.     else {
  9.         int *arr = big->arr;
  10.         printf("%d", arr[--i]); // 輸出第一位數字
  11.         for(--i; i>=0; --i) // 輸出其他數字
  12.             printf("%0*d", BIG_LEN, arr[i]);
  13.     }
  14. }

 

這裡給一份測試碼

 

Code Snippet
  1. ////
  2. int main()
  3. {
  4.     char s1[2000], s2[2000];
  5.     Big A, B;    
  6.     while(scanf("%s%s", s1,s2)==2){
  7.  
  8.         big_from_string(&A, s1); printf("\nA = "), big_print(&A);
  9.         big_from_string(&B, s2); printf("\nB = "), big_print(&B);
  10.     }
  11.     return 0;
  12. }

加上 stdio.h、stdlib.h、string.h,上面四段合起來可做一份測試了。

 

 

[4] 大數加法

 

 

(1) 結果位數 : 結果是 Max(A的位數 , B 的位數) + 1,最後的 +1 不一定會有,如果最後還有進位的話就有 +1 。

(2) 偷用 memcpy : 試想一下 123456 + 14 的時候,其實真正在計算的只有 3 位數,56+14,還有一個 carry + 4,其他的是用 memcpy 直接照抄。

(3) 特殊值判斷 : +0 的話就直接 memcpy 。

 

Code Snippet
  1. // C = A + B
  2. void big_add(Big * BigA , Big * BigB , Big * BigC)
  3. {
  4.     int i, carry;
  5.     int len_a = BigA->dig, len_b = BigB->dig;
  6.     int *A = BigA->arr, * B = BigB->arr;
  7.     int *C = BigC->arr;
  8.  
  9.     // 特殊值判斷
  10.     if(len_a==0 && len_b==0) { // A = B = 0
  11.         BigC->dig = 0;
  12.         return ;
  13.     } else if(len_a==0) { // A = 0, C = A + B = B
  14.         memcpy((void*)C, (void*)B, BSIZE(len_b));
  15.         BigC->dig = len_b;
  16.         return ;
  17.     } else if(len_b==0) { // B = 0, C = A + B = A
  18.         memcpy((void*)C, (void*)A, BSIZE(len_a));
  19.         BigC->dig = len_a;
  20.         return ;
  21.     }
  22.  
  23.     // 保證 len_a >= len_b
  24.     if(len_a < len_b) {
  25.         SWAP(int , len_a, len_b);
  26.         SWAP(int * , A, B);
  27.     }
  28.     A[len_a]=0, ++len_a;
  29.  
  30.     // 加到 len_b
  31.     for(carry=i=0; i<len_b; ++i) {
  32.         C[i]  = A[i] + B[i] + carry;
  33.         carry = C[i]/BIG_CAP;
  34.         C[i]%=BIG_CAP;
  35.     }
  36.  
  37.     // 若有進位,將連續之 BIG_CAP-1 全改 0
  38.     if(carry)
  39.         while(A[i]==BIG_CAP2)
  40.             C[i]=0, ++i;
  41.     C[i] = A[i] + carry, ++i; // 補上最後一個進位
  42.  
  43.     // 剩下之 A 全丟到 C
  44.     if(i<len_a)
  45.         memcpy( (void*)(C+i), (void*)(A+i), BSIZE(len_a-i));
  46.  
  47.     // 調整 len_a, 寫回結果
  48.     BigC->dig = len_a - !C[len_a-1];
  49. }

 

 

[5] 大數減法

 

 

(1) 結果位數 : 結果是 0~Max(A的位數 , B 的位數),這裡最後是用 polling 方式去查。

(2) 偷用 memcpy : 和上面一樣

(3) 特殊值判斷 : -0 的話就直接 memcpy 。

(4) 我們不考慮結果為負的時候之處理。


Code Snippet
  1. // C = A - B , 逐次借位
  2. void big_sub(Big * BigA , Big * BigB , Big * BigC)
  3. {
  4.     int i, borrow;
  5.     int len_a = BigA->dig, len_b = BigB->dig;
  6.     int *A = BigA->arr, * B = BigB->arr;
  7.     int *C = BigC->arr;
  8.  
  9.     // 特殊值判斷
  10.     if(len_a==0 && len_b==0) { // A = B = 0
  11.         BigC->dig = 0;
  12.         return ;
  13.     } else if(len_a==0 || len_b > len_a) { // A = 0, C = A - B = -B
  14.         // 直接設 0 不處理 < 或以十補數方式處理 >
  15.         BigC->dig = 0;
  16.         return ;
  17.     } else if(len_b==0) { // B = 0, C = A - B = A
  18.         memcpy((void*)C, (void*)A, BSIZE(len_a));
  19.         BigC->dig = len_a;
  20.         return ;
  21.     }
  22.  
  23.     // 先做 len_b 個減法
  24.     for(borrow=i=0; i<len_b; ++i) {
  25.         C[i] = A[i] - B[i] - borrow;
  26.         if(C[i] < 0) C[i]+=BIG_CAP, borrow=1;
  27.         else borrow=0;
  28.     }
  29.  
  30.     // A < B , 溢位處理
  31.     if(len_a<=len_b && borrow){
  32.         BigC->dig = 0;
  33.         return ;
  34.     }
  35.  
  36.     // 把借位處理完, 遇到 0 直接給 BIG_CAP-1
  37.     if(borrow)
  38.         while(A[i]==0)
  39.             C[i]=BIG_CAP2, ++i;
  40.     C[i]=A[i]-borrow, ++i; // 補上最後借位
  41.  
  42.     // 剩下的 A 全給 C
  43.     if(i < len_a)
  44.         memcpy( (void*)(C+i), (void*)(A+i), BSIZE(len_a-i));
  45.  
  46.     //回頭查剩幾位, 寫回 C
  47.     for(i=len_a-1 ; i>0 && C[i]==0 ; --i);
  48.     if(i==0) BigC->dig = (C[0]!=0);
  49.     else BigC->dig = i+1;
  50. }

 

 

[6] 大數乘法

 

 

大數乘法這是效率低的版本 O(n^2),另外有 O(n^log(3)) 、O(nlogn) 版本。

 

(1) 結果位數 : 正確結果應該是 (A位數 + B位數 -1 ) ~ (A 位數 + B 位數 + 1),那個 -1 蠻多人略過、覺得不對,舉個例,2 * 4 = 8,就是有 -1 的情況。

(2) 這裡不能用 memcpy 加速。

(3) 特殊值判斷 : 有 0 設 0 ,有 1 設 A/B。

(4) 這裡沒處理位數不夠放之情況,但註解裡會提出要注意。

(5) 拿掉 carry 變數。


Code Snippet
  1. //////////////////////////////////////////////////////////////////////////
  2. // 大數乘大數, 逐次進位, 無 ov
  3.  
  4. void big_mul(Big * BigA , Big * BigB , Big * BigC)
  5. {
  6.     int i, j,ij;
  7.     int len_a = BigA->dig, len_b = BigB->dig;
  8.     int *A = BigA->arr, *B = BigB->arr, *C=BigC->arr;
  9.     int len_c = len_a + len_b + 1;
  10.  
  11.     // 特殊值判斷
  12.     if(!len_a || !len_b) { // A==0 , B==0, C = A * B = 0
  13.         BigC->dig = 0;
  14.         return ;
  15.     } else if(len_a==1 && A[0]==1) { // A==1, C = A * B = B
  16.         memcpy( (void*)C, (void*)B , BSIZE(len_b));
  17.         BigC->dig = len_b;
  18.         return ;
  19.     } else if(len_b==1 && B[0]==1) { // B==1, C = A * 1 = A
  20.         memcpy( (void*)C, (void*)A , BSIZE(len_a));
  21.         BigC->dig = len_a;
  22.         return ;
  23.     }
  24.     // 做 arr-size , len_c 之防呆
  25.  
  26.     // 開始做乘法
  27.     memset((void*)C, 0, BSIZE(len_c));
  28.     for(i = 0; i < len_a ; ++i) {
  29.         for( j=0 ; j<len_b ; ++j) {
  30.             ij = i+j;
  31.             C[ij] += A[i] * B[j];
  32.             C[ij+1] += C[ij]/BIG_CAP;            
  33.             C[ij]%=BIG_CAP;
  34.         }
  35.     }
  36.     // 修正 len_c
  37.     if(C[len_c-1]==0) --len_c;
  38.     if(C[len_c-1]==0) --len_c;
  39.     BigC->dig = len_c;
  40. }


 

[7] 大數乘整數


大數乘整數習慣上還是用一次性進位。

 

Code Snippet
  1. void big_mul_int(Big * BigA , int numB , Big * BigC)
  2. {
  3.     int i, carry;
  4.     int *A = BigA->arr, * C = BigC->arr;
  5.     int len_a = BigA->dig;
  6.  
  7.     // 特殊值判斷
  8.     if( !len_a || !numB) { // 設0
  9.         BigC->dig = 0;
  10.         return ;
  11.     } /* 其他 A, numB = 1 之判斷 */
  12.  
  13.     // 開始做乘法
  14.     for( carry = i = 0; i<len_a; ++i){
  15.         C[i] = A[i] * numB + carry; // (!!!) 注意這行
  16.         carry = C[i] / BIG_CAP;
  17.         C[i] %= BIG_CAP;
  18.     }
  19.  
  20.     // 再調整進位
  21.     while(carry) {
  22.         C[i] = carry % BIG_CAP;
  23.         carry /= BIG_CAP;
  24.         ++i;
  25.     }
  26.     BigC->dig = i;
  27. }

 

注意到上面的算法其實會有溢位的可能性,假設為 10 進位時,當 numB = 2147483647 就溢位了,會溢位就不叫大數運算了。但相對的,這種算法在針對 numB 是小量的,比如說在 BIG_CAP 以內,且 numB 也在 BIG_CAP 內的話,其實是不會有溢位問題的。

上面這段碼是拿來做大數除大數用的乘法,因在過程中的設計的確用這段碼就行了。要考慮不溢位的話其實很簡單,再將 numB 轉成一個小量的大數 (BIG_CAP = 10 的時候了不起陣列大小也才 10 而已) 即可。

 

Code Snippet
  1. void big_mul_int(Big * BigA , int numB , Big * BigC)
  2. {
  3.     int i, j, ij, carry;
  4.     int *A = BigA->arr, * C = BigC->arr;
  5.     int len_a = BigA->dig;
  6.     int len_b, len_c;
  7.     int B[20];
  8.  
  9.     // 特殊值判斷
  10.     if( !len_a || !numB) { // 設0
  11.         BigC->dig = 0;
  12.         return ;
  13.     } /* 其他 A, numB = 1 之判斷 */
  14.  
  15.     // 將整數 numB 轉到 array 裡去
  16.     for(len_b=0; numB; ++len_b, numB/=BIG_CAP)
  17.         B[len_b] = numB%BIG_CAP;
  18.     len_c = len_a + len_b +1;
  19.  
  20.     // 開始做乘法
  21.     memset((void*)C, 0, BSIZE(len_c));
  22.     for(i = 0; i < len_a ; ++i) {
  23.         for(carry = j=0 ; j<len_b ; ++j) {
  24.             ij = i+j;
  25.             C[ij] += A[i] * B[j];
  26.             C[ij+1] += C[ij]/BIG_CAP;            
  27.             C[ij]%=BIG_CAP;
  28.         }
  29.     }
  30.     // 修正 len_c
  31.     if(C[len_c-1]==0) --len_c;
  32.     if(C[len_c-1]==0) --len_c;
  33.     BigC->dig = len_c;    
  34. }

 

[7] 大數除整數

 

實際上考慮直式除法便行。

 

Code Snippet
  1. /////////////////////////////////////////////////
  2.  
  3. //
  4. // 大數除整數, 無溢位問題
  5. // BigA : 被除數, numB : 除數 , C : 商
  6. // 傳回 : 成功傳回餘數(>=0), 失敗傳回負值
  7. //
  8. int big_div_int(Big * BigA , int numB, Big * BigC)
  9. {
  10.     int * A = BigA->arr, * C = BigC->arr;
  11.     int i, len_a = BigA->dig;
  12.     int r, len_c; //餘數
  13.  
  14.     // 特殊值判斷
  15.     if(numB<=0) { // 除零錯誤, 除到負數
  16.         BigC->dig=0;
  17.         return -1;
  18.     } else if(!len_a) {// 分子為 0
  19.         BigC->dig = 0;
  20.         return 0; // 餘數設 0
  21.     } else if(numB==1) { // 分母為 1
  22.         memcpy( (void*)BigC, (void*)BigA, len_a);
  23.         BigC->dig = len_a;
  24.         return 0; // 餘數為 0
  25.     }
  26.  
  27.     // 開始做除法
  28.     r=0;
  29.     for(i=len_a-1; i>=0; --i) {
  30.         r = r *BIG_CAP + A[i];
  31.         C[i] = r / numB;
  32.         r %= numB;
  33.     }
  34.  
  35.     for(i=len_a-1; i>0 && C[i]==0; --i);
  36.  
  37.     BigC->dig = i + !(i==0 && C[0]==0);
  38.     //     if(i==0 && C[0]==0) BigC->dig = 0;
  39.     //     else BigC->dig=i+1;
  40.  
  41.     return r; // 傳回餘數
  42. }

 

再來是大數除大數,不過這該換一篇了...

arrow
arrow
    全站熱搜

    edisonx 發表在 痞客邦 留言(0) 人氣()