close

 

結構體與部份函式 續上篇。

 

[0] 前言

 

 

(1) 這篇提的大數除大數效能也還不是最好的，但應算「可用」。之前看過的一種作法是：被除數一直去減除數，商慢慢加1。

(2) 這裡提的方法其實可以實作 3  種出來，但最好的我沒時間做出來。

 

 

[1] 直式除法說明

 

 

首先我們觀察位數關係。


９９９９／１０　＝　９９９．．．９    ---->  4 位數 / 2 位數 = 3 位數

１０００／９９　＝　１０　．．．１０ ---->   4 位數 / 2 位數 = 2 位數

 

所以我們可以得到的結果是： A / B ，得到的結果，最多是 (A 位數 - B 位數 + 1) 位數。

我們再想一下 7163456(被除數) / 123(除數) ，一開始的位數會怎麼填？


　　　　　　　？？？？？
　　　┌────────────
１２３│　７１６３４５６　

 

由於 123 是 3 位數，所以我們一開始也在被除數的前三位數看做一組，做試乘。但如果第一次算完的時候想一下情況會變怎樣

 


　　　　　　　５？？？？
　　　┌────────────
１２３│　７１６３４５６　
　　　　　６１５
　　　──────────────
　　　　　１０１３

　

 這時候就變成了 4 位數為一組做除法。而後面我們也不確定到底該是有幾位數做除法，簡單一點的方法是，我們都假設是 123(除數，3位數) 的位數 + 1，統一一次以 4 位 數為一組做試乘。如此一來，就必須在被除數(7163456)　最前面補一個 0。


　　　　　　　？？？？？
　　　┌────────────
１２３│０７１６３４５６　　

 

下面才是真正開始做除法的細節。

我們令被除數 7163456 為 A，除數 123 為 B。

 

－－－－－－－－－－ （１）第一次迭代　Ａ〔７：４〕 －－－－－－－－－－－－－

 


　　　　　　　Ｎ？？？？
　　　┌───────────
１２３│０７１６３４５６　


 

（１．１）　

 


Ｎ　從１開始猜到９，先找一個數滿足 ＞＝ Ａ〔７：４〕之條件


　Ｂ　　　　Ｎ　　　　ＢＮ　　　　Ａ〔７：４〕
１２３　＊　１　＝　１２３　　＜　　０７１６
１２３　＊　２　＝　２４６　　＜　　０７１６
１２３　＊　３　＝　３６９　　＜　　０７１６
１２３　＊　４　＝　４９２　　＜　　０７１６
１２３　＊　５　＝　６１５　　＜　　０７１６
１２３　＊　６　＝　７３８　　＞　　０７１６

上面的 BN 就代表 B * N 之意思。這裡得到，

Ｎ　＝６
ＢＮ＝７３８


（１．２）


Ｎ猜到６時，ＢＮ已滿足大於等於Ａ〔７：４〕條件，故停下來。


（１．３）


ＢＮ＝７３８　＞　Ａ〔７：４〕，做回覆，

Ｎ　＝　Ｎ－１

　　＝　６－１

　　＝　５

ＢＮ＝　Ｂ　＊　Ｎ　

　　＝　０１２３　＊　５　

　　＝　０６１５


表格變如下


　　　　　　　５？？？？
　　　┌────────────
１２３│０７１６３４５６
　　　│０６１５
　　　──────────────

　　

（１．４）

 

這步是重點，注意運算對象。


Ａ〔７：４〕＝ Ａ〔７：４〕－　ＢＮ

　　　　　　＝　０７１６　－　０６１５　
　
　　　　　　＝　０１０１


表格

　　　　　　　５？？？？
　　　┌───────────
１２３│０７１６３４５６
　　　│０６１５
　　　───────────
　　　　０１０１


再整理一下


　　　　　　　５？？？？
　　　┌───────────
１２３│０１０１３４５６
　　　│
　　　

 

 

－－－－－－－－－－－ （２）第二次迭代　Ａ〔６：３〕 －－－－－－－－－

 


　　　　　　　５Ｎ？？？
　　　┌─────────
１２３│０１０１３４５６

 


（２．１）　


Ｎ　從１開始猜到９，先找一個數滿足 ＞＝ Ａ〔６：３〕之條件


