計算圓周率

編輯: 逍遙路 關鍵詞: 高中數(shù)學 來源: 高中學習網(wǎng)


古今中外,許多人致力于圓周率的研究與計算。為了計算出圓周率的越來越好的近似值,一代代的數(shù)學家為這個神秘的數(shù)貢獻了無數(shù)的時間與心血。十九世紀前,圓周率的計算進展相當緩慢,十九世紀后,計算圓周率的世界紀錄頻頻創(chuàng)新。整個十九世紀,可以說是圓周率的手工計算量最大的世紀。進入二十世紀,隨著計算機的發(fā)明,圓周率的計算有了突飛猛進。借助于超級計算機,人們已經(jīng)得到了圓周率的2061億位精度。歷史上最馬拉松式的計算,其一是德國的Ludolph Van Ceulen,他幾乎耗盡了一生的時間,計算到圓的內(nèi)接正262邊形,于1609年得到了圓周率的35位精度值,以至于圓周率在德國被稱為Ludolph數(shù);其二是英國的William Shanks,他耗費了15年的光陰,在1874年算出了圓周率的小數(shù)點后707位?上,后人發(fā)現(xiàn),他從第528位開始就算錯了。把圓周率的數(shù)值算得這么精確,實際意義并不大,F(xiàn)代科技領域使用的圓周率值,有十幾位已經(jīng)足夠了。如果用Ludolph Van Ceulen算出的35位精度的圓周率值,來計算一個能把太陽系包起來的一個圓的周長,誤差還不到質(zhì)子直徑的百萬分之一。以前的人計算圓周率,是要探究圓周率是否循環(huán)小數(shù)。自從1761年Lambert證明了圓周率是無理數(shù),1882年Lindemann證明了圓周率是超越數(shù)后,圓周率的神秘面紗就被揭開了。現(xiàn)在的人計算圓周率, 多數(shù)是為了驗證計算機的計算能力,還有,就是為了興趣。

圓周率的計算方法

古人計算圓周率,一般是用割圓法。即用圓的內(nèi)接或外切正多邊形來逼近圓的周長。Archimedes用正96邊形得到圓周率小數(shù)點后3位的精度;劉徽用正3072邊形得到5位精度;Ludolph Van Ceulen用正262邊形得到了35位精度。這種基于幾何的算法計算量大,速度慢,吃力不討好。隨著數(shù)學的發(fā)展,數(shù)學家們在進行數(shù)學研究時有意無意地發(fā)現(xiàn)了許多計算圓周率的公式。下面挑選一些經(jīng)典的常用公式加以介紹。除了這些經(jīng)典公式外,還有很多其他公式和由這些經(jīng)典公式衍生出來的公式,就不一一列舉了。

1、 Machin公式


這個公式由英國天文學教授John Machin于1706年發(fā)現(xiàn)。他利用這個公式計算到了100位的圓周率。Machin公式每計算一項可以得到1.4位的十進制精度。因為它的計算過程中被乘數(shù)和被除數(shù)都不大于長整數(shù),所以可以很容易地在計算機上編程實現(xiàn)。

Machin.c 源程序

還有很多類似于Machin公式的反正切公式。在所有這些公式中,Machin公式似乎是最快的了。雖然如此,如果要計算更多的位數(shù),比如幾千萬位,Machin公式就力不從心了。下面介紹的算法,在PC機上計算大約一天時間,就可以得到圓周率的過億位的精度。這些算法用程序?qū)崿F(xiàn)起來比較復雜。因為計算過程中涉及兩個大數(shù)的乘除運算,要用FFT(Fast Fourier Transform)算法。FFT可以將兩個大數(shù)的乘除運算時間由O(n2)縮短為O(nlog(n))。

2、 Ramanujan公式


1914年,印度數(shù)學家Srinivasa Ramanujan在他的論文里發(fā)表了一系列共14條圓周率的計算公式,這是其中之一。這個公式每計算一項可以得到8位的十進制精度。1985年Gosper用這個公式計算到了圓周率的17,500,000位。

1989年,David & Gregory Chudnovsky兄弟將Ramanujan公式改良成為:

這個公式被稱為Chudnovsky公式,每計算一項可以得到15位的十進制精度。1994年Chudnovsky兄弟利用這個公式計算到了4,044,000,000位。Chudnovsky公式的另一個更方便于計算機編程的形式是:

3、AGM(Arithmetic-Geometric Mean)算法

Gauss-Legendre公式:

初值:

重復計算:

最后計算:

這個公式每迭代一次將得到雙倍的十進制精度,比如要計算100萬位,迭代20次就夠了。1999年9月Takahashi和Kanada用這個算法計算到了圓周率的206,158,430,000位,創(chuàng)出新的世界紀錄。

