1、系统仿真与 matlab综合试题 .0M/M/N 排队系统的模拟仿真 1摘 要 .11. 问题分析 .22. 模型假设 .23. 符号说明 .34. 模型准备 .34.1 排队系统的组成和特征 34.1.1 输入过程 .44.1.2 排队规则 .44.1.3 服务过程 .44.1.4 排队系统的主要指标 .54.2 输入过程与服务时间的分布 .54.2.1 负指数分布 .54.2.2 泊松分布 .54.3 生灭过程 .65. 标准 M/M/N 模型 .85.1 多服务台模型准备 .85.2 多服务台模型建立 .95.2.1 服务利用率 .95.2.2 平均排队长 .95.2.3 平均队长 10
2、5.2.4 平均等待时间 106. 程序设计 116.1 动画流程图 116.2 M/M/N 流程图 127. 程序运行实例介绍 137.1 动画实例讲解 137.2M/M/N 排队系统实例讲解 .148. 程序实现难点和模型评价 178.1 程序实现难点 178.2 模型评价 179. 参考文献 1710. 附录 .1710.1 动画实现的核心程序 .1710.2 M/M/N 模型计算主要程序 .221M/M/N 排队系统的模拟仿真摘 要排队是在日常生活中经常遇到的事，由于顾客到达和服务时间的随机性，使得排队不可避免。因此，本文建立标准的 M/M/N 模型，并运用 Matlab 软件，对 M
3、/M/N 排队系统就行了仿真，从而更好地深入研究排队问题。问题一，基于顾客到达时间服从泊松分布和服务时间服从负指数分布，建立了标准的 M/M/N 模型。运用 Matlab 软件编程，通过输入服务台数量、泊松分布参数以及负指数分布参数，求解出平均队长、服务利用率、平均等待时间以及平均排队长等重要指标。然后，分析了输入参数与输出结果之间的关系。得出当服务台数增加时，几个参数都会变小的结论。问题二，为了更加清晰地反映出实际排队过程。本文通过运用 Matlab 软件编程，制作了 M/M/1 排队过程的动画仿真，通过输入泊松分布参数以及负指数分布参数来模拟不同情况下的排队过程。通过仿真动画，可以看到明显
4、的等待和排队过程。问题三，为了清晰地展示程序执行的效果以及程序功能的使用方法。本文特意制作了程序运行指南，并做了程序运行实例分析。通过详细地介绍，使读者能更好地理解 M/M/N 模型以及如何使用该仿真程序。最后，对建立的 M/M/N 模型做了评价，并提出了一些改进的思路。同时，指出了程序实现的难点等问题。关键词： M/M/N 排队系统 泊松分布 负指数分布 动画模拟仿真21. 问题分析排队论 （Queuing Theory）也称随机服务系统理论，就是为解决有关排队问题而发展的一门学科。它研究的内容有下列三部分： 1性态问题，即研究各种排队系统的概率规律性，主要是研究队长分布、等待时间分布和忙期
5、分布等，包括了瞬态和稳态两种情形。 2最优化问题，又分静态最优和动态最优，前者指最优设计。后者指现有排队系统的最优运营。3排队系统的统计推断，即判断一个给定的排队系统符合于哪种模型，以便根据排队理论进行分析研究。 其过程如下图：本文需要解决的问题：1建立顾客到达时间服从泊松分布、服务时间服从负指数分布的 M/M/N 排队模型，并利用 Matlab 软件实现输入参数的键入以及输出参数的显示。2运用 Matlab 软件编程制作 M/M/1 排队系统的动态仿真模拟动画，并拥有输入参数的键入功能。3制作程序运行指南，并结合程序运行实例对程序功能作深入分析。4对本文建立的标准 M/M/N 排队模型作评价
6、。2. 模型假设针对本问题，建立如下合理的假设：1顾客源是无穷的；2排队长度没有限制；33到达系统的顾客按先到先服务原则依次进入服务；4服务员在仿真过程中没有休假；5顾客到达时排成一队，当有服务台空闲时进入服务状态；6单位时间内到达的顾客数量服从泊松分布；7顾客所需的服务时间服从负指数分布；8各服务台工作是相互独立且平均服务时间相同。3. 符号说明符号 说明 单位顾客到达时间参数 人数/分顾客服务时间参数 人数/分P出现某种状态的概率 ps服务利用率 L平均排队长 人s平均队长 人Ws平均逗留时间 分钟q平均等待时间 分钟4. 模型准备4.1 排队系统的组成和特征 一般的排队过程都由输入过程、
7、排队规则、服务过程三部分组成，现分述如下： 44.1.1 输入过程 输入过程是指顾客到来时间的规律性，可能有下列不同情况： 1. 顾客的组成可能是有限的，也可能是无限的。2. 顾客到达的方式可能是一个个的，也可能是成批的。3. 顾客到达可以是相互独立的，即以前的到达情况对以后的到达没有影响，否则是相关的。4. 输入过程可以是平稳的，即相继到达的间隔时间分布及其数学期望、方差等数字特征都与时间无关，否则是非平稳的。 4.1.2 排队规则 排队规则指到达排队系统的顾客按怎样的规则排队等待，可分为损失制，等待制和 混合制三种。 1. 损失制（消失制） 。当顾客到达时，所有的服务台均被占用，顾客随即离
8、去。 2. 等待制。当顾客到达时，所有的服务台均被占用，顾客就排队等待，直接受完服务才离去。例如出故障的机器排队等待维修就是这种情况。 3. 混合制。介于损失制和等待制之间的是混合制，即既有等待又有损失。有队列长度有限和排队等待时间有限两种情况，在限度以内就排队等待，超过一定限度就离去。 排队方式还分为单列、多列和循环队列。 4.1.3 服务过程 1. 服务机构。主要有以下几种类型：单服务台；多服务台并联（每个服务台同时为不同顾客服务） ；多服务台串联（多服务台依次为同一顾客服务） ；混合型。 2. 服务规则。按为顾客服务的次序采用以下几种规则： 1) 先到先服务，这是通常的情形。2) 后到先
9、服务，如情报系统中，最后到的情报信息往往最有价值，因而常被优先处理。3) 随机服务，服务台从等待的顾客中随机地取其一进行服务，而不管到达的先后。4) 优先服务，如医疗系统对病情严重的病人给予优先治疗。 54.1.4 排队系统的主要指标1. 平均队长：指系统内顾客数（包括正被服务的顾客与排队等待服务的顾客）的数学期望。2. 平均排队长:指系统内等待服务的顾客数的数学期望。3. 平均逗留时间：顾客在系统内逗留时间（包括排队等待的时间和接受服务的时间）的数学期望。4. 平均等待时间：指一个顾客在排队系统中排队等待时间的数学期望。5. 平均忙期：指服务机构连续繁忙时间（顾客到达空闲服务机构起，到服务机
10、构再次空闲止的时间）长度的数学期望。4.2 输入过程与服务时间的分布排队系统中的事件流包括顾客到达流和服务时间流。由于顾客到达的间隔时间和服务时间不可能是负值，因此，它的分布是非负随机变量的分布。最常用的分布有泊松分布、确定型分布，指数分布和爱尔朗分布。由于本文只用到了泊松分布和负指数分布，因此只对这两种分布加以说明。4.2.1 负指数分布指数分布是单参数 的非对称分布，记作 ，概率密度函数为：)(Exp10,ttfet它的数学期望为 ，方差为 。121指数分布是唯一具有无记忆性的连续型随机变量，即有，在排队论、可靠性分析中有广泛应用。本文将用sXPtstXP|负指数分布来产生顾客的服务时间。
11、4.2.2 泊松分布泊松分布与指数分布有密切的关系。当顾客平均到达率为常数 的到达间隔服从指数分布时，单位时间内到达的顾客数 K 服从泊松分布，即单位时间内到达 k 位顾客的概率为62!keP记作 Poisson() 。泊松分布在排队服务、产品检验、生物与医学统计、天文、物理等领域都有广泛应用。 本文将用泊松分布来产生单位时间内到达的顾客数目。 4.3 生灭过程在排队论中，如果 表示时刻 系统中的顾客数，则 就构成了一个tNt 0,tN随机过程。如果用“生”表示顾客的到达， “灭”表示顾客的离去，则对许多排队过程来说， 就是一类特殊的随机过程生灭过程。0,t定义 1 设 为一个随机过程。若 的
12、概率分布具有以下性质： t1. 假设 ，则从时刻 起到下一个顾客到达时刻止的时间服从参数为ntNt的负指数分布， 。 n 2,102. 假设 ，则从时刻 起到下一个顾客离去时刻止的时间服从参数为tt的负指数分别， 。n ,n3. 同一时刻只有一个顾客到达或离去。则称 为一个生灭过程。0,tN当系统运行相当时间而到达平衡状态后，对任一状态 ，即n其中 表示系统中一共有 名顾客，单位时间内2,1ntPtpn n进入该状态的平均次数和单位时间内离开该状态的平均次数应该相等，这就是系统在统计平衡下的“流入流出”原理。根据这一原理，可得到任一状态下的平衡方程如下： ppnnnn 112231 11200
13、121：7有上述平衡方程，可求得 pppnn011012301201210：因此，记 32,11021 nnnC则平稳状态的分布为 42,0npn由概率分布的要求 510n可以得到 610nCp即系统空闲状态的概率。注意只有当级数 收敛时才有意义。1n85. 标准 M/M/N 模型5.1 多服务台模型准备设顾客单个到达，相继到达时间间隔服从参数为 的负指数分布，系统中共有个服务台，每个服务台的服务时间相互独立，且服从参数为 的负指数分布。s 当顾客到达时，若有空闲的服务台则马上接受服务，否则便排成一个队列等待，等待时间为无限。记 为系统达到平稳状态后的队长 N 的概率分布，注2,10nNPpn
14、意到对个数为 s 的多服务台系统，有72,10,n和 8,1,2snsn记服务强度 ，s 9ss，则当 时，由式（3） 、 （4） 、 （5） 、 （6） ，可以得到1s10,!,21,!snsssnCn 故91,!,21,00snspnn其中 1210!1ssn为系统空闲的概率。因此，系统中有任意多个顾客的概率都可以由式（11）和式（12）得到。从而，就可以计算出反映该系统性能的各种重要指标。5.2 多服务台模型建立5.2.1 服务利用率在本文中，简单的把服务利用率定义为系统处于非空闲状态的概率。因此，由公式（12） ，可以得到服务利用率为 1310!1ssnPS5.2.2 平均排队长对于多
15、服务台排队系统，由公式（11）和公式（12）可以得到系统到达平稳状态的平均排队长 。Lp14!1201001 sdssnssnnss snssnp 105.2.3 平均队长对于到达平稳状态的多服务台排队系统，平均队长 等于平均排队长与正Ls在接受服务的顾客的平均数之和，即 平 均 数正 在 接 受 服 务 的 顾 客 的平 均 排 队 长平 均 队 长 Lps记系统中正在接受服务的顾客的平均数为 ，则s 151!10 0100010sn sssssnnsnp因此，平均队长 16Lps5.2.4 平均等待时间对于顾客在系统中平均等待时间的建模，可以先计算顾客在系统中的逗留时间，然后再减去顾客在系
16、统中的服务时间，这样就可以得到顾客的平均等待时间。顾客在系统中的逗留时间可以看作一个服从 的负指数分布，即-17ettTP因此，平均逗留时间 81Ws又因为顾客在系统中的逗留时间 等于等待时间 与接受服务的时间 之和，TqV即 19Vq故平均等待时间 201sq116. 程序设计6.1 动画流程图 开始输入顾客到达的泊松分布参数、服务时间的负指数分布参数判断时间是否等于第i 个人的到达时间？NOYESNO进入系统接受服务YES离开是否空闲？服务完毕？ NOYES结束？ENDYESNO126.2M/M/N 流程图输入参数判断输入是否正确？NO代入公式计算,并显示结果YES结束？ENDYESNO开
17、始137. 程序运行实例介绍7.1 动画实例讲解当运行程序时就会出现下面的窗口：如果要观看 M/M/1 排队系统的仿真动画，可以现在左下角输入每分钟到达人数和每分钟服务人数。然后，点击“观看动画”按钮，就可以观看仿真动画。例如，设每分钟到达人数为 0.35、每分钟服务人数为 0.4.程序运行一段时间后，我们可以看到如下的画面：14由上图，可以得到等待人数、总接待人数以及离开人数。如果要重新设置参数，可以按左下角的“重置”键。当“重置”键按下后，我们可以看到如下画面：然后，重新设置每分钟到达人数以及每分钟服务人数，并点击“观看动画” 。7.2M/M/N 排队系统实例讲解在窗口的“请设置参数”下方
18、，可以设置服务台数、每分钟到达人数以及每分钟服务人数。然后，点击“开始”按钮，在“模型仿真结果”下方就可以得到平均队长、服务利用率、平均等待时间以及平均排队长度。15如果要重新设置参数，可以点击开始键旁边“重置” ，点击后可得到如下画面：可以看见“M/M/N 排队系统计算 ”下面所有的框都被置为 0。如果输入参数出错时，程序会弹出窗口自动提示，如下图：16并且在切换到提示时，并不会影响前面动画的执行。因此，不必担心输入参数出错而影响其他程序的运行。点击“返回”即可回到前一画面，并且可以看见“M/M/N 排队系统计算”所有的框都被置为 0。点击后的画面如下图：注：本程序在执行时是以一秒钟代替一分
19、钟178. 程序实现难点和模型评价8.1 程序实现难点1. 本段程序的主要难点在制作动画中，在制作动画时不是很容易把泊松分布和负指数分布产生的数据和实际的时间联系起来。最后，利用了 clock 函数解决了问题。2. 在 GUI 设计时，切换画面时总是会有一定的闪动，此问题暂时还没有解决。8.2 模型评价本文建立的是标准的 M/M/N 排队系统模型，而在实际生活中并不是这样的。在实际的生活中往往是一种有损失制的排队系统，当人们在排队的时候看见排队的人数很多就会部分人就会离开，而从模型并没有考虑此情况。因此，模型还有待改进。9. 参考文献1运筹学 教材编写组, 运筹学（第三版） ，北京：清华大学出
20、版社M,20052齐欢 王小平系，统建模与仿真 北京：清华大学出版社M，2004。3姜启源 谢金星，数学模型(第三版) ，北京：高等教育出版社M，2003。4赵广元,MATLAB 与控制系统仿真实践，北京：北京航空航天大学出版社M,20095施晓红 周佳，精通 GUI 图形界面编程，北京：北京大学出版社M，200310. 附录10.1 动画实现的核心程序function Mypush_Callback(hObject, eventdata, handles)% hObject handle to Mypush (see GCBO)% eventdata reserved - to be def
21、ined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)axes(handles.Myaxe);18N=10000;% geta=str2double(get(findobj(tag,MYGETA),String);% serveb=str2double(get(findobj(tag,MYSERVEB),String);geta=str2double(get(handles.MYGETA,String);serveb=str2double(get(handles
22、.MYSERVEB,String);geta=1/geta;serveb=1/serveb;reset=0;%动画标志nownum=0;%初始人数servenum=0;%接受服务的总人数outnum=0;% state=0;%0服务台无人，1有人gettime0=ceil(poissrnd(geta,1,N);%产生到达时间gettime=cumsum(gettime0(1,:);servetime0=ceil(exprnd(serveb,1,N);%产生服务时间servetime=cumsum(servetime0(1,:);leavetime(1)=gettime(1)+servetime
23、0(1);for j=2:Nleavetime(j)=leavetime(j-1)+servetime0(j);%等于前一个离开时间加其服务时间endwaittime(1)=0;for j=2:Nwaittime(j)=-gettime(j)+leavetime(j-1);%等待时间endx=10;y=0;h=plot(x,y,.);hold onx0=10;y0=8;h1=plot(x0,y0,.);hold onxl=10;yl=16;h2=plot(xl,yl,.);hold on%服务台19x1=8 12 12 8;%指定x坐标的值y1=15 15 17 17;%指定y坐标的值X1=x
24、1 x1(1);%首尾相连Y1=y1 y1(1);%首尾相连plot(X1,Y1);fill(X1,Y1,y);text(8,18,服务台);hold on%等待处x2=8 12 12 8;%指定x坐标的值y2=7 7 9 9;%指定y坐标的值X2=x2 x2(1);%首尾相连Y2=y2 y2(1);%首尾相连plot(X2,Y2);fill(X2,Y2,r);text(8,10,等待处);hold on%等待人数x3=14 16 16 14;%指定x坐标的值y3=7 7 9 9;%指定y坐标的值X3=x3 x3(1);%首尾相连Y3=y3 y3(1);%首尾相连plot(X3,Y3);fil
25、l(X3,Y3,g);text(13,10,等待人数);hold on%总接待人数x4=3 5 5 3;%指定x坐标的值y4=15 15 17 17;%指定y坐标的值X4=x4 x4(1);%首尾相连Y4=y4 y4(1);%首尾相连plot(X4,Y4);fill(X4,Y4,b);text(1,18,总接待人数);hold on%入口x5=9 11 11 9;%指定x坐标的值y5=0 0 2 2;%指定y坐标的值X5=x5 x5(1);%首尾相连Y5=y5 y5(1);%首尾相连plot(X5,Y5);fill(X5,Y5,g);text(11.5,1,入口);20hold on%出口x6
26、=17 19 17 19;%指定x坐标的值y6=15 15 17 17;%指定y坐标的值X6=x6 x6(1);%首尾相连Y6=y6 y6(1);%首尾相连plot(X6,Y6);fill(X6,Y6,g);text(18,18,出口);hold on%离开人数x7=14 16 16 14;%指定x坐标的值y7=15 15 17 17;%指定y坐标的值X7=x7 x7(1);%首尾相连Y7=y7 y7(1);%首尾相连plot(X7,Y7);fill(X7,Y7,g);text(12,18,离开人数);hold onaxis(0 20 0 20)axis squaregrid offset(h
27、,EraseMode,xor,MarkerSize,18)%异或实现动画set(h1,EraseMode,xor,MarkerSize,18)set(h2,EraseMode,xor,MarkerSize,18)i=1;j=1;t=1;temp=clock;nowtime=temp(6)+temp(5)*60+temp(4)*60*60+temp(3)*24*60*60;nowtime0=nowtime;%记录开始仿真的时间while 1 geta=str2double(get(handles.MYGETA,String);serveb=str2double(get(handles.MYSER
28、VEB,String);if(geta=0|serveb=0)%按下重置键h_axes=findobj(gcf,type,axes); %获得当前图中所有坐标的句柄h_children_axes=allchild(h_axes); %获得坐标的子对象的句柄delete(h_children_axes); %删除所有坐标子对象句柄，达到清空目的break;end% if(geta=-1|serveb=-1)%按下退出键% break;21% endif(reset=1)%把坐标置到（10,0）入口处x=10;y=0;reset=0;endtime=clock;curtime=time(6)+ti
29、me(5)*60+time(4)*60*60+time(3)*24*60*60;%顾客到来if(i=gettime(i)y=y+0.005;set(h,XData,x,YData,y)if(y=8)reset=1;fill(X3,Y3,g)nownum=nownum+1;text(15,8,num2str(nownum);i=i+1;% tempnum=1;% nowtime=curtime;%更新时间% plot(x,y,o);% hold onendendend%顾客去服务if(j0)if(curtime-nowtime0=gettime(j)+waittime(j)x0=10.000;y
30、0=y0+0.005;set(h1,XData,x0,YData,y0)if(y0=16)nownum=nownum-1;servenum=servenum+1;% state=1;fill(X3,Y3,g)text(15,8,num2str(nownum);fill(X4,Y4,b);text(3.7,16,num2str(servenum);x0=10;y0=8;j=j+1;22endendendend%顾客离开if(t0)% if(state0)if(curtime-nowtime0=leavetime(t)xl=xl+0.005;yl=16;set(h2,XData,xl,YData,
31、yl)% if(xl=20)% t=t+1;% xl=10;% endif(xl=18)t=t+1;xl=10;outnum=outnum+1;fill(X7,Y7,g);text(14.5,16,num2str(outnum);% state=0;% endendendendendif(isize(gettime)if(jsize(waittime)if(tsize(leavetime)break;endendenddrawnowend10.2M/M/N 模型计算主要程序function MYstart_Callback(hObject, eventdata, handles)% hObje
32、ct handle to MYstart (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB23% handles structure with handles and user data (see GUIDATA)serve_num=str2double(get(handles.MY_taishu,String);%服务台数get_mmn=str2double(get(handles.MY_get,String);%到达时间参数serve_mmn=str2double(get(handles
33、.MY_serve,String);%服务时间参数p=get_mmn/serve_mmn;%服务强度ps=p/serve_num;for i=1:1if(ps=1|serve_num0)for j=1:ip1=p1*p;endendmmn=mmn+p1/factorial(i);endp2=1;for i=1:serve_nump2=p2*p;%求p的serve_num次方endmmn1=p2/(factorial(serve_num)*(1-ps);p0=1/(mmn+mmn1);%求平均排队长LP、平均队长LS、平均等待时间WS、服务利用率PSELP=p0*p2*ps/(factorial
34、(serve_num)*(1-ps)*(1-ps);LS=LP+p;24WS=LS/get_mmn;PSE=1-p0;%保留3位小数LP=round(LP*1000)/1000;LS=round(LS*1000)/1000;WS=round(WS*1000)/1000;PSE=round(PSE*1000)/1000;%显示set(handles.MY_duichang,string,num2str(LS);set(handles.MYpaiduichang,string,num2str(LP);set(handles.MY_waittime,string,num2str(WS);set(handles.MY_servetime,string,num2str(PSE);end25
