艱苦環境裡的大整數精確計算

曾慶潭2011-07-16於紐西蘭
Ching-Tang Tseng
ilikeforth@gmail.com
http://forthfortnight.blogspot.com

相當的年紀才能有相當豐富的回憶，過往的時日還不能醉生夢死，回憶才有意義。這幾天，這裡天氣變得很冷了，冷到坐下來就能大量的回想到過去的事情，因為，不少老友都找上門來了，自然就勾起了我的許多回憶。
先是彰化的林伯龍醫師，透過此處的網文，捎來了訊息，彼此都感到非常意外。我也立刻回信，告訴他，不會忘記將近30年前，中華民國FORTH語言協會草創初期，他曾經給過我們的金錢資助，並請他問候也住彰化的陳信良先生。當年，他們都捐助過協會不少錢，協會卻是一個地下組織，我拿了錢就能辦事，搞的也不是什麼偉大革命，只是推廣FORTH。
從小，我也經歷過許多的欺壓與霸凌，但對方後來的下場都不太好。打過我的人，有的還都死掉了。我也常聽說，欺負過我的人，後來的生活都很悽涼、潦倒。可是，曾經幫助過我的人，再度與我聯絡上時，看來都有了滿意的人生，林醫師就是其中之一。我不是鬼神，也不會詛咒別人，蒙受任何打擊時，內心只會要求自己更堅強的站起來，更要充實自己，也有著想當世界天皇的偉大志向，其他甚麼天皇我都不會放在眼裡。
30幾年前，職業軍人敢在社會上搞地下組織，是屬於違反軍紀的重大犯罪行為，嚴重的，可以撤職查辦，我卻胆大妄為到在社會上私募款項。曾經被也是同行的退伍軍人，一狀告回當時的服務單位，他們把持著國營公司主管的權位，不捐錢也就罷了，還陷害忠良，差一點就令我落得與郭冠英同樣的下場。幸好，承辦告案的魏元勳副所長深明大義，放我一馬，對我算是有著不殺之恩。後來，他出題考我，能不能用FORTH把圓周率π算到小數點後面1000位數？我算不出來，於是FORTH被判定成沒有用。我終身記得此事，便立定志向要自己設法解決，雖然多年以前，我就能解這個問題，但我要等到電腦軟硬體的發展都足夠成熟之後，才願意實際使用FORTH解決問題。能夠解決，才值得在網文中公佈，算是了了一股心願，從此之後，也才可以不必再談π的事情。
我有終身反省的習慣，2009-5-15刊出的：『FORTH解大數問題』談到π，那時，我用一般FORTH程式解問題，可以獲得正確答案，但是沒有特別意義。前一篇網文公佈了我想將π算到幾位數都可以了，就是為了回報魏副所長的不殺之恩，用的不僅只是實際程式，還是我親自設計出來的程式語言系統。以前，我若堅持已見，想做這種發展，會被視同為不務正業，結果也可能會是撤職查辦。我的本業被上級視同為掃放射性的工作，而且全單位都沒人愛做，我就只好忍氣吞聲的一直幹著吃放射性的工作，難怪七早八早就被國軍淘汰了，也不夠資格領終身俸，但沒有怨言，現在只能搞FORTH。

我還有另外兩個必須反省之後，透過網文向讀者交待的問題，但不必計較會有多少讀者來點閱。2010-01-16刊出的網文：『高級數學函數』，文內提到，當年一月八日出刊的『科學美國人』(Scientific American)裡面的一篇文章，論及一個十進制為232位數的密碼性數字，可以因式分解成兩個均為116位數的質數。該文中，我寫下誓言，『它日我必會用到這份資料驗證程式』，現在，實現這個諾言，程式要用大整數除法算對結果才有價值，因為，大數的加、減、乘法都很容易設計，除法不容易。
列示這個程式，代表我對所有研究問題認真的態度，不能設計出算對結果的系統前，我長期被大數無法正確輸出與輸入系統的問題所困擾，請教國際網站的負責人，也得不到答覆，他們反而要我幫忙做另外一件事情。最後，只好自行解決問題，完成這件工作後，令我憶起，從小都是這樣子自立自強長大的，我經常鼓勵中國同胞努力站起來，我自己更應該努力站起來。