4、Borwein四次迭代式:

初值:

重復計算:  

最后計算:

這個公式由Jonathan Borwein和Peter Borwein于1985年發(fā)表,它四次收斂于圓周率。

5、 Bailey-Borwein-Plouffe算法

 

這個公式簡稱BBP公式,由David Bailey, Peter Borwein和Simon Plouffe于1995年共同發(fā)表。它打破了傳統(tǒng)的圓周率的算法,可以計算圓周率的任意第n位,而不用計算前面的n-1位。這為圓周率的分布式計算提供了可行性。1997年,F(xiàn)abrice Bellard找到了一個比BBP快40%的公式:

圓周率的計算歷史

時間

紀錄創(chuàng)造者

小數(shù)點后位數(shù)

前2000

古埃及人

1

前1200

中國

1

前500

圣經(jīng)

1

前250

Archimedes

3

263

劉徽

5

480

祖沖之

7

1429

Al-Kashi

14

1593

Romanus

15

1596

Ludolph Van Ceulen

20

1609

Ludolph Van Ceulen

35

1699

Sharp

71

1706

John Machin

100

1719

De Lagny

127(112位正確)

1794

Vega

140

1824

Rutherford

208(152位正確)

1844

Strassnitzky & Dase

200

1847

Clausen

248

1853

Lehmann

261

1853

Rutherford

440

1874

William Shanks

707(527位正確)

20世紀后

紀錄創(chuàng)造者

所用機器

小數(shù)點后位數(shù)

1946

 

Ferguson

 

620

1947

1

Ferguson

 

710

1947

9

Ferguson & Wrench

 

808

1949

 

Smith & Wrench

 

1,120

1949

 

Reitwiesner et al

ENIAC

2,037

1954

 

Nicholson & Jeenel

NORC

3,092

1957

 

Felton

Pegasus

7,480

1958

1

Genuys

IBM 704

10,000

1958

5

Felton

Pegasus

10,021

1959

 

Guilloud

IBM 704

16,167

1961

 

Shanks & Wrench

IBM 7090

100,265

1966

 

Guilloud & Filliatre

IBM 7030

250,000

1967

 

Guilloud & Dichampt

CDC 6600

500,000

1973

 

Guilloud & Bouyer

CDC 7600

1,001,250

1981

 

Miyoshi & Kanada

FACOM M-200

2,000,036

1982

 

Guiloud

 

2,000,050

1982

 

Tamura

MELCOM 900II

2,097,144

1982

 

Tamura & Kanada

HITACHI M-280H

4,194,288

1982

 

Tamura & Kanada

HITACHI M-280H

8,388,576

1983

 

Kanada, Yoshino & Tamura

HITACHI M-280H

16,777,206

1983

10

Ushiro & Kanada

HITACHI S-810/20

10,013,395

1985

10

Gosper

Symbolics 3670

17,526,200

1986

1

Bailey

CRAY-2

29,360,111

1986

9

Kanada & Tamura

HITACHI S-810/20

33,554,414

1986

10

Kanada & Tamura

HITACHI S-810/20

67,108,839

1987

1

Kanada, Tamura & Kubo et al

NEC SX-2

134,217,700

1988

1

Kanada & Tamura

HITACHI S-820/80

201,326,551

1989

5

Chudnovskys

CRAY-2 & IBM-3090/VF

480,000,000

1989

6

Chudnovskys

IBM 3090

525,229,270

1989

7

Kanada & Tamura

HITACHI S-820/80

536,870,898

1989

8

Chudnovskys

IBM 3090

1,011,196,691

1989

11

Kanada & Tamura

HITACHI S-820/80

1,073,741,799

1991

8

Chudnovskys

 

2,260,000,000

1994

5

Chudnovskys

 

4,044,000,000

1995

8

Takahashi & Kanada

HITACHI S-3800/480

4,294,967,286

1995

10

Takahashi & Kanada

 

6,442,450,938

1997

7

Takahashi & Kanada

 

51,539,600,000

1999

4

Takahashi & Kanada

 

68,719,470,000

1999

9

Takahashi & Kanada

HITACHI SR8000

206,158,430,000

圓周率的最新計算紀錄

1、新世界紀錄

