实验六蒙特卡罗方法.ppt

,实验六:蒙特卡罗方法实验,面积、体积计算问题 冰淇淋锥的体积计算 思考题与练习题 二维填充图绘制方法,rand 产生一个0到 1之间均匀随机数 rand(m,n) 产生mn个0到 1之间均匀随机数,X=rand(10000,1); hist(X),蒙特卡罗方法利用随机试验做近似计算,10000个随机数较均匀地分布在各个小区间上，随机变量X落入小区间的概率仅与小区间长度有关,而与小区间位置无关,unifrnd(a,b) 产生一个a到 b之间均匀随机数,X=unifrnd(2,10) Y=2+(10-2)*rand,例1 计算两条抛物线y = x2, x = y2所围面积.,在正方形0,10,1区域投入2000个均匀随机点 则随机点落入抛物线所围区域的概率为所求面 积与正方形面积之比,function S=area1(N) if nargin=0,N=2000;end X=rand(N,1);Y=rand(N,1); II=find(Y=X.2); m=length(II);S=m/N; x1=0:0.01:1;x2=1:-0.01:0; y1=sqrt(x1);y2=x2.2; fill(x1,x2,y1,y2,c),S = 0.3333,定积分数值计算方法quad()的使用格式,quad(F,a,b) 返回被积函数F(X)从a 到 b的定积分值,F是被积函数名 构成的字符串.相关命令:dblquad() 重积分计算,例2.计算定积分,fun=inline(sqrt(x)-x.2); S=quad(fun,0.01,1) t=0:0.01:1;y=fun(t); fill(0,t,0,y,c) syms u, S0=int(sqrt(u)-u2,0,1),S = 0.3327 S0 = 1/3,例3. 计算下面两条曲线所围区域面积,function S=area2(N) if nargin=0,N=2000;end X=2*rand(N,1)-1; Y=2*rand(N,1); II=find(Y=abs(X); m=length(II);S=4*m/N; x1=0:0.01:1;y1=x1; x2=1:-0.01:-1;y2=1+sqrt(1-x2.2); x3=-1:-0.01:0;y3=-x3; fill(x1,x2,x3,y1,y2,y3,c),S=2.5460,例4. 计算两个半径为1的直交圆柱面所围成体积,x2 + y2 = 1 , x2 + z2 = 1,function V=mlab4(N) if nargin=0,N=2000;end P=rand(N,3); x=P(:,1);y=P(:,2);z=P(:,3); II=find(x.2+y.2=1 meshz(x,y,zz) colormap(0 0 1) axis off view(119,34),V = 16/3,D:,实验:冰淇淋锥的体积计算,&,x=2*rand-1产生 1到1之间的随机数 y=2*rand-1产生 1到1之间的随机数 z=2*rand;产生0到2之间的随机数,冰淇淋锥含于体积 = 8 的六面体,2,2,由于rand 产生0 到1之间的随机数,所以,N个点均匀分布于六面体中,锥体中占有m个,则锥体与六面体体积之比近似为 m : N,function V=icecream(N) if nargin=0,N=10000;end P=rand(N,3); X=2*P(:,1)-1;Y=2*P(:,2)-1;Z=2*P(:,3); R2=X.2+Y.2; R=sqrt(R2); II= find(Z=R,h=2*pi/100;t=0:h:2*pi;r=0:.05:1; x=r*cos(t);y=r*sin(t); z1=sqrt(x.2+y.2); z2=1+sqrt(1+eps-x.2-y.2); mesh(x,y,z1),hold on mesh(x,y,z2) colormap(0 0 1) axis off view(0,-18),思考题与练习题,1. 如何用极坐标变换处理重积分,二维多边形填充图 fill() 使用格式,fill(x,y,c) 用c所指定的颜色对多边形填充.其中,多边形的顶点由x,y确定. c 指定颜色不允许省略;x,y确定的点必须形成封闭的多边形.,2. 圆 x2 + (y 2)2 = 4内挖去小圆x2 + (y 1)2 =1后图形,h=2*pi/100; t=-pi/2:h:3*pi/2; x=cos(t);y=1+sin(t); X=2*cos(t);Y=2+2*sin(t); n=length(x); x=x(n:-1:1);y=y(n:-1:1); fill(x,X,y,Y,c) axis off,大圆,小圆,38/wlxt/ncourse/mathlab/web/shiyan/,function V=icecream(N) if nargin=0,N=10000;end P=rand(N,3); X=2*P(:,1)-1;Y=2*P(:,2)-1;Z=2*P(:,3); R2=X.2+Y.2; R=sqrt(R2); M= sum(Z=R,syms r t x=r*cos(t);y=r*sin(t); z1=sqrt(x.2+y.2); z2=1+sqrt(1+eps-x.2-y.2); f=r*(z2-z1) f=int(f,t,0,2*pi); f=int(f,r,0,1),function V=icecream(N) if nargin=0,N=10000;end m=0; for i=1:N X=-1+2*rand;Y=-1+2*rand;Z=2*rand; R2=X2+Y2; R=sqrt(R2); if Z=R end end V=8*