2009B.rar_可靠性_蒙特卡洛 matlab_蒙特卡洛 可靠_蒙特卡洛 排队_蒙特卡洛随机

clc lamdawsh=1.0484; lamdashwm=2.8167; lamdabnzhd=0.4*3.7742; lamdabnzhsh=0.6*3.7742; lamdaqgy=1.0333; bed=zeros(1,79); quenewsh=0; queneshwm=0; quenebnzhd=0; quenebnzhsh=0; queneqgy=0; lq=0; ls=0; N=50; n1=0; waittime=0; for day=1:N quenewsh=quenewsh+1 queneshwm=queneshwm+1 quenebnzhd=quenebnzhd+1 quenebnzhsh=quenebnzhsh+1 queneqgy=queneqgy+1 bed=bed-1 week=mod(day,7) renshuwsh=poissrnd(lamdawsh) if renshuwsh>0 renshuwsh=rand(1,renshuwsh); renshuwsh=sort(renshuwsh); renshuwsh=renshuwsh(end:-1:1); quenewsh=[quenewsh renshuwsh] end renshushwm=poissrnd(lamdashwm) if renshushwm>0 renshushwm=rand(1,renshushwm); renshushwm=sort(renshushwm); renshushwm=renshushwm(end:-1:1); queneshwm=[queneshwm renshushwm] end renshubnzhd=poissrnd(lamdabnzhd) if renshubnzhd>0 renshubnzhd=rand(1,renshubnzhd); renshubnzhd=sort(renshubnzhd); renshubnzhd=renshubnzhd(end:-1:1); quenebnzhd=[quenebnzhd renshubnzhd] end renshubnzhsh=poissrnd(lamdabnzhsh) if renshubnzhsh>0 renshubnzhsh=rand(1,renshubnzhsh); renshubnzhsh=sort(renshubnzhsh); renshubnzhsh=renshubnzhsh(end:-1:1); quenebnzhsh=[quenebnzhsh renshubnzhsh] end renshuqgy=poissrnd(lamdaqgy) if renshuqgy>0 renshuqgy=rand(1,renshuqgy); renshuqgy=sort(renshuqgy); renshuqgy=renshuqgy(end:-1:1); queneqgy=[queneqgy renshuqgy] end while max([length(quenewsh),length(queneshwm),length(quenebnzhd),length(quenebnzhsh),length(queneqgy)])>1&min(bed)<1 index=find(bed<1); n1=n1+1; if length(quenewsh)>1 waittime=waittime+floor(quenewsh(2)) quenewsh(2)=[]; r=rand; % the rand number if r<2/55 t=4; %time in hospital elseif r<14/55 t=5; elseif r<24/55 t=6; elseif r<33/55 t=7; elseif r<44/55 t=8; elseif r<49/55 t=9; elseif r<52/55 t=10; else t=11; end bed(index(1))=t; else if (week==0|week==6)&(length(quenebnzhd)>1|length(quenebnzhsh)>1) if length(quenebnzhsh)>1 waittime=waittime+floor(quenebnzhsh(2)) quenebnzhsh(2)=[]; r=rand; % the rand number if r<17/83 t=2; %time in hospital elseif r<70/83 t=3; else t=4; end if week==0 bed(index(1))=t+1; else bed(index(1))=t+2; end else waittime=waittime+floor(quenebnzhd(2)) quenebnzhd(2)=[]; r=rand; % the rand number if r<21/71 t=2; %time in hospital elseif r<57/71 t=3; else t=4; end if week==0 bed(index(1))=t+1; else bed(index(1))=t+2; end end elseif (week==1|week==2)&length(quenebnzhd)>1 waittime=waittime+floor(quenebnzhd(2)) quenebnzhd(2)=[]; r=rand; % the rand number if r<21/71 t=2; %time in hospital elseif r<57/71 t=3; else t=4; end if week==1 bed(index(1))=t+1; else bed(index(1))=t+2; end else if length(queneqgy)==1&length(queneshwm)==1 break; end if length(queneqgy)==1 shwm=1; else if length(queneshwm)==1 shwm=0; else if queneshwm(2)==max([queneshwm(2),queneqgy(2)]) shwm=1; else shwm=0; end end end if shwm>0 queneshwm(2) waittime=waittime+floor(queneshwm(2)) queneshwm(2)=[]; r=rand; if r<0.05 t=6; % observation time elseif r<0.12 t=7; elseif r<0.28 t=8; elseif r<0.38 t=9; elseif r<0.57 t=10; elseif r<0.72 t=11; elseif r<0.82 t=12; elseif r<0.91 t=13; else t=14; end if week==0 %sunday bed(index(1))=t+2; elseif week==1|week==6 bed(index(1))=t+3; else if r<0.62 bed(index(1))=t+2; else bed(index(1))=t+3; end end else waittime=waittime+floor(queneqgy(2)) queneqgy(2)=[]; r=rand; if r<4/38 t=6; % observation time elseif r<12/38 t=7; elseif r<25/38 t=8; elseif r<33/38 t=9; elseif r<35/38 t=10; else t=14; end if week==0 %sunday bed(index(1))=t+2; elseif week==1|week==6 bed(index(1))=t+3; else if r<0.59 bed(index(1))=t+2; else bed(index(1))=t+3; end end end end end end lq=lq+length(quene)-1 ls=ls+length(quene)-1+length(find(bed>0)) end lq=lq/N ls=ls/N waittime=waittime/n1