電子科技大學(xué)數(shù)學(xué)實(shí)驗(yàn)實(shí)驗(yàn)五 Monte Carlo方法計(jì)算Viviani體體積

1、1/12實(shí)驗(yàn)五實(shí)驗(yàn)五:蒙特卡羅方法實(shí)驗(yàn)蒙特卡羅方法實(shí)驗(yàn) 面積、體積計(jì)算問(wèn)題面積、體積計(jì)算問(wèn)題Viviani體的體積計(jì)算體的體積計(jì)算思考題與練習(xí)題思考題與練習(xí)題2/12蒙特卡羅方法蒙特卡羅方法隨機(jī)投點(diǎn)試驗(yàn)求近似解隨機(jī)投點(diǎn)試驗(yàn)求近似解引例引例. 給定曲線(xiàn)給定曲線(xiàn)y =2 x2 和曲線(xiàn)和曲線(xiàn)y3 = x2，曲線(xiàn)的交點(diǎn)為曲線(xiàn)的交點(diǎn)為：P1( 1，1 )、P2( 1，1 )。曲線(xiàn)圍成平面有限區(qū)域，用曲線(xiàn)圍成平面有限區(qū)域，用蒙特卡羅方法計(jì)算區(qū)域面積。蒙特卡羅方法計(jì)算區(qū)域面積。 P=rand(10000,2);x=2*P(:,1)-1;y=2*P(:,2);II=find(y=x.2);M=length(I

2、I);S=4*M/10000plot(x(II),y(II),g.)S = 2.11361532)2(11322dxxxS3/12例例5.14 計(jì)算計(jì)算 其中其中D為為y= x 2與與y2 = x 所圍所圍 D的邊界曲線(xiàn)交點(diǎn)為的邊界曲線(xiàn)交點(diǎn)為：(1，-1)，(4，2)，被積函數(shù)在求被積函數(shù)在求積區(qū)域內(nèi)的最大值為積區(qū)域內(nèi)的最大值為1616。積分值是三維體積，該三維積分值是三維體積，該三維圖形位于立方體區(qū)域圖形位于立方體區(qū)域0 x 4，1 y 2，0 z 16內(nèi)，立方體區(qū)域的體積為內(nèi)，立方體區(qū)域的體積為1 192。 data=rand(10000,3); x=4*data(:,1); y=-1+3

3、*data(:,2); z=16*data(:,3); II=find(x=y.2&x=y+2&zsqrt(x.2+y.2)&z=1&u=x.2+y.2+z.2); M=length(II); V=8*M/100005/12實(shí)驗(yàn)實(shí)驗(yàn):蒙特卡羅方法計(jì)算蒙特卡羅方法計(jì)算Viviani體體積體體積( , , )|02, 11,02x y zxyz x=2*rand產(chǎn)生產(chǎn)生0到到2之間的隨機(jī)數(shù)之間的隨機(jī)數(shù)y=2*rand-1產(chǎn)生產(chǎn)生-1到到1之間的隨機(jī)數(shù)之間的隨機(jī)數(shù)z=2*rand產(chǎn)生產(chǎn)生0到到2之間的隨機(jī)數(shù)之間的隨機(jī)數(shù)含于體積含于體積 = 8 的六面體的六面體NmV8

4、由于由于rand 產(chǎn)生產(chǎn)生0 到到1之間的隨機(jī)數(shù)之間的隨機(jī)數(shù),所以所以N個(gè)點(diǎn)均勻分布于六面體中個(gè)點(diǎn)均勻分布于六面體中, Viviani體中體中占有占有m個(gè)個(gè),則錐體與六面體體積之比近似為則錐體與六面體體積之比近似為 m : N222222&2xyzRxyRxR ，2222xyzR22xyRx 求球體求球體 被圓柱面被圓柱面所割下部分的體積所割下部分的體積 ( 稱(chēng)為維維安尼稱(chēng)為維維安尼 (Viviani) 體體 ). 解解 由所求立體的對(duì)稱(chēng)性，只要求出在第由所求立體的對(duì)稱(chēng)性，只要求出在第 一卦限內(nèi)的部分體積一卦限內(nèi)的部分體積, ,再乘以再乘以4, ,即得所求立體的體即得所求立體的體 yx

5、zR(圖圖21-32), ,而曲頂?shù)姆匠虨槎數(shù)姆匠虨?222.zRxy 所以所以 2224d,DVRxy 222223042.323RRx xRx xVRxy dydxR 0y 22xyRxxy 平面內(nèi)由平面內(nèi)由和和所確定的區(qū)域所確定的區(qū)域 D 積積. 在第一卦限內(nèi)的立體是一個(gè)曲頂柱體在第一卦限內(nèi)的立體是一個(gè)曲頂柱體, 其底為其底為 其中其中 22( ,)0,.Dx yyxyRx4)2(222RyRx8/12syms x y;f=sqrt(4-x2-y2);y1=-sqrt(2*x-x2);y2=sqrt(2*x-x2);S1=int(f,y,y1,y2);S2=int(S1,x,0,2)

6、V=double(S2)P=rand(10000,3);X=2*P(:,1);Y=2*P(:,2)-1;Z=2*P(:,3);II=find(X-1).2+Y.2=1&Z=sqrt(4-X.2-Y.2);V1=8*length(II)/10000Error=V-V1實(shí)驗(yàn)參考程序?qū)嶒?yàn)參考程序 蒙特卡羅方法計(jì)算體積蒙特卡羅方法計(jì)算體積9/12實(shí)驗(yàn)任務(wù)一實(shí)驗(yàn)任務(wù)一：記錄：記錄L次實(shí)驗(yàn)的實(shí)驗(yàn)數(shù)據(jù)及誤差次實(shí)驗(yàn)的實(shí)驗(yàn)數(shù)據(jù)及誤差實(shí)驗(yàn)任務(wù)二實(shí)驗(yàn)任務(wù)二：修改實(shí)驗(yàn)程序修改實(shí)驗(yàn)程序MonteC計(jì)算計(jì)算L次實(shí)驗(yàn)數(shù)次實(shí)驗(yàn)數(shù)據(jù)均值及均值誤差據(jù)均值及均值誤差( mean 計(jì)算平均值計(jì)算平均值 )序號(hào)序號(hào) 1 2 3

