数学软件作业圆周率的计算.doc

2007数学应用软件设计内容 ——怎样计算圆周率π的值 实验报告 姓名：。。 学号：；；；； 实验目的： 利用所学的数学应用软件，解决下列问题： 1. 用反正切函数的幂级数展开式结合有关公式求π，若要精确到以40位、50位数字，试比较简单公式和Machin公式所用的项数； 2. 用数值积分计算π，分别用梯形法和Simpson法精确到10位数字，用Simpson法精确到15位数字； 3. 用Monte Carlo 法计算π，除了加大随机数，在随机数一定时可重复算若干次后求平均值，看能否求得5位精确数字？ 4. 设计方案用计算机模拟Buffon实验； 5.利用学习过的知识(或查阅资料)，提出其他计算π的方法。 实验的基本理论与方法： 1.利用反正切函数的幂级数表达式： 利用积分公式两边积分得： 在上式中令x=1得： 用Mathematica计算In[1] k=1000; S1=N[4*Sum[(-1)^(n-1)/(2n-1),{n,1,k}],18] [Out2] 3.14059265383979293 In[3] k=10000; [Out4] 3.14149265359004324 In[5] k=15000; [Out6] 3.14152598692320065 In[7] k=20000 [Out8] 3.14154265358982449 由上述计算可以知道，用反正切函数的幂级数展开来计算π的值，运算速度比较慢，结果精度也不高,因此此方法一般不选取。 下面对精确到相同精度时简单公式和Machin公式所用的项数进行比较： 简单公式： Machin公式： 再用Mathematica计算： 当精确到40位时： 对于简单公式： Clear[k,n,S] In[9] k =62; s2 = N[4*Sum[(-1)^(n - 1)*(1/2)^(2n - 1)/(2n - 1) + (-1)^(n - 1)*(1/3)^(2n - 1)/(2n - 1), {n, 1, k}], 40] [Out10] 3.141592653589793238462643383279502884197 对于Machin公式： In[11] k = 28; s3 = N[4*Sum[ 4*(-1)^(n - 1)*(1/5)^(2n - 1)/(2n - 1) - (-1)^(n - 1)*(1/239)^(2n - 1)/(2n - 1), {n, 1, k}], 40] [Out12] 3.141592653589793238462643383279502884197 当精确到50位时： 对于简单公式： In[13] k = 78; s2 = N[4*Sum[(-1)^(n - 1)*(1/2)^(2n - 1)/(2n - 1) + (-1)^(n - 1)*(1/3)^(2n - 1)/(2n - 1), {n, 1, k}], 40] [Out14] 3.1415926535897932384626433832795028841971693993750 对于Machin公式： In[15] k = 35; s3 = N[4*Sum[ 4*(-1)^(n - 1)*(1/5)^(2n - 1)/(2n - 1) - (-1)^(n - 1)*(1/239)^(2n - 1)/(2n - 1), {n, 1, k}], 50] [Out16] 3.1415926535897932384626433832795028841971693993751 由上述比较可以看出：在精确到相同精度时，用Machin公式计算的结果比较准确，并且所需的项数也比较少，当精确到40位时简单公式和Machin公式分别需要62和28位，而当精确到50位时，分别为78和35位。所以用Machin公式计算圆周率速度较快，结果也较精确，并且程序也较简单，是一种比较好的方法。 用数值积分计算π： 将[0,1]区间n等分，取xk=k/n, yk= 1/ (1+xk2) 用Mathematica计算： 梯形法（精确到10位）： In[17] y[x_]:=4/(1+x^2); n=20000; S4=N[