　Ｂ　　　　Ｎ　　　　ＢＮ　　　　Ａ〔７：４〕
１２３　＊　１　＝　１２３　　＜　　１０１３
１２３　＊　２　＝　２４６　　＜　　１０１３
１２３　＊　３　＝　３６９　　＜　　１０１３
１２３　＊　４　＝　４９２　　＜　　１０１３
１２３　＊　５　＝　６１５　　＜　　１０１３
１２３　＊　６　＝　７３８　　＜　　１０１３
１２３　＊　７　＝　８６１　　＜　　１０１３
１２３　＊　８　＝　９８４　　＜　　１０１３
１２３　＊　９　＝　１１０７　＞　　１０１３

這裡得到，

Ｎ　＝９
ＢＮ＝１１０７


（２．２）


Ｎ猜到９時，ＢＮ已滿足大於等於Ａ〔６：３〕條件，故停下來。


（２．３）


ＢＮ＝１１０７　＞　Ａ〔６：３〕，做回覆，

Ｎ　＝　Ｎ－１

　　＝　９－１

　　＝　８

ＢＮ＝　Ｂ　＊　Ｎ　

　　＝　０１２３　＊　８

　　＝　０９８４


表格變如下

 

　　　　　　　５８？？？
　　　┌───────────
１２３│０１０１３４５６
　　　│　０９８４　　　
　　　───────────


　　

（２．４）


Ａ〔６：３〕＝ Ａ〔６：３〕－　ＢＮ

　　　　　　＝　１０１３　－　０９８４　
　
　　　　　　＝　００２９


表格


　　　　　　　５８？？？
　　　┌───────────
１２３│０１０１３４５６
　　　│　０９８４　　　
　　　───────────
　　　　　００２９


再整理一下


　　　　　　　５８？？？
　　　┌───────────
１２３│０００２９４５６
　　　│

 

－－－－－－－－－－－－ （n）第n次迭代　Ａ〔3：0〕 －－－－－－－－－－－－

....

 

[2] 重點流程 / 虛碼

 

說明示範兩個步驟應該可以知道規則是什麼，就不再畫下去。重點流程大致如下。


(0) 

 

判斷除數 B 是否為 0，是的話傳回 0 代表失敗。

 

(1) 

 

找被除數 A 之使用位數，a_len，記得在高位補一個零，所以 a_len 找到後還要再加1；

找除數 B 之使用位數，b_len；

我們還需要一個暫存器 BN，存 除數(B)*猜測值(N)，其長度為 b_len+1；

再令商 C 之長度為 c_len，c_len = a_len - b_len+1

 

(2) 猜謎遊戲

 

直接看虛碼

 

Code Snippet
  1. /* 依序猜測每位數之 商，C[c_len:0] */
  2. for ( c_len = a_len - b_len ; c_len >=0 ; ++c_len) {
  3.  
  4.     // 1. 先將 BN 清成 0
  5.  
  6.     // 2. N 從 1 猜到 9, BN = B * N,
  7.     // 找出 BN >= A[c_len:c_len+r_len] 之最小 N
  8.     // (大數*=整數)
  9.  
  10.     /* 猜謎遊戲結束 ，準備填商 C[c_len] 和調整被除數 A */
  11.  
  12.  
  13.     // 3.1.
  14.     // 如果猜完 BN 還是小於 A[c_len:c_len+r_len]，
  15.     // 就代表商是 9，填入 C[c_len] = 9，
  16.     // 直接進行 A[c_len:c_len+r_len]-=BN
  17.  
  18.  
  19.     // 3.2
  20.     // 如果找到的 BN 剛好等於 A[c_len:c_len+r_len]，
  21.     // 不需做回覆，C[c_len] = N，不要呼叫減法，
  22.     // 直接將 C[c_len] 填入
  23.  
  24.     // 3.3
  25.     // 如果找到的 BN 大於 A[c_len:c_len+r_len]，
  26.     // 需要做回覆，N=N-1，BN=B*N(大數乘整數)，
  27.     // 再用 A[c_len:c_len+r_len]-=BN
  28.     // (大數-=大數)
  29. }

 

