計算機模擬方法

1、 計算機模擬方法一四人追逐實驗模擬如圖1,在正方形ABCD的四個頂點各有一個人。設(shè)在初始時刻時，四人同時出發(fā)勻速以沿順時針走向下一個人。如果他們始終對準下一個人為目標行進，最終結(jié)果會如何。作出各自的運動軌跡。 解：該問題可以通過計算機模擬來實現(xiàn)。這需要將時間離散化。設(shè)時間間隔為，時刻表示時間 設(shè)第個人時刻的位置坐標為： 對前面3個人表達式為： 其中 對第4個人表達式為： 其中 Matlab實現(xiàn)程序run.m如下：%模擬運動n=240;x=zeros(4,n);y=zeros(4,n);dt=0.05; %時間間隔v=10; %速度x(1,1)=100; y(1,1)=0; %第1個人初始坐標x

2、(2,1)=0; y(2,1)=0; %第2個人初始坐標x(3,1)=0; y(3,1)=100; %第3個人初始坐標x(4,1)=100; y(4,1)=100; %第4個人初始坐標for j=1:n-1 for i=1:3 d=sqrt(x(i+1,j)-x(i,j)2+(y(i+1,j)-y(i,j)2); %第i個人和第i+1個人距離 cosx=(x(i+1,j)-x(i,j)/d; %求cos值 sinx=(y(i+1,j)-y(i,j)/d; %求sin值 x(i,j+1)=x(i,j)+v*dt*cosx; %求新x坐標 y(i,j+1)=y(i,j)+v*dt*sinx; %求

3、新y坐標 end %考慮第1,2,3人運動一步 d=sqrt(x(1,j)-x(4,j)2+(y(1,j)-y(4,j)2); %第4個人和第1個人距離 cosx=(x(1,j)-x(4,j)/d; %求cos值 sinx=(y(1,j)-y(4,j)/d; %求sin值 x(4,j+1)=x(4,j)+v*dt*cosx; %求第4點新x坐標 y(4,j+1)=y(4,j)+v*dt*sinx; %求第4點新y坐標 end%plot(x,y)for j=1:nplot(x(1,j),y(1,j),x(2,j),y(2,j),x(3,j),y(3,j),x(4,j),y(4,j) %作點圖ho

4、ld on %保持每次作圖,實現(xiàn)各次圖行迭加end執(zhí)行結(jié)果見圖1 圖1 模擬結(jié)果圖形二、電梯問題隨機模擬 設(shè)有個人在一樓進入電梯，樓上共有層。設(shè)每個乘客在任何一層樓出電梯的可能性相同，求直到電梯中的人下完為止，電梯需停次數(shù)的數(shù)學期望。并對進行計算機模擬驗證。 計算機模擬算法思想： 樓上層的序號記為，定義整型數(shù)組,并將其個元素的初始值置0。計算機隨機產(chǎn)生1到的個隨機整數(shù)，若該數(shù)為，則令，表示第層電梯有人要下，從而電梯在該層要停；當有多人在第層要下時，仍用表示。統(tǒng)計數(shù)組中為1的元素個數(shù)，即為該次模擬實驗電梯要停的次數(shù)?？偣材M實驗次，求平均值就是電梯所停次數(shù)的模擬值。 Matlab模擬程序dian

5、ti.mN=5000; %模擬次數(shù)n=28; %電梯層數(shù)r=14; %電梯開始進的人數(shù)s=0;x=zeros(n,1); for k=1:N %模擬N次 s1=0; for i=1:n x(i)=0; end for j=1:r %對每個人是否下電梯進行模擬 i=1+floor(rand(1,1)*n); x(i)=1; %第i層有人下 end s1=sum(x); %該次模擬中總共要下的人數(shù) s=s+s1; %統(tǒng)計各次模擬中下的人數(shù) end eq=s/N; %模擬平均值 ei=n*(1-(1-1/n)r); fprintf('模擬電梯要停次數(shù)的理論值%6.2f,模擬值為%6.2fn&

6、#39;,ei,eq);某兩次運行結(jié)果如下：11.15,11.18。與理論數(shù)值11.17很接近。三、 礦工選門問題模擬實驗 一名礦工陷入一個有三扇門的礦井中，第一扇門通到一個隧道，走2小時可以到達安全區(qū)；第二扇門通到另一個隧道，走3小時回到礦井中；第三扇門通到又一個隧道，走5小時回到礦井中。假定礦工總是等可能的選擇三扇門的任何一扇門，求到達安全區(qū)的平均時間計算機模擬的Matlab程序： N=50000; sum=0;for index=1:N flag=1; Time=0;while flag=1 %一次試驗直到走出為止door=1+floor(3*rand(1,1); %隨機選擇1個門if(

7、door=1) Time=Time+2; flag=0;%選擇1號門elseif(door=2) Time=Time+3; %選擇2號門elseif(door=3) Time=Time+5; %選擇3號門 end end sum=sum+Time; end %end for index aver=1.0*sum/N; fprintf('平均時間為%5.2f小時.n',aver);某三次模擬結(jié)果為：10.04，10.14，9.91。與理論值10小時很接近。 基于時間步長的計算機仿真，就是按照時間流逝的順序，一步一步地對系統(tǒng)進行仿真，在整個仿真過程中，通常將時間步長固定不變。 在進

8、行系統(tǒng)仿真過程中，可以把整個仿真過程分為許多相等的時間間隔，時間步長的長度可根據(jù)實際問題分別取為秒、分、小時等。程序中按此步長前進的時鐘就是仿真的時鐘。四 容器水流問題仿真實驗 某容器內(nèi)盛入鹽水池有3000m3水，其中含鹽3kg，以每分鐘10 m3的速率向水池中注入含鹽率為0.5kg/ m3的鹽水。同時以每分鐘7 m3的速率從水池流出攪拌均勻的鹽水。求1小時后水池中的含鹽量。并做出水池含鹽量的變化曲線。分析：1 理論求解設(shè)初始時刻為0，t時刻水池含鹽量為kg。t時刻水池的水量為 m3經(jīng)過時間流入水池的鹽量為 t時刻水池的含鹽率為：經(jīng)過時間流出水池的鹽量為 則水池含鹽量滿足關(guān)系： 計算機仿真設(shè)仿

9、真的時間步長為分鐘，時刻時間為，水池含鹽量為kg由前面分析可得遞推關(guān)系為： 終止時間為1小時=60分鐘，對應步數(shù) Matlab程序為：dt=0.01; %分鐘n=60/dt; %終止步數(shù)t=zeros(1,n+1);x=zeros(1,n+1);x(1)=3;t(1)=0;for i=1:nt(i+1)=t(i)+dt; x(i+1)=x(i)+5*dt-x(i)*7*dt/(3000+3*i*dt); endfprintf('時間%5.1f分鐘,水池含鹽量為%6.3f kgn',t(n+1),x(n+1);plot(t,x);grid on;輸出結(jié)果為：時間 60.0分鐘,水

10、池含鹽量為283.306 kg.仿真值和理論值十分接近。圖形為： 圖2.16 水池含鹽量曲線圖五 庫存問題仿真實驗 在物資供應過程中，由于到貨與銷售不可能做到同步同量，因此總要保持一定的庫存量。如果庫存過多，就會造成積壓浪費和保管費的增大；如果庫存過少，會造成缺貨。如何選擇庫存和缺貨策略，是一個值得研究的問題。庫存問題有多種類型，通常比較復雜，這里我們通過計算機仿真求解一種較為簡單的情況。某自行車商店的倉庫管理人員采取一種簡單定貨策略，當庫存量降低到P輛自行車就向廠家訂貨，每次訂貨Q量，如果某一天的需求量超過了庫存量，商店就有缺貨損失和信譽損失；但如果庫存量過多，將會導致資金積壓和保管費增加。

11、若現(xiàn)在有六種策略見表2.12,試選擇一種總費用最少的策略。 表2.12 訂貨策略方案編號訂貨點P量訂貨量Q輛112515021252503150250417525051753006200300已知條件為：（1） 從發(fā)出貨到收到貨需要隔3天。如周一訂貨則周四才能到貨。（2） 每輛自行車每天保管費為0.75元，每輛自行車每天缺貨損失為1.80元，每次訂貨費為75元。（3） 每天自行車的市場需求量服從參數(shù)為50輛的Poisson分布。（4） 原始庫存為120輛，并假設(shè)第一天沒有到貨也沒有發(fā)出訂貨。 該問題采用理論求解比較麻煩，我們采用計算機仿真，計算每天的庫存量，銷售量，計算一段時間如n=200天的

12、總費用及平均費用。然后選取總費用最少的策略。 計算機仿真分析： 設(shè)第天早上（還沒有營業(yè)前）商店擁有量為輛，白天銷售后到晚上（停止營業(yè)后）擁有量為輛(可為負，負表示該天缺貨量)，市場需求量為輛，每天到貨量為輛（沒有到貨則為0）。第天缺貨量為輛，庫存量為輛。保管費用元，缺貨損失費元，訂貨費。 則有如下關(guān)系： 1）初始值第一天早上商店擁有量，到貨量. 2）隨機產(chǎn)生每天服從Poisson分布的市場需求量, 3）白天銷售后到晚上（停止營業(yè)后）擁有量為滿足： 4）缺貨量計算式為： 5）庫存量計算式為： 6）到貨量及訂貨費計算： 初始置所有， 計算當庫存量，則,7）第天早上商店擁有量為滿足： 8) 第天保管

13、費計算式為： 9) 第天缺貨損失費計算式為： 10)總費用計算11)平均每天費用 Matlab程序：n=200; %總天數(shù)x=zeros(1,n+3);%早上擁有量y=zeros(1,n+3);%晚上擁有量a=zeros(1,n+3); %到貨量d=zeros(1,n+3); %市場需求量u=zeros(1,n+3); %缺貨量s=zeros(1,n+3); %庫存量c=zeros(1,n+3); %保管費r=zeros(1,n+3); %缺貨損失費v=zeros(1,n+3); %訂貨費P=120; %訂貨點Q=150; %訂貨量x(1)=120; %初始早上擁有量a(1)=0; %初始到貨

14、量d=poissrnd(50,1,n); %產(chǎn)生每天的市場需求for i=1:n y(i)=x(i)-d(i)+a(i); %第i天晚上擁有量 if(y(i)<0) u(i)=-y(i); end %當天缺貨量 if(y(i)>0) s(i)=y(i); end %當天庫存量 if(s(i)<P) a(i+3)=Q; %后3天將有訂貨到達 v(i)=75 ; %當天訂貨費 end x(i+1)=s(i); %第二天早上擁有量為前一晚庫存量 c(i)=0.75*s(i); %當天庫存費 r(i)=1.8*u(i) ;%當天缺貨損失費 endCost=sum(c(1:n)+sum(r(1:n)+sum(v(1:n); %總費用Aver=Cost/n; %平均費用fprintf('當P=%2d,Q=%2d,n=%3d天,總費用%6.2f,平均費用%6.2fn',P,Q,n,Cost,Aver);某次仿真實驗結(jié)果為：當P=125,Q=150,n=200天,總費用35630.10,平均費用178.15當P=125,Q=250,n=200天,總費用48503.70,平均費用242.52當P