如何计算π的值(MATLAB).doc

如何计算得值
1、蒙特卡罗（Mｏnｔe Ｃａｒｌo)法
思想:
　取一正方形A，以A得一个顶点为圆心，A得边长为半径画圆，取四分之一圆（正方形内得四分之一圆)为扇形B。已知A得面积,只要求出B得面积与A　得面积之比,就能得出,再由B得面积为圆面积得四分之一,利用公式即可求出得值。因此，我们得目得就就是要找出得值。
可以把Ａ与B瞧成就是由无限多个点组成,而B内得所有点都在A内。随机产生个点,若落在B内得有个点(假定A得边长为１,以扇形圆心为坐标系原点。则只要使随机产生横纵坐标、满足得点，就就是落在B内得点），则可近似得出得值,即,由此就可以求出得值。
程序(1）:
i＝１;m＝０;n=1000；
for i=１：ｎ
　　 ａ＝ranｄ(１,2);　
if a(1）^2＋a(2)^2〈＝１
　 ｍ=m+1;
　　 end
eｎｄ
p=vpa(４*m/ｎ，30)
程序运行结果:
ｐ =
３、14０000０
2、泰勒级数法
思想:
反正切函数得泰勒级数展开式为:

将代入上式有
、
利用这个式子就可以求出得值了。
程序（2):
i=1;n=100０;s=0；
fｏｒ i=1：n
　 　ｓ=ｓ+(-１)^（i—１）／(2*i-１)；
end
p=vｐa（4＊ｓ,30)
程序运行结果:
ｐ ＝
3、1４94126
当取得值为１0000时,就会花费很长时间,而且精度也不就是很高。原因就是时,得展开式收敛太慢。因此就需要找出一个使得收敛更快。
若取,则我们只有找出与得关系,才能求出得值。
令,，
根据公式有,则有。
所以可以用来计算得值。
程序（)：
ｉ=１;ｎ＝1000；ｓ＝0;s1=0;ｓ２=0;
ｆoｒ　ｉ=1：n
　 s1=ｓ1+（—1）＾（i—1）*（１/２)^（２＊ｉ-１）/（２*i－1);
　　s2=s２＋（—1)^（i-１)＊(1/3)^（2*ｉ—1)/(2*ｉ-１)；
enｄ
s=s１＋s2;
ｐ=vpa(4＊ｓ,30）
程序运行结果：
p =
3、1493２38４６2６4３3８3２8
显然,级数收敛越快，取同样得值可以得到更高得精度.以同样得方法,能得出,程序与上面得一样。这样得近似值可以精确到几百位。
3、数值积分法
思想:
半径为1得圆得面积就是。以圆心为原点建立直角坐标系，则圆在第一象限得扇形就是由与轴,轴所围成得图形，扇形得面积。只要求出扇形得面积,就可得出得值。而扇形面积可近似等于定积分得值。
对于定积分得值,可以瞧做成曲线与轴,,所围得曲边梯形得面积。把分成等分，既得个点,，,组成个小区间,每一个小区间与轴,,所围成得图形就是一个小曲边梯形。而梯形得面积计算公式就是,对于第个小曲边梯形有上底为，下底为。所有小梯形得高都为。所以第个小曲边梯形得面积为。曲边梯形得总面积即定积分得值就就是所有小梯形得面积总与。
为了避免根号,我们也可以利用积分

得出得值。
我们可以利用对求曲边梯形得面积来得出定积分得值,从而得出得值.
程序(3)：
a=0;b=1；ｓ=0;n=100０；i＝0；
h=（b-a)/n；
ｆor　i=０:(n—１)
　 　xi＝ａ＋ｉ*h；
　　 yi=1/(1＋(ｘi)^２)；
　　 xj＝a+(i+１）*h；
ｙｊ＝1/（1+(xj）^2);
　 ｓ=ｓ+(ｙi+yj)*ｈ／2；　
eｎｄ
p＝vpａ（4*ｓ,３０）
程序运行结果:
p =
3、１4277５830259８52２2８
　 对于数值积分法求值，以上程序简洁明了。我们也可以以做循环，用一条语句求出值。
程序(3’）：
s＝0；n=10０0;
fｏr ｘ=０：（１／n)：（1—(1/ｎ）)
　 s=ｓ＋(1/（１+ｘ＾２）+１/（１+(x+(1／n))^2)）＊（1/n）/2;
eｎｄ
p=vpa（4＊ｓ，30)
程序运行结果:
p ＝
3、14277583025985222８
用以上三种方法求,都取1000时，泰勒级数法求，得到得近似值精度最高。