[matlab]Monte Carlo模拟学习笔记
理論基礎(chǔ):大數(shù)定理,當(dāng)頻數(shù)足夠多時,頻率可以逼近概率,從而依靠概率與$\pi$的關(guān)系,求出$\pi$
? ? ? ?? 所以,rand在Monte Carlo中是必不可少的,必須保證測試數(shù)據(jù)的隨機(jī)性。
?
用蒙特卡洛方法進(jìn)行計算機(jī)模擬的步驟:
[1] 設(shè)計一個邏輯框圖,即模擬模型.
[2] 根據(jù)流程圖編寫程序,模擬隨機(jī)現(xiàn)象.可通過具有各種概率分布的模擬隨機(jī)數(shù)來模擬隨機(jī)現(xiàn)象.
[3] 分析模擬結(jié)果,計算所需要結(jié)果.
?
ex1.投針試驗求$\pi$
%蒲豐投針實驗的計算機(jī)模擬 format long; %設(shè)置15位顯示精度 a=1; l=0.6; %兩平行線間的寬度和針長 figure; axis([0,pi,0,a/2]); %初始化繪圖板 set(gca,'nextplot','add'); %初始化繪圖方式為疊加 counter=0; n=1120; %初始化計數(shù)器和設(shè)定投針次數(shù) x=unifrnd(0,a/2,1,n); phi=unifrnd(0,pi,1,n); %樣本空間Ω frame=moviein(n); %建立一個1120列的大矩陣 for i=1:nif x(i)<l*sin(phi(i))/2 %滿足此條件表示針與線的相交plot(phi(i),x(i),'b.');counter=counter+1; %統(tǒng)計針與線相交的次數(shù)frame(:,counter)=getframe; %描點并取幀end end fren=counter/n; pihat=2*l/(a*fren); %用頻率近似計算π disp(counter); disp(pihat);
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.
在我方某前沿防守地域,敵人以一個炮排(含兩門火炮)為單位對我方進(jìn)行干擾和破壞.為躲避我方打擊,敵方對其陣地進(jìn)行了偽裝并經(jīng)常變換射擊地點.經(jīng)過長期觀察發(fā)現(xiàn),我方指揮所對敵方目標(biāo)的指示有50%是準(zhǔn)確的,而我方火力單位,在指示正確時,有1/3的射擊效果能毀傷敵人一門火炮,有1/6的射擊效果能全部毀傷敵人火炮.現(xiàn)在希望能用某種方式把我方將要對敵人實施的20次打擊結(jié)果顯現(xiàn)出來,確定有效射擊的比率及毀傷敵方火炮的平均值。
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)
?
轉(zhuǎn)載于:https://www.cnblogs.com/elpsycongroo/p/7215619.html
總結(jié)
以上是生活随笔為你收集整理的[matlab]Monte Carlo模拟学习笔记的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 收红莲表情,有的留名字。我有资格哦
- 下一篇: .NET调用JAVA的WebServic