0 BIGVARIABLE A1 1000 ALLOT
0 BIGVARIABLE B1 1000 ALLOT
0 BIGVARIABLE C1 1000 ALLOT

: MAIN BASIC
10 LET B{ A1 =
S" 1230186684,5301177551,3049495838,4962720772,8535695953" FIRST-SECTION
S" 3479219732,2452151726,4005072636,5751874520,2199786469" NEXT-SECTION
S" 3899564749,4277406384,5925192557,3263034537,3154826850" NEXT-SECTION
S" 7917026122,1429134616,7042921431,1602221240,4792747377" NEXT-SECTION
S" 9408066535,1419597459,8569021434,13"
NEXT-SECTION
S>BIG }B
20 LET B{ B1 =
S" 3347807169,8956898786,0441698482,1269081770,4794983713" FIRST-SECTION
S" 7685689124,3138898288,3793878002,2876147116,5253174308" NEXT-SECTION
S" 7737814467,999489"
NEXT-SECTION
S>BIG }B
30 LET B{ C1 = A1 / B1 }B
40 PRINT " A1 = : "
50 RUN A1 BIG.
60 PRINT " B1 = : "
70 RUN B1 BIG.
80 PRINT " C1=A1/B1= : "
90 RUN C1 BIG.
100 END
;

執行結果如下：

MAIN
A1 = :
232 digits
12301866845301177551304949583849627207728535695953
34792197322452151726400507263657518745202199786469
38995647494277406384592519255732630345373154826850
79170261221429134616704292143116022212404792747377
94080665351419597459856902143413
B1 = :
116 digits
33478071698956898786044169848212690817704794983713
76856891243138898288379387800228761471165253174308
7737814467999489
C1=A1/B1= :
116 digits
36746043666799590428244633799627952632279158164343
08764267603228381573966651127923337341714339681027
0092798736308917 ok

驗證無誤。

另一件必須向讀者交待的事情，是有關於算出超大質數(Prime number)之事，我在2009-02-10『質數』網文中，提供過算出質數之一般FORTH程式，它並不是什麼偉大的創作，只適合用來當作小型範例程式，當時(2008-09-18)梅森質數大搜索組織(GIMPS, Great Internet Mersenne Prime Search)活動的結果，已經算出了2的43112609次方減1為新發現的大質數，此數大於1200萬位數(10^12978188)，現在該組織可能已經發現更大的數字了。
我重提此事的目的，是想告訴大家，尤其是全體中國人，找梅森質數的數學計算方法非常簡單，但一般程式語言系統不易設計出那樣的程式，必須是擁有無窮位數大整數計算功能的程式語言系統，才能用來輕易的設計出程式。我不想批評別種程式語言在這方面的能力，只問自己設計出來的ABC FORTH數學計算系統，又能做到什麼程度？現在則是展示本領的時候了。
因為我已設計出了輕巧的工具，可以實現這套尋找大質數方法的程式，並快速執行後，得到測試結果。也因此能通順、詳細的解釋這套技術。我們先行詳細解釋一些相關術語。
任何數只要可以表示為M(P) = 2 ^ P – 1，就被稱為梅森數(Mersenne number，梅森是一位法國數學家，1588~1648)，如果P不是質數，則M(P)也不會是質數，若M(P)為質數，則P必定為質數，但是，P為質數，M(P)卻不一定為質數。
梅森當初提出上述研究時，宣稱：當P = 2,3,5,7,13,17,19,31,67,127,257等11個質數時，M(P)都是質數，而除了這11個質數以外，比257小的其他質數P，計算出來的M(P)都不是質數。
可是，後人發現M(67)不是質數，M(257)也不是質數，但M(61)、M(89)、M(107)反而都是質數。
當P為較大的質數時，能夠計算出來的M(P)就會非常的大，因此，就必須經過更大數字的大量計算測試，才能確定該M(P)是否為質數？早期沒有電腦可用的時代，這種工作，實在難以實現。後來，法國數學家盧卡斯(Edouard Lucas)，於1856年與1878年，分別提出過發現了用來測試M(P)是否為質數的方法，1930年代，萊默(D. N. Lehmer)又研究出改良的測試方法，這種測試方法，後來就被稱為盧卡斯萊默測試法(Lucas Lehmer test, LLT)，可以高效率測試出梅森數是否為一質數？目前，此法是被普遍用來發現梅森數中更大質數的最佳方法。國際網路梅森質數大搜索組織，就是專門從事於這種活動的一個組織。
盧卡斯萊默測試法，可以用下列程式語言式的敘述表達清楚，意義能夠直讀，我就不用詳細解釋。