7、 4 5 6 7 數(shù)據(jù)數(shù)據(jù)誤差誤差 L 8 16 32 64 128 256均值均值誤差誤差10/12圓柱面在圓柱面在 平面的投影用極坐平面的投影用極坐標(biāo)表示為標(biāo)表示為theta=(-180:10:1800)*pi/180;r=0:0.1:1;X=r*(1+cos(theta);Y=r*sin(theta);Z=sqrt(4-X.2-Y.2);meshz(X,Y,Z)axis offcolormap(0 0 0)view(130,30)Viviani體繪制程序體繪制程序xy0.5 (1cos ),0.5 sin ,0,2 xRyR 11/12思考題與練習(xí)題思考題與練習(xí)題5. 說(shuō)明說(shuō)明L次蒙特卡

8、羅實(shí)驗(yàn)平均值計(jì)算次蒙特卡羅實(shí)驗(yàn)平均值計(jì)算Viviani的體積誤的體積誤差與實(shí)驗(yàn)次數(shù)之間關(guān)系差與實(shí)驗(yàn)次數(shù)之間關(guān)系1.蒙特卡羅方法計(jì)算面積和計(jì)算體積方法有何差異？蒙特卡羅方法計(jì)算面積和計(jì)算體積方法有何差異？2.說(shuō)明蒙特卡羅方法計(jì)算說(shuō)明蒙特卡羅方法計(jì)算Viviani的體積誤差與哪些因的體積誤差與哪些因素有關(guān)。素有關(guān)。3.概率論中的貝努里大數(shù)定律在本次實(shí)驗(yàn)中體現(xiàn)如何概率論中的貝努里大數(shù)定律在本次實(shí)驗(yàn)中體現(xiàn)如何4.敘述概率論中著名獨(dú)立同分布大數(shù)定律，并以這一敘述概率論中著名獨(dú)立同分布大數(shù)定律，并以這一大數(shù)定律解釋實(shí)驗(yàn)的數(shù)據(jù)變換規(guī)律。大數(shù)定律解釋實(shí)驗(yàn)的數(shù)據(jù)變換規(guī)律。12/125.下面程序繪出二維圖形填下面程序繪出二維圖形填充圖充圖(右圖右圖)。分析每條語(yǔ)句功。分析每條語(yǔ)句功能能,給程序中語(yǔ)句寫(xiě)注記給程序中語(yǔ)句寫(xiě)注記x1=-1:0.