摘要:本系列文章將深入探討大數階乘之計算的各種演算法和實現,從最簡單的大家最容易想到的演算法,到使用各種大數乘法的版本,包括硬乘法、分治法、快速數論變換(FNT)和快速傅立葉變換(FFT)的版本,甚至還有使用組合語言寫成的迷你版,使用斯特林公式的極速版。這些系列文章不僅分析和講述演算法思想,而且給出一些版本的部分或全部代碼,並探討代碼的最佳化方法。除了重點討論自己寫的程式外,該系列文章還會對網上的一些計算大數階乘的程式和論文以及一些學術刊物上的論文。 大數階乘的計算是一個有趣的話題,從中學生到大學教授,許多人都投入到這個問題的探索和研究之中,並發表了他們自己的研究成果。如果你用階乘作關鍵字在google上搜尋,會找到許多此類文章,另外,如果你使用google學術搜尋,也能找到一些計算大數階乘的學術論文。但這些文章和論文的深度有限,並沒有給出一個高速的演算法和程式。
我和許多對大數階乘感興趣的人一樣,很早就開始編製大數階乘的程式。從2000年開始寫第一個大數階乘程式算起,到現在大約己有6-7年的時光,期間我寫了多個版本的階乘計算機,在階乘計算機的演算法探討和程式的編寫和最佳化上,我花費了很大的時間和精力,品嘗了這一過程中的種種甘苦,我曾因一個新演算法的實現而帶來速度的提升而興奮,也曾因費了九牛二虎之力但速度反而不及舊的版本而懊惱,也有過因解一個bug而通宵敖夜的情形。我寫的大數階乘的一些程式碼片段散見於互連網絡,而演算法和構想則常常縈繞在我的腦海。自以為,我對大數階乘計算機的演算法探索在深度和廣度上均居於先進水平。我常常想,應該寫一個關於大數階乘計算的系列文章,一為整理自己的勞動成果,更重要的是可以讓同行分享我的知識和經驗,也算為IT界做一點兒貢獻吧。
我的第一個大數階乘計算機始於2000年,那年夏天,我買了一台電腦,開始在家專心學習VC,同時寫了我的第一個VC程式,一個仿製windows介面的計算機。該計算機的特點是高精度和高速度,它可以將四則運算的結果精確到6萬位以內,將三角、對數和指數函數的結果精確到300位以內,也可以計算開方和階乘等。當時,我碰巧看到一個叫做實用計器的軟體。值得稱頌的是,該計算機的作者薑邊竟是一個高中生,他的這個計算機功能強大,獲得了各方很高的評價。該計算的功能之一是可以計算9000以內階乘的精確值,且速度很快。在佩服之餘,也激起我寫一個更好更快的程式的決心,經過幾次改進,終於使我的計算機在做階乘的精確計算時 (以9000!為例),可比實用計算機快10倍。而當精確到30多位有效數字時,可比windows內建的計算機快上7500倍。其後的2001年1月,我在csdn上看到一個貼子,題目是“有誰可以用四行代碼求出1000000的階乘”,我回複了這個貼子,給出一個相對簡潔的代碼,這是我在網上公布的第一個大數階乘的程式。這一時期,可以看作是我寫階乘計算機的第一個時期。
我寫階乘計算機的第二個時期始於2003年9月,在那時我寫了一組專門計算階乘的程式,按運行速度來分,分為三個層級的版本,初級版、中級版和進階版。初級版本的演算法許多人都能想到,中級版則採用大數乘以大數的硬乘法,進階版本在計算大數乘法時引入分治法。期間在csdn社區發了兩個貼子,“擂台賽:計算n!(階乘)的精確值,速度最快者2000分送上”和“擂台賽:計算n!(階乘)的精確值,速度最快者2000分送上(續)”。其進階算的版本完成於2003年11月。此後,郭先強於2004年5月10日也發表了系列貼子,“擂台:超大整數高精度快速演算法”、“擂台:超大整數高精度快速演算法(續)”和“擂台:超大整數高精度快速演算法(續2)”, 該貼重點展示了大數階乘計算機的速度。這個貼子一經發表即引起了熱列的討論,除了我和郭先強先生外,郭雄輝也寫了同樣功能的程式,那時,大家都在持續改進自己的程式,看誰的程式更快。初期,郭先強的稍稍領先,中途郭子將apfloat的代碼應用到階乘計算機中,使得他的程式勝出,後期(2004年8月後)在我將程式作了進一步的改進後,其速度又稍勝於他們。在這個貼子中,arya提到一個開放源碼的程式,它的大數乘法採用FNT+CRT(快速數論變換+中國剩餘定理)。郭雄輝率先採用apflot來計算大數階乘,後來郭先強和我也參於到apfloat的學習和改進過程中。在這點上,郭先強做得非常好,他在apfloat的基礎上,成功地最佳化和改時演算法,並應用到大數階乘計算機上,同時他也將FNT演算法應用到他的<超大整數高精度快速演算法庫>中,並在2004.10.18正式推出V3.0.2.1版。此後,我在2004年9月9日也完成一個採用FNT演算法的版本,但卻不及郭先強的階乘計算機快。那時,郭先強的程式是我們所知的運算速度最快的階乘計算機,其速度超過久負盛名的數學軟體Mathematica和Maple,以及通用高精度演算法庫GMP。
我寫階乘計算機的第三個時間約開始於2006年,在2005年8月收到北大劉楚雄老師的一封e-mail,他提到了他寫的一個程式在計算階乘時比我們的更快。這使我非常吃驚,在詢問後得知,其核心部分使用的是ooura寫的FFT函數。ooura的FFT代碼完全公開,是世界上啟動並執行最快的FFT程式之一,從這點上,再次看到了我們和世界先進水平的差距。佩服之餘,我決定深入學習FFT演算法,看看能否寫出和ooura速度相當或者更快的程式,同時一個更大計劃開始形成,即寫一組包括更多演算法的階乘計算機,包括使用FFT演算法的終極版和使用無窮級數的stirling公式來計算部分精度的極速版,除此之外,我將重寫和最佳化以前的版本,力爭使速度更快,代碼更優。這一計劃的進展並不快,曾一度停止。
目前,csdn上blog數量正在迅速地增加,我也萌生了寫blog的計劃,藉此機會,對大數階乘之計算作一個整理,用文字和代碼詳述我的各個版本的演算法和實現,同時也可能分析一些我在網上看到的別人寫的程式,當然在這一過程中,我會繼續編寫未完成的版本或改寫以前己經實現的版本,爭取使我公開的第一份代碼都是精品,這一過程可能是漫長的,但是我會儘力做下去。