然．後．怎麼好像有點怪怪的？上面虛碼竟然有兩個怪怪的？網路上找的、這份 blog 講的大數演算法，大多都只是長這樣

 

Code Snippet
  1. void big_add( int * A , int * B , int *C); // C = A + B
  2. void big_sub( int * A , int * B , int *C); // C = A - B
  3. void big_mul( int * A , int * B , int *C); // C = A * B

 

沒辦法指定，A 從第 c_len 位置，連續 r_len 長度，和 B 做運算耶！

是的，正常來講，為了複用性， prototype 應該還要再加一個參數 : Size。如果是一般在調用的時候就長這樣

 

Code Snippet
  1. void big_add(int *A , int *B, int *C , int Size);
  2. /* do something */
  3.  
  4. const int BigSize = 200; /* VC 之 C 編譯只能用 #define 或 enum */
  5. int BigA[BigSize], BigB[BigSize], BigC[BigSize];
  6. /* do something */
  7.  
  8. big_add(A, B, C, BigSize); /* C = A + B (1) */
  9. big_add(A, B, A, BigSize); /* A = B + A, A+=B (2) */
  10.  
  11. big_add(A+a, B+b, C, op_size); /* (3) */
  12.  
  13. /*
  14. (3) C[0:op_size-1] = A[a:a+op_size-1] + B[b:b+op_size-1]
  15. */

 

在上述程式碼裡之 (1) 是常見的沒問題；(2) 的作法是想複用 C = A + B，達成 A = A + B 的作用，這手法在 C 語言裡一些場合裡面常見到，但筆者很不喜歡這樣用，一方面不是每個函式這樣塞進去都可以正確執行(甚至不少時候會出問題)，另一方面是，其實直接再寫一個 A+=B ，裡面的程式碼會有更多優化的空間。

但相對的開發時程就拉長，這是對 C language 最嚴重的致命傷！今天要寫 A+=B，可能就代表著接寫下還要寫 A-=B, A*=B，大數 op 大數一份，大數 op 整數又一份，這樣下去沒完沒了，還是快跳 C++ 吧 ，不過 C++ 不熟的話，也只是從這個地嶽，跳到另一個 C++ 地嶽而已。

離題了，繼續講(3) 吧。

(3) 的作法對老手而言是常見，其代表意義從虛碼看比較清楚，

 

Code Snippet
  1. for(i = 0 ; i<op_size ; ++i)
  2.     C[i] = A[i+a] + B[i+b];


 

是的，沒看錯，所以如果當時在設計時，Prototype 有多一個 size，
可以省很多事。總之經過分析後，這裡決定不再改之前寫的 prototype，
而是直接再額外寫三個 function 出來。

然後筆者在 C 裡，函式本身有自我賦值(也就是賦合運算子)的，
習慣多一個 s(self)代表，C++ 就直接搞 operator +=。

 

Code Snippet
  1. // A*=numB
  2. void big_mul_int(Big * A , Big numB, int Size);
  3.  
  4. // A-=B
  5. void big_sub_parts(Big * A , Big * B , int Size);
  6.  
  7. // compare A and B
  8. int compare_parts(Big * A , Big * B , int Size);

 

當然，有興趣也可以把之前寫過的 fucntion prototype 全改過再調用。基於系列文之設計架構，以下是一份完整的大數除法範例。

 

