function f81(n)
x=unifrnd(-2*pi,2*pi,1,n); y=(1-x.^2).*sin(3*x); max(y)
x=-2*pi:0.001:2*pi; y=(1-x.^2).*sin(3*x); plot(x,y) xlabel('x'); ylabel('y'); 运行结果: >>f81(1000) ans =
32.3293
>> f81(10000) ans =
32.4002
>> f81(100000) ans =
32.4006
做出函数的图像,并且标出最高点的值
结果分析:可以看到,蒙特卡罗法求出的最大值接近于32.4,而从图中可以看出最大值是32.33,求出的结果比较符合。 II、(2)使用均值估计法 分析:由于x1=x2+10,所以可以消元,使其变为两个自变量x2和x3。x2,x3在它们被允许的范围内生成多个随机的数值,利用max函数可以近似地求出结果。然后做出图像,进行结果的比较。 程序:
function f82(n)
x2=unifrnd(10,20,1,n); x1=10+x2;
x3=unifrnd(-10,20,1,n); fori=1:n
if -x1(i)+2*x2(i)+2*x3(i)>=0 if x1(i)+2*x2(i)+2*x3(i)<=72 y(i)=(x1(i))*(x2(i))*(x3(i)); end end end max(y)
x2=10:0.1:20; x3=-5:21/100:16;
[X,Y]=meshgrid(x2,x3); err1 = X+2*Y<10; err2 = 3*X+2*Y>62; X(err1) = nan; Y(err2) = nan; Z=X.*Y.*(X+10); surf(X,Y,Z) 运行结果: >>f82(1000) ans =
3.3889e+03 >> f82(10000) ans =
3.4357e+03 >> f82(100)
ans =
3.3726e+03 >> f82(100000) ans =
3.4441e+03
结果分析:可以看到,蒙特卡罗法求出的最大值接近于3400,而从图中可以看出最大值是3437,求出的结果比较符合。 II、(3)使用蒙特卡罗法
分析:x,y在它们被允许的范围内生成多个随机的数值,利用max函数可以近似地求出结果。然后做出图像,进行结果的比较。 程序:
function f83(n)
x=unifrnd(-1.5,1.5,1,n); y=unifrnd(-1.5,1.5,1,n);
z=(x.^2+2*(y.^2)+x.*y).*exp(-x.^2-y.^2); max(z)
x=-1.5:0.1:1.5; y=-1.5:0.1:1.5;
[X,Y]=meshgrid(x,y);
Z=(X.^2+2*(Y.^2)+X.*Y).*exp(-X.^2-Y.^2);
surf(X,Y,Z) 运行结果: >>f83(1000) ans =
0.8105 >>f83(10000) ans =
0.8117
作出函数图,并且标出最大值
结果分析:可以看到,蒙特卡罗法求出的最大值接近于0.81,而从图中可以看出最大值是0.8025,求出的结果比较符合。
3.实验总结和实验感悟
这次蒙特卡洛法令我印象比较深刻,特别是可以利用多次模拟实验的方法来求圆周率,这是我以前没有接触过的。蒙特卡洛法可以理解成一种思想,就是多次随机的实验来求近似值。不过这种方法比较适合电脑模拟,模拟次数足够高才可以保证误差不过大,而且某些可以直接求解的问题并不需要用蒙特卡罗法来做。