无理数的近似计算

1、西安交通大学数学实验报告实验题目: 无理数的近似计算班级 姓名 学号 填写日期 2013年6月7日0 / 6一、实验问题问题1.对P89页示例4：利用定积分0111+x2dx=4计算圆周率的近似值利用矩形、和抛物线形方法求近似值，并与书上的方法进行比较，看哪种方法收敛速度快?问题2.设有一制作均匀的冰淇淋可以看成由圆锥面z=x2+y2和球面x2+y2+(z-1)2=1围成，采用取随机数的方法，用蒙特卡罗法计算这个冰淇淋的体积。二、问题分析问题1.矩形公式 ：abfxdxi=1nfxixi抛物线形公式：abfxdx=i=1nx2i-1x2i+1fxdxb-a6ni+1n(y2i-1+4y2i+y

2、2i+1)问题2.根据“冰淇淋”的形状容易发现由平面x=-1,x=1,y=-1,y=-1,z=0,z=2围成的正方体恰好能将其包围在其中。假设冰淇淋的体积为V，正方体的体积为U试想：可以在正方体区域内等概率地随意各处取点，所取点落入冰淇淋的概率应该是冰淇淋的体积与正方体的体积之比，即VU。若总取点数为n，其中落入冰淇淋中的次数为cs，则有csn=VU,于是V=Ucsn，其中U=8。三、程序设计问题1.矩形公式：n=input('请输入等分积分区间数');i=0:1/n:1;s=0;for k=1:length(i)-1 s=s+(1/(1+(i(k)2)*1/n;endpai=

3、4*s 抛物线形公式：n=input('请输入等分积分区间数(n必须为偶数)');i=0:1/n:1;s=0;for k=1:(length(i)-1)/2 s=s+1/(1+(i(2*k-1)2)+4*(1/(1+(i(2*k)2)+1/(1+(i(2*k+1)2);endpai=4/(3*n)*s问题2.cs=0;n=input('请输入随机取点数');for i=1:n; a=rand(1,3); x=2*a(1)-1;y=2*a(2)-1;z=2*a(3); if z-sqrt(x.2+y2)>0 & x2+y2+(z-1)2<1 c

4、s=cs+1; endendV=8*cs/n四、实验结果和结论问题1. 三种方法收敛速度比较如下等分区间数n圆周率的近似值矩形公式法抛物线形公式法书上所用方法103.2399259889071593.1415926139392153.142425985001098203.1911759869541283.1415926529697863.141800986893093503.1615259869232543.1415926535872543.1416259869230041003.1515759869231273.1415926535897553.1416009869231235003.143

5、5919869231283.1415926535897973.14159298692312410003.1425924869231223.1415926535897923.14159273692313150003.1417926469231213.1415926535898053.141592656923131结论:抛物线形公式法收敛速度最快，书上所用方法次之，矩形公式法收敛速度最慢。问题2. 令随机取点数取为100000000，运行结果如下五、实验总结和体会完成一次数学实验需要各方面的因素，其中最重要的因素是编写出一套没有谬误的程序，而要编写好程序，就需要对matlab中必要的命令了解且能应

6、用。对于matlab命令的不够熟悉，在实验过程中造成了一些麻烦。例如问题1中，圆周率的近似值在书上精确到了小数点后14位，但我运行程序时却只能显示到小数点后4位，这就使通过抛物线形公式法计算出的圆周率近似值从等分区间数为10时到等分区间数为5000时都是相同的，无从比较其收敛速度。于是只好求助于百度知道，得知可以通过matlab菜单栏这样设置：file=>Preferences=>Command Window=>Text display，在Numeric format选择long。如果对vpa命令足够熟悉，就可以在此应用vpa命令来控制结果的精度。而在问题2中，我最大的收获在于其对问题新颖的观察角度和解决方法给了我一种新的思想。其次，还学会了“&”符号的应用。六、实验扩展1.通过将所求数用积分式表示，再转化为多项式求和的方法不仅可以解决求的问题，也可以推广到其他能够用积分式表示的无理数的求解问题上，如求无理数e和欧拉数。2.蒙特卡罗法为估计无理数的值提供了一种新思路，但是其误差太大，只有随机取点数非常大时，才能得到一个比较接近真实值的解。例如实验中取了1亿个随机点，程序运行了很长时间，才只能