Code Snippet
  1. #include <stdio.h>
  2. #include <stdlib.h>
  3. #include <string.h>
  4.  
  5. #define BIG_SIZE 100
  6. #define BIG_CAP  10     /* 10,100,1000,100000   */
  7. #define BIG_LEN  1      /* 10^BIG_LEN = BIG_CAP */
  8. #define BIG_CAP2 9      /* BIG_CAP - 1          */
  9.  
  10. #define BSIZE(n) ( sizeof(int) * (n))
  11.  
  12. typedef struct tagBig{
  13.     int arr[BIG_SIZE]; // 大數陣列
  14.     int dig;           // 使用位數
  15. }Big;
  16.  
  17. // 輸入
  18. void big_from_string(Big * big , const char * str)
  19. {
  20.     int i, j, len, pwr;
  21.     int * arr = big->arr;
  22.     // 去前導0
  23.     while(*str=='0' && *str!='\0') ++str;
  24.     // 轉入大數
  25.     len = strlen(str);
  26.     arr[0] = 0;
  27.     for(j=0, i=len-1; i>=0; ++j, arr[j]=0)
  28.         for(pwr=1; pwr!=BIG_CAP && i>=0; pwr*=10, --i)
  29.             arr[j]+=(str[i]-'0')*pwr;
  30.     // 計算位數
  31.     big->dig = j;
  32. }
  33.  
  34. // 輸出
  35. void big_print(Big * big)
  36. {
  37.     int i = big->dig;
  38.  
  39.     if(!i) printf("0");
  40.     else {
  41.         int *arr = big->arr;
  42.         printf("%d", arr[--i]); // 輸出第一位數字
  43.         for(--i; i>=0; --i) // 輸出其他數字
  44.             printf("%0*d", BIG_LEN, arr[i]);
  45.     }
  46. }
  47.  
  48. // 除法用，大數比較
  49. int big_compare_part(int * A , int * B , int size)
  50. {
  51.     int i;
  52.     for(i=size-1; i && A[i]==B[i]; --i);
  53.     return A[i] - B[i];
  54. }
  55.  
  56. // 除法用，大數*整數，不需考慮一次進位之 ov 問題
  57. void big_mul_int_part(int * A , int numB, int * C, int size)
  58. {
  59.     int i;
  60.     for(i = 0 ; i<size ; ++i)
  61.         C[i] = numB * A[i];
  62.     for(i = 0 ; i<size-1; ++i) {
  63.         C[i+1] += C[i] / BIG_CAP;
  64.         C[i]%=BIG_CAP;
  65.     }
  66. }
  67.  
  68. // 除法用，大數減法，A-=B
  69. void big_subs(int * A , int * B , int size)
  70. {
  71.     int i, borrow=0;
  72.     for(i=0; i<size; ++i){
  73.         A[i]-=(B[i]+borrow);
  74.         if(A[i] < 0) borrow=1, A[i]+=BIG_CAP;
  75.         else borrow = 0;
  76.     }
  77. }
  78.  
  79. //
  80. // 大數除大數
  81. // BigA : 被除數，執行完後變餘數
  82. // BigB : 除數
  83. // BigC : 餘數
  84. // return : 成功(1), 失敗(0)
  85. int big_div(Big * BigA , Big * BigB , Big * BigC)
  86. {
  87.     int i, N, cmp;
  88.     int * NB; // 存放 N * B
  89.     int * A = BigA->arr, * B = BigB->arr;
  90.     int * C = BigC->arr;
  91.  
  92.     // 算長度
  93.     int a_len = BigA->dig ; // 被除數位數, 高位補0, 故要+1
  94.     int b_len = BigB->dig ; // 除數位數
  95.     int c_len = a_len - b_len + 1; // 商數位數
  96.     int r_len = b_len + 1;  // 暫存位數
  97.  
  98.     if(b_len==0) {  // 除零錯誤
  99.         BigC->dig=0;
  100.         return 0;
  101.     } else if(b_len > a_len) { // B > A, 0
  102.         BigC->dig = 0;
  103.         return 1;
  104.     }
  105.  
  106.     NB = (int*)malloc(sizeof(int) * r_len);
  107.     BigC->dig = c_len, A[a_len++]=0; // A 前面補0
  108.  
  109.     /* 尋找 C[c_len:0] 之猜商遊戲 */
  110.     for( ; c_len >=0 ; --c_len){
  111.         for(N=1; N<BIG_CAP ; ++N) {
  112.             big_mul_int_part(B, N, NB, r_len);
  113.             cmp = big_compare_part(A+c_len, NB, r_len);
  114.             if(cmp<=0) break; // A >= NB 時結束
  115.         }
  116.         /* 依比較結果 cmp 調整 */
  117.  
  118.         if(cmp==0){ // 相等, 將商填入, A[...] 填入
  119.             memset( (void*)(A+c_len), 0, BSIZE(r_len));
  120.             C[c_len] = N;
  121.         } else if(cmp < 0) { // 猜太大, 變小, 回覆
  122.             C[c_len ] = --N;
  123.             big_mul_int_part(B, N, NB, r_len); // NB = N * B
  124.             big_subs(A + c_len, NB, r_len); // A -=NB
  125.         } else { // 唯一可能剛好情況, N=BIG_CAP, 不需再做回覆
  126.             big_subs(A+c_len, NB, r_len); // A-= NB
  127.             C[c_len] = BIG_CAP2;
  128.         }
  129.     }
  130.     free( (void*)NB);
  131.     // 更新 c_len
  132.     for(i=BigC->dig; i>0 && C[i]==0; --i);
  133.     if(i==0) BigC->dig = ( (C[i]==0) ? 0 : 1 );
  134.     else BigC->dig = i+1;
  135.  
  136.     // 更新 a_len, a_len < b_len
  137.     for(i=b_len; i>0 && A[i]==0 ; --i);
  138.     if(i==0) BigA->dig = ( (A[0]==0) ? 0 : 1);
  139.     else BigA->dig = i+1;
  140.  
  141.     return 1; 
  142. }
  143.  
  144. /////////////////////////////////////////
  145. int main()
  146. {
  147.     char * s1[] = {"7123456", "0", "7428571",
  148.         "999999999999999999999999999999999999999999999999999999999999"};
  149.     char * s2[] = {"123", "120", "7123745",
  150.         "1"};
  151.     Big A, B, C, R;
  152.     int i;
  153.     ///
  154.  
  155.     for(i=0; i<sizeof(s1) / sizeof(s1[0]) ; ++i){
  156.         big_from_string(&A, s1[i]);
  157.         big_from_string(&B, s2[i]);
  158.  
  159.         R = A ; // 除法會變更除數成餘數，故將除數先備份
  160.         big_print(&A), printf(" / "), big_print(&B),  printf(" = ");
  161.         big_div(&R, &B, &C);
  162.         big_print(&C), printf(" ... ");
  163.         big_print(&R), puts("\n");
  164.     }
  165.     getchar();
  166.     return 0;
  167. }

 

 