\ Determine if M ( P ) = 2 ^ P – 1 is a prime?
\ Lucas-Lehmer test
Lucas-Lehmer( P )
var S = 4
var M = 2 ^ P – 1
repeat P – 2 times:
S = ( ( S * S ) – 2 ) MOD M
if S = 0 return PRIME else return COMPOSITE

但是，讀者請注意，敘述中的S及M都必須能夠被處理成無限位數的大整數，上列程式才有效，P則暫時還可以被處理成21億以下的一般單整數。不是所有現成的程式語言都能辦到此事，能辦到的大概也不會是免費的軟體。我改寫出一個精簡交談式的ABC FORTH程式，簡單到就如下列所示，重點顯現在可以快速測試出能夠告訴大家特性的實際數字結果。

INTEGER I
INTEGER P

0 BIGVARIABLE S 10008 ALLOT
0 BIGVARIABLE M 10008 ALLOT

0 BIGVARIABLE BIG0
1 BIGVARIABLE BIG1
2 BIGVARIABLE BIG2

: LUCAS BASIC
10 PRINT " Input P for (2^P-1), but P must be > 2 !!! "
20 INPUTI P
30 RUN CR CR
40 LET B{ S = I>BIG ( 4 ) }B
50 LET B{ M = BIG2 ^ P - BIG1 }B
70 FOR I = 1 TO P - 2
80 LET B{ S = ( ( S * S ) - BIG2 ) MOD M }B
90 NEXT I
100 IF B{ S = BIG0 }B THEN 130
110 PRINT " For P = " ; P ; " : (2^P-1) is a composite. "
120 GOTO 150
130 PRINT " For P = " ; P ; " : (2^P-1) is a prime. It is: "
140 RUN M BIG.
150 RUN CR CR
160 END
;

程式再經簡單的改寫、調整，就能自動的執行從所有梅森數中找出梅森質數的工作，以下為只計算到2000位數時的印出結果：

For P = 3 : (2^P-1) is a prime. It is:
1 digits
7
For P = 5 : (2^P-1) is a prime. It is:
2 digits
31
For P = 7 : (2^P-1) is a prime. It is:
3 digits
127
For P = 13 : (2^P-1) is a prime. It is:
4 digits
8191
For P = 17 : (2^P-1) is a prime. It is:
6 digits
131071
For P = 19 : (2^P-1) is a prime. It is:
6 digits
524287
For P = 31 : (2^P-1) is a prime. It is:
10 digits
2147483647
For P = 61 : (2^P-1) is a prime. It is:
19 digits
2305843009213693951
For P = 89 : (2^P-1) is a prime. It is:
27 digits
618970019642690137449562111
For P = 107 : (2^P-1) is a prime. It is:
33 digits
162259276829213363391578010288127
For P = 127 : (2^P-1) is a prime. It is:
39 digits
170141183460469231731687303715884105727
For P = 521 : (2^P-1) is a prime. It is:
157 digits
68647976601306097149819007990813932172694353001433
05409394463459185543183397656052122559640661454554
97729631139148085803712198799971664381257402829111
5057151
For P = 607 : (2^P-1) is a prime. It is:
183 digits
53113799281676709868958820655246862732959311772703
19231994441382004035598608522427391625022652292856
68889329486246501015346579337652707239409519978766
587351943831270835393219031728127
For P = 1279 : (2^P-1) is a prime. It is:
386 digits
10407932194664399081925240327364085538615262247266
70480531911235040360805967336029801223944173232418
48424216139542810077913835662483234649081399066056
77320762924129509389220345773183349661583550472959
42054768981121169367714754847886696250138443826029
17323488853111608285384165850282556046662248318909
18801847068222203140521026698435488732958028878050
869736186900714720710555703168729087