圓周率的最新計算紀錄由兩位日本人Daisuke Takahashi和Yasumasa Kanada所創(chuàng)造。他們在日本東京大學的IT中心,以Gauss-Legendre算法編寫程序,利用一臺每秒可執(zhí)行一萬億次浮點運算的超級計算機,從日本時間1999年9月18日19:00:52起,計算了37小時21分04秒,得到了圓周率的206,158,430,208(3*236)位十進制精度,之后和他們于1999年6月27日以Borwein四次迭代式計算了46小時得到的結(jié)果相比,發(fā)現(xiàn)最后45位小數(shù)有差異,因此他們?nèi)⌒?shù)點后206,158,430,000位的?值為本次計算結(jié)果。這一結(jié)果打破了他們于1999年4月創(chuàng)造的68,719,470,000位的世界紀錄。

2、最后20位

圓周率小數(shù)點后206,158,430,000位的最后20位為:

22144 96687 55157 30964

3、 p小數(shù)點后2000億位中各數(shù)字出現(xiàn)的次數(shù):

0 : 20000030841 1 : 19999914711

2 : 20000136978 3 : 20000069393

4 : 19999921691 5 : 19999917053

6 : 19999881515 7 : 19999967594

8 : 20000291044 9 : 19999869180

4、一些有趣的數(shù)字序列在p小數(shù)點后出現(xiàn)的位置

數(shù)字序列

出現(xiàn)的位置

01234567891

26,852,899,245

41,952,536,161

99,972,955,571

102,081,851,717

171,257,652,369

01234567890

53,217,681,704

148,425,641,592

432109876543

149,589,314,822

543210987654

197,954,994,289

98765432109

123,040,860,473

133,601,569,485

150,339,161,883

183,859,550,237

09876543210

42,321,758,803

57,402,068,394

83,358,197,954

10987654321

89,634,825,550

137,803,268,208

152,752,201,245

27182818284

45,111,908,393

PC機上的計算

1、PiFast

目前PC機上流行的最快的圓周率計算程序是PiFast。它除了計算圓周率,還可以計算e和sqrt(2)。PiFast可以利用磁盤緩存,突破物理內(nèi)存的限制進行超高精度的計算,最高計算位數(shù)可達240億位,并提供基于Fabrice Bellard公式的驗算功能。

2、 PC機上的最高計算記錄

最高記錄

12,884,901,372位

時間

2000年10月10日

記錄創(chuàng)造者

Shigeru Kondo

所用程序

PiFast ver3.3

機器配置

Pentium III 1G, 1792M RAM,WindowsNT4.0,40GBx2(IDE,FastTrak66)

計算時間

1,884,375秒 (21.8天)

驗算時間

29小時

圓周率小數(shù)點后1000位

1415926535 8979323846 2643383279 5028841971 6939937510 5820974944 5923078164 0628620899 8628034825 3421170679 8214808651 3282306647 0938446095 5058223172 5359408128 4811174502 8410270193 8521105559 6446229489 5493038196 4428810975 6659334461 2847564823 3786783165 2712019091 4564856692 3460348610 4543266482 1339360726 0249141273 7245870066 0631558817 4881520920 9628292540 9171536436 7892590360 0113305305 4882046652 1384146951 9415116094 3305727036 5759591953 0921861173 8193261179 3105118548 0744623799 6274956735 1885752724 8912279381 8301194912 9833673362 4406566430 8602139494 6395224737 1907021798 6094370277 0539217176 2931767523 8467481846 7669405132 0005681271 4526356082 7785771342 7577896091 7363717872 1468440901 2249534301 4654958537 1050792279 6892589235 4201995611 2129021960 8640344181 5981362977 4771309960 5187072113 4999999837 2978049951 0597317328 1609631859 5024459455 3469083026 4252230825 3344685035 2619311881 7101000313 7838752886 5875332083 8142061717 7669147303 5982534904 2875546873 1159562863 8823537875 9375195778 1857780532 1712268066 1300192787 6611195909 2164201989

圓周率的探索者

 

Archimedes (BC287 - BC212)

 祖沖之 (430 - 501)

 Ludolph van Ceulen (1540 - 1610)

 John Machin(1680 - 1751)

 Johann Heinrich Lambert (1728 - 1777)

 Adrien-Marie Legendre (1752 - 1833)

 Johann Carl Friedrich Gauss (1777 - 1855)

 Carl Louis Ferdinand von Lindemann (1852 - 1939)

 Srinivasa Ramanujan(1887 - 1920)

 ENIAC (1946)

 David & Gregory Chudnovsky

 David H. Bailey

 Jonathan M. Borwein

 Peter Borwein

 Simon Plouffe

 Fabrice Bellard (1973)


 Yasumasa Kanada

從左至右:Eugene Salamin,Yasumasa Kanada, David H. Bailey,William Gosper


本文來自:逍遙右腦記憶 http://www.yy-art.cn/gaozhong/196179.html

相關閱讀:高一數(shù)學公式 圓的標準方程和一般方程_高中數(shù)學公式