[3] 猜商策略

 

好了，上面的 code 其實效率很糟，今天是 10 進位所以最多算 10 次，若是 萬進位、65536 進位，最多不就要算 65536 次了嗎？是的，有幾種改善方式。一種是整個把大數 * 整數拿掉，直接以加、減法方式慢慢做。

 

去掉一些拉哩拉喳的東西，大致如下。

 

Code Snippet
  1. //
  2. // 除法用, 大數自身加法, A+=B
  3. void big_adds(int * A , int * B , int size)
  4. {
  5.     int i,c=0;
  6.     for(i=0; i<size; ++i){
  7.         A[i]+=B[i]+c;
  8.         c = A[i]/BIG_CAP;
  9.         A[i]%=BIG_CAP;
  10.     }
  11. }
  12. //
  13. // 大數除大數
  14. // BigA : 被除數，執行完後變餘數
  15. // BigB : 除數
  16. // BigC : 餘數
  17. // return : 成功(1), 失敗(0)
  18. int big_div2(Big * BigA , Big * BigB , Big * BigC)
  19. {
  20.     // 前面全略...
  21.      
  22.     NB = (int*)malloc(sizeof(int) * r_len);
  23.     BigC->dig = c_len, A[a_len++]=0; // A 前面補0
  24.  
  25.     /* 尋找 C[c_len:0] 之猜商遊戲 */
  26.     for( ; c_len >=0 ; --c_len){
  27.         memset( (void*)NB, 0, sizeof(int)*r_len);
  28.         for(N=1; N<BIG_CAP ; ++N) {
  29.             big_adds(NB, B, r_len); // NB+=B, 取代乘法
  30.             cmp = big_compare_part(A+c_len, NB, r_len);
  31.             if(cmp<=0) break; // A >= NB 時結束
  32.         }
  33.         /* 依比較結果 cmp 調整 */
  34.  
  35.         if(cmp==0){ // 相等, 將商填入, A[...] 填入
  36.             memset( (void*)(A+c_len), 0, BSIZE(r_len));
  37.             C[c_len] = N;
  38.         } else if(cmp < 0) { // 猜太大, 變小, 回覆
  39.             C[c_len ] = --N;
  40.             big_subs(NB, B, r_len); // NB-=B, 取代乘法
  41.             big_subs(A + c_len, NB, r_len); // A -=NB
  42.         } else { // 唯一可能剛好情況, N=BIG_CAP, 不需再做回覆
  43.             big_subs(A+c_len, NB, r_len); // A-= NB
  44.             C[c_len] = BIG_CAP2;
  45.         }
  46.     }
  47.     free( (void*)NB);
  48.     // 後面一堆也略 ... 
  49.     return 1;
  50. }

 