我以不建立質數表的方式，讓上列程式就直接逐次加1的往上計算到2000為止，所以可以確定其間沒有其它的質數是梅森數了。我寫出程式的目的，不在與前人比高下，只想用程式證明，除了上列數字外，確實沒有其他的梅森質數了。我沒有用我的電腦繼續算下去，只算出了刊出這篇文章所需要的數字就夠了，讀者大約很難在其他的場合，見到一個名符其實的386位數的真正質數了。我將系統暫時設定成可以得到一個萬位數以下的質數，因此，就目前全世界已經公開的47個梅森質數而言，我的系統不改設定，可以算得第32個以下的梅森質數。修改設定，可以提高記錄，但我只有艱苦環境裡的電腦，不必與有錢人爭長短。我的電腦在這個程式執行期間，CPU的冷卻風扇不斷的起動運轉，這是紐西蘭冬日裡，環境溫度接近0度C時的罕見現象，可見它們跑得有多累！
這些數據都是有用的大數據，用來測試系統時可以派得上用場。我從科學美國人期刊上取得116位數的質數，就能檢驗出系統的輸出輸入，會受到微軟作業系統內不良設計的嚴重影響。大家都花了那麼多錢向微軟購買軟體，受到這種氣，卻不能告它來索賠，只能透過這種技術網文，將問題披露出來。我也不是光會在社會上取得資源而不回饋的人，上述資料，就回饋給全世界，比科學美國人期刊能夠提供的數字大三倍位數的質數，更大位數的資料，我也能夠產生。
我在2009-02-10：『質數』網文中，詳細介紹的程式及其執行結果，對找出超大質數的研究也有幫助。如果將那些質數建成檔案，再利用本文內的程式，兩者配合起來執行，就能更快速的找出更大的梅森質數，不必以每次加1的方式，逐次的往上計算。兩者配合起來的找法，可以輕易的執行到2的二十一億次方的計算，可以超過第47個梅森質數的限量，但我對安排硬體來執行可能可以打破記錄的事情，沒有興趣。

談完了技術，可以話點家常。

