理论基础:大数定理,当频数足够多时,频率可以逼近概率,从而依靠概率与$pi$的关系,求出$pi$
所以,rand在Monte Carlo中是必不可少的,必须保证测试数据的随机性。
用蒙特卡洛方法进行计算机模拟的步骤:
[1] 设计一个逻辑框图,即模拟模型.
[2] 根据流程图编写程序,模拟随机现象.可通过具有各种概率分布的模拟随机数来模拟随机现象.
[3] 分析模拟结果,计算所需要结果.
ex1.投针试验求$pi$
%蒲丰投针实验的计算机模拟 format long; %设置15位显示精度 a=1; l=0.6; %两平行线间的宽度和针长 figure; axis([0,pi,0,a/2]); %初始化绘图板 set(gca,'nextplot','add'); %初始化绘图方式为叠加 counter=0; n=1120; %初始化计数器和设定投针次数 x=unifrnd(0,a/2,1,n); phi=unifrnd(0,pi,1,n); %样本空间Ω frame=moviein(n); %建立一个1120列的大矩阵 for i=1:nif x(i)
ex2.依然求$pi$
n=10000000; a=2; m=0; for i=1:nx=rand*a; y=rand*a;if ( (x-a/2)^2+(y-a/2)^2 <= (a/2)^2 )m=m+1;end end disp(['投点法近似计算的π为: ',num2str(4*m/n)]);
ex3.
在我方某前沿防守地域,敌人以一个炮排(含两门火炮)为单位对我方进行干扰和破坏.为躲避我方打击,敌方对其阵地进行了伪装并经常变换射击地点.经过长期观察发现,我方指挥所对敌方目标的指示有50%是准确的,而我方火力单位,在指示正确时,有1/3的射击效果能毁伤敌人一门火炮,有1/6的射击效果能全部毁伤敌人火炮.现在希望能用某种方式把我方将要对敌人实施的20次打击结果显现出来,确定有效射击的比率及毁伤敌方火炮的平均值。
p=0.5;m=2000; efreq=zeros(1,m);efreq2=zeros(1,m); randnum1 = binornd(1,p,1,m); randnum2 = unidrnd(6,1,m);k1=0;k2=0;k3=0; for i=1:mif randnum1(i)==0k1=k1+1;elseif randnum2(i)<=3k1=k1+1;elseif randnum2(i)==6k3=k3+1;elsek2=k2+1;endendefreq(i)=(k2+k3)/i;efreq2(i)=(k2+2*k3)/i; end num=1:m;plot(num,efreq,num,efreq2)