[4] 二分搜尋法

 

 

另一種猜商策略是用二分搜尋法，若是 65536 進位時，猜一個商最多需要 16 次。為留參，這份  code 有些長，盡可能縮。只有部份變動，有興趣從 for 開始看便可。

 

Code Snippet
  1. // 大數除大數
  2. // BigA : 被除數，執行完後變餘數
  3. // BigB : 除數
  4. // BigC : 餘數
  5. // return : 成功(1), 失敗(0)
  6. int big_div3(Big * BigA , Big * BigB , Big * BigC)
  7. {
  8.     int i, N, cmp, low, up;
  9.     int * NB; // 存放 N * B
  10.     int * A = BigA->arr, * B = BigB->arr, * C = BigC->arr;
  11.     
  12.     // 算長度
  13.     int a_len = BigA->dig, b_len=BigB->dig;
  14.     int c_len = a_len - b_len + 1 , r_len = b_len+1;
  15.  
  16.     if(!b_len) {    BigC->dig=0;return 0;}
  17.     else if(b_len > a_len) { BigC->dig = 0;return 1;}
  18.  
  19.     if(! (NB = (int*)malloc(BSIZE(r_len)))) return 0; // fail
  20.     BigC->dig = c_len, A[a_len++] = 0;
  21.  
  22.     /* 尋找 C[c_len:0] 之猜商遊戲 */
  23.     for( --c_len; c_len >=0 ; --c_len){
  24.         low = 0 , up = BIG_CAP-1;
  25.         while(low<=up){
  26.             big_mul_int_part(B , N = (low + up) >>1 , NB, r_len);
  27.             if(!(cmp = big_compare_part(A+c_len, NB, r_len))) break;
  28.             else if(cmp>0)  low = N+1; // 猜太小
  29.             else up =  N-1;
  30.         }
  31.         /* 依比較結果 cmp 調整 */
  32.         if(cmp==0){ // 相等, 將商填入, A[...] 填入
  33.             memset( (void*)(A+c_len), 0, BSIZE(r_len));
  34.             C[c_len] = N;
  35.         } else if(cmp < 0) { // 猜太大, 變小, 回覆
  36.             C[c_len ] = --N;
  37.             big_subs(NB, B, r_len); // NB-=B
  38.             big_subs(A + c_len, NB, r_len); // A -=NB
  39.         } else {
  40.             big_subs(A+c_len, NB, r_len); // A-= NB
  41.             C[c_len] = N;
  42.         }
  43.     }
  44.     free( (void*)NB);
  45.     // 更新 c_len
  46.     for(i=BigC->dig-1; i>0 && C[i]==0; --i);
  47.     if(i==0) BigC->dig = ( (C[i]==0) ? 0 : 1 );
  48.     else BigC->dig = i+1;
  49.  
  50.     // 更新 a_len, a_len < b_len
  51.     for(i=b_len; i>0 && A[i]==0 ; --i);
  52.     if(i==0) BigA->dig = ( (A[0]==0) ? 0 : 1);
  53.     else BigA->dig = i+1;
  54.  
  55.     return 1;
  56. }

 

 

 

[5] 有更好的猜商策略嗎？

 

 

上面的二分法和比較，其實寫在一起，效率會較高，乘的時候先乘高位數部份，可以同時做比較，細節留下不贅述！

 

另一種猜商方式，假設是 56789012345 除以  123456，capacity 為 10 進位，商可以先「大致猜一下」，示意如下。

 

　　　　　　　　　　　　　　？？？？？？
　　　　　　　┌────────────────
　１２３４５６│０５６７８９０１２３４５
　　　　　　　│
 

直接先把被除數的前兩位合併成一個數字，變　０５，

除數取出最高位，變１，而５／１＝５，發現用 5 太大，但往下降 1 變 4 後就剛好。

 

