2022年pi计算实验报告（18页）.docx

精选word文档 下载可编辑　　实验报竹　　班级：02　　姓名：张海洋　　学号：9　　圆周率（选做）　　要求：先查资料看看古人是怎样计算 n的，再对n的各种计算 方法进行研究和讨论（收敛速度等），并给出不同算法算出 n　　的小数点后第10000位的数字是什么，你觉得该数字应该是多少　　第一部分:　　圆周率简介　　圆周率是指平面上圆的周长与直径之比 （ratio of the circumfere nee of a circle to the　　diameter）用符号n （读音：p a）i表示中国古代有圆率、圆率、周等名称它是一 个常数（约等于）它是一个无理数，即无限不循环小数在日常生活中，通常都用代 表圆周率去进行近似计算而用十位小数便足以应付一般计算即使是工程师或物理学 家要进行较精密的计算，充其量也只需取值至小数点后几百个位计算圆周率的方法　　“历史上一个国家所算得的圆周率的准确程度， 可以作为衡量这个国家当时数学发展　　水平的指标　　历史上最马拉松式的计算，其一是德国的 Ludolph Van Ceulen他几乎耗尽了一生的　　时间，计算到圆的内接正262边形，于1609年得到了圆周率的35位精度值，以至于圆 周率在德国被称为Ludolph数；其二是英国的 William Shanks,他耗费了 15年的光阴，在 1874年算出了圆周率的小数点后707位。

可惜，后人发现，他从第528位开始就算错了　　把圆周率的数值算得这么精确，实际意义并不大现代科技领域使用的圆周率值，有十 几位已经足够了如果用Ludolph Van Ceulen算出的35位精度的圆周率值，来计算一个 能把太阳系包起来的一个圆的周长，误差还不到质子直径的百万分之一以前的人计算 圆周率，是要探究圆周率是否循环小数自从 1761年Lambert证明了圆周率是无理数，　　1882年Lin dema nn证明了圆周率是超越数后，圆周率的神秘面纱就被揭开了　　在中国，公元263年前后，刘徽提出著名的 “割圆术”求出了比较精确的圆周率　　他发现：当圆内接正多边形的边数不断增加后，多边形的周长会越来越逼近圆周长，而 多边形的面积也会越来越逼近圆面积于是，刘徽利用正多边形面积和圆面积之间的关 系，从正六边形开始，逐步把边数加倍：正十二边形、正二十四边形，正四十八边形……， 一直到正三七二边形，算出圆周率等于三点一四一六，将圆周率的精度提高到小数点 后第四位在刘徽研究的基础上，祖冲之进一步地发展，经过既漫长又烦琐的计算，一 直算到圆内接正24576边形，而得到一个结论： V n V 同时得到n 的两个近似分　　2　　2、　　数：约率为22 7;密率为355 113。

他算出的n的8位可靠数字，不但在当时是最 精密的圆周率，而且保持世界记录九百多年以致于有数学史家提议将这一结果命名为　　“祖率” 现在的人计算圆周率，多数是为了验证计算机的计算能力，还　　有，就是为了兴趣　　第二部分：　　古人计算圆周率方法　　古人计算圆周率，一般是用割圆法即用圆的内接或外切正多边形来逼近圆 的周长Archimedes用正96边形得到圆周率小数点后3位的精度；刘徽用正3072 边形得到5位精度；Ludolph Van Ceulen用正262边形得到了 35位精度这种基 于几何的算法计算量大，速度慢，吃力不讨好随着数学的发展，数学家们在进 行数学研究时有意无意地发现了许多计算圆周率的公式 下面挑选一些经典的常　　116 arctg 亏14arg tg 239　　1　　16 arctg 亏　　1　　4arg tg 239　　Machin公式　　2 n 1　　3 5 7　　arctg(X)X； X5 X7 L　　(1) n 1 x　　2n　　这个公式由英国天文学教授 John Machin于1706年发现他利用这个公式计算到了 100位的圆周率　　Machin公式每计算一项可以得到位的十进制精度。

因为它的计算过程中被乘数和被除数都不大于长整　　数，所以可以很容易地在计算机上编程实现　　还有很多类似于 Machin公式的反正切公式在所有这些公式中， Machin公式似乎是最快的了虽然　　如此，如果要计算更多的位数，比如几千万位， Machin公式就力不从心了下面介绍的算法，在 PC　　机上计算大约一天时间，就可以得到圆周率的过亿位的精度这些算法用程序实现起来比较复杂因 为计算过程中涉及两个大数的乘除运算，要用 FFT(Fast Fourier Transform)算法FFT可以将两个大数的　　乘除运算时间由 0(n2)缩短为0(nlog(n))　　Ramanujan 公式　　Ramanujan 公式　　9801　　1914年，印度数学家Srinivasa Ramanujan在他的论文里发表了一系列共 14条圆周率的计算公式， 这是　　其中之一这个公式每计算一项可以得到 8位的十进制精度1985年Gosper用这个公式计算到了圆周　　率的 17,500,000 位　　3、　　AGM(Arithmetic-Geometric Mea n)算法　　Gauss-Legendre 公式：　　初值：　　a+6重复计算：　　a+6　　STWHHE HFW　　这个公式每迭代一次将得到双倍的十进制精度，比如要计算 100万位，迭代20次就够了'　　1999年9月Takahash和Kanada用这个算法计算到了圆周率的 206,158,430,000位，创出 新的世界纪录。