前兩天，一對老朋友夫婦來訪，順便找我協修一些修不好的電器，肉眼能夠瞧出來的問題，大概都難不倒我，所以東西很快就都修好了。
第一個問題電器，是被拆到螺絲頭完全毀損、變型而無法受力之深槽內的螺絲，老友確定已經拆不下來了，所以帶來找我協修。我找遍了自己花了不少錢買來的所有工具，確實也沒有適當的螺絲起子，可以派得上用場，因為那根螺絲是一個故意將螺絲頭設計成三線凹槽式的螺絲，目的就在防止一般使用者，任意自行拆開電器來維修，發生觸電意外。
後來，我突然想到，一般人動手工作的方式與我不同，我因經常自己動手做事，知道一些工作技巧，也許那根螺絲換我直接拆，就能拆得下來，因此，就取用了一根1字型的長起子直接拆卸，果然不錯，深槽內螺絲頭全毀成圓洞的螺絲，我直接拆出來了，電路斷掉了的電器，當然也就修好了。
第二個問題電器，是一個新型獨立的液晶顯示幕，不能正常工作，系統告知傳輸線有一條是斷訊狀態，老友不敢擅拆機背板來查線，所以拿來找我協修。我只看了看整組機子的情況，就告訴他，這種東西敢出廠銷售，多股訊號傳輸線發生機內斷線的可能機率相當低，反倒是訊號傳輸接頭又接了轉換接頭，出問題的可能性較高，不如先拆驗接頭。直接拆開一看，又是果然不錯，D型多孔接頭的一個凸針，根本沒插入對應插孔，全被壓歪了。以小起子併同尖嘴鉗，小心扳回銅針後，看來還能使用，修好了這個大概並不便宜我都還沒買過的液晶顯示幕。
做這些事情，會令我想起自己曾經有過的艱苦人生。
我拆螺絲時，知道扳轉起子的第一時間，必須壓得夠緊，並只施以小角度的衝擊式扭力，感到螺絲確實鬆開時，必須立刻解除壓緊之力，改成單純的只施以扭轉之力，特別是金屬螺絲安裝於塑膠機殼的情況尤然，否則螺絲要不是拆壞螺頭都還拆不下來，就是發生滑牙現象，從此之後都難再拆裝了。英文有句俗語，能幫助經常搞不清楚扭轉方向的人，記住『旋緊』或『扭鬆』的方向，他們說：『Right ~y tight ~y, Left ~y loose』，~y要拖長了發[i]的音，是『右旋變緊，左旋變鬆』的意思。
我的表達能力很好，這一段話是標準起子使用技巧的經典說明，如果您能讀懂這一段話，終身牢記，身體力行，保證您受用不盡。至於話語的意義，您應自己體會，我大學學的是機械工程，是一名具有工程師執照的淘汰老兵，當然能夠詳細解釋出道理，但我會告訴想學的人，就像學FORTH程式語言一樣，您只宜拿著起子直接動手照著做，多做熟練之後，手腦才能配合得起來，光知道道理沒有用。大家同時拿著起子拆螺絲，同做這麼簡單的一件工作，為什麼會拆出不一樣的結果？因為我出身卑微，有過艱苦的人生。
從小到井邊打水、挑水回家、生火煮飯，其實都還稱不上艱苦，但現代人都不必做了。
我沒上學之前，就陪著家人到海邊防風林撿柴火，當個小樵童，經常要躲避森林巡守員的為難，這也不算什麼，但聽到穿林的松風，環顧四下無人的荒野森林時，心中卻油然生起強烈的意志，立志將來要讓家人不必這樣過日子。
我年紀雖小，找柴迅速，每次都是我陪著媽媽最先回到產業道路旁約好的集合地點，我們要先藏好媽媽的一擔柴，與我能背的一竹簍子柴，等了一陣子還不見鄰居來會合，我媽媽就會要我去路的盡頭把風，她才高喊：
『花~ ~ 草嬸媽哦 ~ ~，好轉去了 ~ ~ 』(要用台語喊，是要花草嬸回去的意思)
這樣的吆喝，一輩子在我腦海裡迴盪。我永遠記得，來到這個世界，媽媽把我養育長大，我負有重任，要改善我們的環境，要改善我們的生活，一生遭受任何打擊，都不得低頭，必須昂首闊步，堅強前行。
從小生活環境不好，經常要倒馬桶，有夠粗鄙吧？那又如何？夜晚穿過黑漆漆的暗巷，端著馬桶，前往公廁後面的糞便池，您以為我心裡能想甚麼？我只想到將來絕不讓我的家人，還要這樣過日子。
我14歲時就上過大夜班，在木器工廠當童工，賺錢貼補家用，清晨四點鐘，下班騎自行車回家，一路上只有我孤獨的一個人，也是黑漆漆的，只有星光照著路面，路過無人看守剛撞死人的鐵路平交道，那時，屍體都只能放在現場好幾天，只蓋上草蓆，您以為我這個童工這麼勇敢，什麼都不怕？但是我從不跟家人提起這些事情，說出來就表示自己將來一定沒出息。
大學畢業時，母親過世了，家人生活還是很苦，我沒能盡任何孝道，我媽媽也從來不曾要求我要回報她。因此，改善環境、改善生活的志向，轉而奉獻社會，這絕對就是我爸爸、媽媽的意思。
有過這些卑下的人生，才會顯示出我拆螺絲與眾不同之處，通常光憑肉眼看不出來。就算成年踏入社會了，我也並非就能過著好日子，一切仍然都要靠著自己，不努力就不可能有好的結果。前面說過，吃公家飯時，我盡是幹著別人不愛幹的掃放射性工作，那是打從進入公家單位開始，就長期被人欺負的結果，後來毅然決然的離職了，還脫離不了不良的影響。
我估計台灣當今也沒幾個能把放射性測量工作做得比我確實的人，十五年前，移民之前，我自購測量放射性用的蓋革─米勒計數器(Geiger-Muller counter)，自行測量家中器具有無輻射？以便決定是否搬到紐西蘭繼續使用。我熟知儀器的性能，就算沒有校準儀器的機會，我也還能量得出有意義的數據。只要儀器能夠反應，您將測量時間加長，精確程度就可增加，這是不變的真理，讀數的意義就靠比較，也能凸顯出來。最後，量出了床上長期使用的兩個枕頭，放射性計量讀數，是環境平均背景讀數的三倍，我知道輻射不算太強，等同於大家每次搭飛機，飛上高空後所必須接受的輻射強度，那時，搭飛機時，沒有不得使用任何儀器的管制，我用同一儀器在搭飛機時測量了幾小時，確實獲得這樣的數據，所以也沒必要將枕頭當輻射物質處理，我太太氣得當場就將兩個枕頭當垃圾丟了。長期在原子爐內工作，不可能每天下班前洗完澡才回家，頭髮就像一個標準的空氣過濾器，什麼東西都能吸附，長期工作下來，能不帶點放射性物質回家？
同樣的，長期使用電腦執行數學計算程式，在現今全世界的數學計算軟體環境中，我怎麼不知道誰會欺負誰？我在學習原子爐爐心設計工程時，是拉計算尺過日子的，那時，買一台進口的小型計算器(Calculator)要幾萬台幣，教我這門課程的鄭季康博士，有一台公家買的裝備，我每星期寫好50張單光紙，留下空格，標示程序，就等著鄭老師來教課時借用他的計算器，利用下課時間，一口氣填完表格，晚上再拉計算尺驗證。這是大四下學期的事情，四下的成績不列入畢業後分發工作的參考，全班同學都不想學太囉嗦的東西了，只剩我一個人努力完成原子爐爐心的工程設計，只有我交了報告，裡面全是大量的數學計算。
今天，中國人的環境，可以說電腦已經普及了，可是，有誰還會在乎與電腦相關的數學計算？大家用的是花大錢買來的鴉片式軟體，吃了就上瘾了，誰還想要自行設計？當老師的，選擇軟體之依據只有一個，看看是否有利可圖？管它什麼公益軟體？我聽到台灣開FORTH月會時，講員吐出來的心聲，學生拿到文憑進入社會工作後，沒人想用老師強教、強用的鴉片式軟體，要不然就盜用。還想正式的用，就得花大錢，而且愈來愈不知道裡面在算些什麼？大家還是希望能夠見到公益性軟體，可是，去那裡找？請問？我怎麼不知道誰在欺負誰！
最近，寫完書後，開始貼一些比較特殊的數學計算網文，有幾天，見到我的網文點閱記錄顯示，德國人的點閱數，比台灣的國人同胞還要高，大家認為我的心裡應該有什麼感想？德國人全國國民平均素質，優於其它國家，有許多技術冠蓋全世界，他們可能比台灣人更了解，在FORTH這麼艱苦的環境中，一個人單打獨鬥，產生出來的創作，比什麼都可貴，所以願意看台灣人不願意看的東西，是中文也無妨。
不單純只是我用起子與別人用起子的方式有所不同，也不單純只是我用電腦與別人用電腦的方式有所不同，更不單純只是我用FORTH與別人用FORTH的方式有所不同。我過過苦日子，從艱苦環境中長大，意志與人不同，心中有點分寸，我寫的全是通順的中文，要想增進全人類的福祉，也是從自己同胞的生活環境作為開始呀！