　　　　　　　　　　　　　　４？？？？？
　　　　　　　┌────────────────
　１２３４５６│０５６７８９０１２３４５
　　　　　　　│０４９３８２４
　　　　　　　－－－－－－－－－－－－－
　　　　　　　　００７４０６６

 

　　　　　　　　　　　　　　４？？？？？
　　　　　　　┌────────────────
　１２３４５６│００７４０６６１２３４５
　　　　　　　│

 

一樣繼續拿　０７　除以　１，得到　７，依此類推，但這次實際拿到的商是 5 了。不過記得的是，最大只能猜到 BIG_CAP-1　便是。拿到初始的商後，接下來的策略變化就個人發輝吧，像是把它拿來當作是二分法的 upper 初值，或是直接拿下來做遞減等等 ( 猜 7 太大 -> 變 6 -> 變 5 )，這些都可拿來當作策略。

 

另一種猜商策略，忘了在哪看到這手法，出處不詳。用到的是牛頓法，求 A / B 時，相當於使用 A * (1/B)， 1/B 怎來的？

 

由於牛頓法一開始只是需要初值，所以即使用 double 勉強接進來也行，算完 double 轉一次大數，稱 x，再來就是迭代公式了

xn+1 = (2 - xn * B) * xn。

 

舉個例子， 求 7123456(A) / 123(B)，先算 1/123

 

(1) 1 / 123，抓個大概，算 x1 = 0.01 就好。

(2)  x2 = (2-x1*B) * x1 = (2 - 0.01 * 123) * 0.01 = 0.0077

(3)  x3 = (2-x2*B) * x2 = (2-0.0077 * 123) * 0.0077 = 8.10733 e-3

(4)  x4 = 8.130017633 e-3

(5) x5 = 8.1300813e-3

 

由於被除數是 7123456 (A , 7位數)，到這裡誤差已小於 1e-7，就可以停下來了，最後再把這兩個答案相乘得 Y，Y 取整數部份，就是商 C；餘數 R = A - B * C。

 

有幾個 note 可以再想想。

 

(1) 如果乘法沒有降到 nlogn ，沒意義。

說到底除法效能還是相依於乘法(吧？)，依上述方式，最耗時的應是在最後 兩個答案相乘(A, 1/B) 的乘法部份，因前面一開始取的精度可以不用那麼高 (偷懶)。

(2) 萬一，B 本身就超過 DBL_MAX 的時候，該怎辦？？

oh 這問題好解決，如果是 100 進位的話，假設 B = 71 23 45 60 ，第一個非零是在 B[4-1]，先用 1/71 = 0.01428，再放到小數後第 4 個百位就行，也就是 0.00 00 00 01 43 當初值。然後小數點的移位的話.. 嗯，看來還有多一個變數做紀錄。

(3) 萬一，實作的是 65536 進位的話，那不是等於要先將大數轉成 10 進位，求完之後再轉成 65536 進位？速度整個慢了？

這問題推起來卡卡，就沒真的實作不敢給結論了。初步想法是，其實 65536 進位也可以用 shift 方式完成，用上述的原理也沒必要再做轉換；再糟的話，就 KK IEEE754吧 Orz。

 

 

 

[6] 效率 Note

 

這篇寫很多了，下次有空再把知道的一些技巧做 note 。目前就我所知的算法效率

 

(1) 大數 + - 大數 : O(n) , 這怎麼做都是 O(n)

(2) 大數 * 大數 : O( nlogn) , FFT

(3) 大數 / % 大數 : O(n logn logn) , 牛頓法 + FFT 乘法

(4) 字串轉 16 / 65536 進位大數 : O(nlogn) , 用連除法效率是 O(n^2)

(5) 大數階乘 : O(???)

 

大數階乘其實筆者也不確定是 O(??) ，比較快的一種做法是先對各數做質因式分解，再做相乘。

 

怪不得不少論壇都說這是一個「永遠的挑戰」，就算 C 寫好，C++ 考慮 copy constructor 所帶來效能問題應也不會太好搞，除了 open soruce ( ex, gmp) 外，目前似乎還沒人在講這部份的東西。

arrow
arrow
    全站熱搜
    創作者介紹
    創作者 edisonx 的頭像
    edisonx

    Edison.X. Blog

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

    E-mail轉寄