　　4、　　Borwe in　　四次迭代式：　　初值：　　'- : 二_ ： 1 重复计算:　　sr泸 1　　飞?厂1+（卜铲）"4　　(1}' 1 -　　最后计算： ■!，　　这个公式由Jonathan Borwein和Peter Borwein于1985年发表，它四次收敛于圆周率　　5、　　Bailey-Borwein-Plouffe 算法　　丄（丄o　　丄（丄　　o16n 8n 1　　2 1 1　　8n 4 8n 5 8n 6　　这个公式简称 BBP公式，由David Bailey, Peter Borwein和Simon Plouffe于1995年共同发 表它打破了传统的圆周率的算法，可以计算圆周率的任意第 n位，而不用计算前面的n-　　1位这为圆周率的分布式计算提供了可行性 1997年，Fabrice Bellard找到了一个比BBP　　32 G 二f快　　32 G 二f　　滋麗一一璽—八 葆一 或I 10?+1 ~ ib i+S " 10n+o 10n+7 10m+9　　第三部分：　　对于n的几种计算的研究和讨论:　　1、　　Comma nd window　　? n=lCOOO;　　y=sq.rt (1 一貨."2);　　h=l n;　　a(ri3=vps r aps (y.) ,11)　　ans=　　3,5　　n=10　　ans　　n=20　　ans　　n=50　　ans　　n=100　　ans　　n=200　　ans　　n=500　　ans　　n=1000　　ans　　n=2000　　ans　　半径为1的圆称为单位圆，它的面积等于 。

只要计算出单位圆的面积，就算出了　　在坐标轴上画出以圆点为圆心，以1为半径的单位圆（如下图），则这个单位圆在第 一象限的部分是一个扇形，而且面积是单位圆的 1 4,于是，我们只要算出此扇形的面　　积，便可以计算出 但计算量较大　　2、　　数值积分法（II）利用公式　　1 1　　4 dx　　01 x2　　Command Window　　Command Window　　X　　X　　Command Window　　Command Window　　X　　X　　? n=50:　　i=a； i r； i； s=0;　　Tou k=i; 1 enetli ti)-1 s=s+(l (l-K((i(i)+i(k-kl)V2) ?nd　　arts=　　3, 14162536692300　　n=10　　ans　　n=20　　ans　　n=50　　ans　　n=100　　ans　　n=200　　ans　　n=500　　ans　　n=1000　　ans　　n=2000 ans=;　　设分点Xi,X2「Xn-1将积分区间[0,1]分成n等分　　所有的曲边梯形的宽度都是h=1 n记yi=f(xi).则第i个曲边梯形的面积A近似地等于梯形 面积，即：　　A=(y(i-1)+yi)h 2　　将所有这些梯形面积加起来就得到：　　A~ 2n[ 2(yi+y2+…yn-i)+y0+yn]　　3、　　利用复化梯形算法求Pi的近似值.　　1 40 1 X　　1 4　　0 1 X2　　=1 2 (Tn+h　　n　　h \　　f a　　i h o　　i 1　　2丿　　Com Window　　yy clear；a-0　　* 1=(b-aj 2* (￡ (a) +f (bD ;　　cr=l;n=l ;　　^uhile er>l. Oe-6 h^(b-a) n.　　5=0；　　for 1 s=s-H(a+i*h-h 2> ;　　end t2=(tHh*s) ?: eu=abs (t2~t 1);　　fprin+f t=K. 6f j r=9t. 6￡\n a 12S er); n=2tn; tl=t2,　　end　　end　　t-3. IOCOOO, r=0.100000　　t=5. 131170J r=C< C31]?G　　t=3.. 138988, r=0. 007912　　：i=3.140942, r=0. 001955　　t=3, 141430, r=0. 00048&　　t=3. 141 啤 r=0. 00C122　　t=3.141B82, r=0.000031　　t=3. 141590, r=C. 000008　　t=3.141592, r=0.000002　　t=3. 141592, r=0. OCOOOO　　?　　4、　　泰勒级数法计算 ：　　3 5　　2k 1　　利用反正切函数的泰勒级数arcta nx x x x　　(1)k1 x。