0. 緣由

去年某位前輩在網路上發表了一段文章，其中討論了圓周率的計算，當看到
方法3時小弟真是佩服的五體投地，深覺自己才疏學淺:

句末這位前輩留一下句話：「您覺得您要被訓練多久？才會有這種想法啊？！」

當時小弟心理想，哇賽！少說也要一百幾十年吧(話說絕世武功都要練個一百幾十年)，
但是事情真的是這樣嗎？

 

1. 考證

後來，作者看了 Dr. Jack W. Crenshaw 在雜誌上發表的一篇文章 - More Strange Integer，
才知道6千年前埃及人就已經知道類似的作法，假如你跟我當初一樣覺得很神奇的話那表示
你跟小弟一樣該回家多唸書，埃及人的作法是:

哇賽，誤差居然只有 0.000008%，神奇吧！而且不只可以這麼做，連也可以如法炮製:

問題是，這些神奇的技巧是怎麼來的？有沒有法則可尋？

 

2. 實踐

假如你用 355/113 這個法子，前輩還是會笑你分子 355 超過 8 位元能表示的範圍，
113 也不是 2 的冪次，用在 8051 這種 8 位元系統上還是不適合，那 201/64 這個完美
的結果是如何找出來的？

Dr. Crenshaw 給了一個很直接的答案，用暴力法！就以為例，取到小數點後第八位:

我們試圖找出近似 x 藉由比值:

用國小就會的交叉相乘，我們可以得到:

接著我們用 C++ 寫一個很簡單的程式，把每個 q 都試試看(只取 2 的冪次)，並且 p 超過
255 的都去掉不要，並且把結果排序:

#include <vector>
#include <stdio.h>
#include <math.h>
#include <algorithm>

#define LIMIT 256

struct Ratio{
   int p, q;
   double err;
};
   

bool less_err(const Ratio& lhs, const Ratio& rhs)
{
   return lhs.err < rhs.err;
}

typedef std::vector<Ratio> RatioVector;

int main(void)
{
   double x; 
   int q;
   RatioVector rv;
   
   x = 3.14159265f;
   
   for(q = 1; q < LIMIT; q<<=1)
   {
       Ratio r;
       r.q = q;
       r.p = int(x * r.q + 0.5);
       double tmp = (double)r.p/r.q;
       r.err = ::abs(x - tmp) * 100.0f;
       if(r.p < LIMIT)
           rv.push_back(r);
   }
   
   std::sort(rv.begin(), rv.end(), less_err);
   
   for(size_t i = 0; i < rv.size(); ++i)
   {
       printf("p = %d, q = %d, error = %f%%\n", rv[i].p, rv[i].q, rv[i].err); 
   }
   
   return 0;
}

 

最後我們看看輸出結果:

果然最佳結果是 201/64！現在各位應該知道201/64並不是亂試試出來的，而是用一個小學就知道
的公式叫電腦去幫我們算(也許前輩提到小學就學過的演算法就是指這個)。各位也可以拿來試
試看，應該也會得到類似的結果。

 

3. 結語

前輩在文中提到很多台灣工程師不事先好好分析就急著寫程式，除了這個毛病外，其實台灣工程師
也很少會找國外的技術刊物來看，各位知道我提的那篇文章是什麼時候刊出的嗎？1995年4月！也就是
12年前國外的雜誌就有寫了，而且這些出版社都很好心幫你把文章收集成一張 CD 讓你網路訂購，
現在您覺得看到問題就馬上硬幹還是應該先看看前人怎麼做的呢? ...(^_^)