1、1,-基本原理、 MATLAB实现及案例分析,排队论模型与蒙特卡罗仿真,2,讲 座 提 纲,一 引例 二 排队现象 三 排队论的研究方法 四 蒙特卡罗仿真原理 五 仿真例子与分析 六 作业,3,一 引例,1 到银行取钱，发现前面有几十个人在排着队，你掉头就走：不能忍受啊！怎么不多开几家银行、再增加几个服务窗口啊！ 假如你是相关人员，你觉得应根据什么来决定是否需要开设新的银行或增加新的服务窗口要知道这次让你心烦具有随机性（偶然性）啊。2 银行一般都有几个服务窗口，过去是顾客每个窗口分别排队等待服务，而现在几乎都改为叫号制，这相当于多个窗口只排一队的服务规则。银行为什么要这么做? 有什么好处？,4
2、,排队是我们日常生活中常见的现象，如：上、下班搭乘公共汽车；顾客到商店购买物品；病员到医院看病；学生去食堂就餐等出现的排队和等待 服务现象。,二 排队现象,5,排队可以是有形的，也可以是无形的。如几个顾客打电话到出租车站要车，如果出租车站无足够车辆，则部分顾客只得在要车处等待，他们分散在不同地方，形成一个无形的排队序列。,排队论就是研究排队现象及其规律的一门学科，是运筹学的一个分支。如同数学的特质那样，排队论研究的内容比我们感觉中的排队现象要广泛得多，它是研究那些本质上都有排队特征的一类现象。具体表现在：,6,排队的不一定是人，也可以是物。例如：生产线上等待加工的原料、半成品；因故障停止运转等
3、待修理的机器等。,7,上述问题虽互不相同，但却都有要求得到某种服务的人或物以及提供服务的人或机构。排队论里把要求服务的对象统称为“顾客”,提供服务的人或机构称为“服务台”或“服务员”。,8,不同的顾客与服务组成了各式各样的服务系统。顾客为了得到某种服务而到达系统、若不能立即获得服务而又允许排队等待，则加入等待队伍，待获得服务后离开系统，见图2-1至图2-3。,图2-1 单服务台排队系统,9,图2-2 单队列S个服务台并联的排队系统,图2-3 S个队列S个服务台的并联排队系统,10,面对拥挤现象，人们总希望尽量设法减少排队，通常的做法是增加服务设施。但是增加设施的数量越多，人力、物力的支出就越大
4、，同时会出现空闲浪费。如果服务设施太少，顾客排队等待的时间就会很长，这样对顾客会带来不良影响。,11,顾客排队时间的长短与服务设施规模的大小，就构成了随机服务系统中的一对矛盾。如何做到既保证一定的服务质量指标，又使服务设施费用经济合理，恰当地解决顾客排队时间与服务设施费用大小这对矛盾。这就是随机服务系统理论排队论所要研究解决的问题。,12,3.1 排队系统的组成与特征排队系统一般有三个基本组成部分：1.输入过程；2.排队规则；3.服务机构。,三 排队论的研究方法,13,输入即顾客的到达，可有下列情况：1）顾客源可能是有限的，也可能是无限的。2）顾客是成批到达或是单个到达。3）顾客到达间隔时间可
5、能是随机的或确定的。4）顾客到达可能是相互独立或关联的。所谓独立就是以前顾客的到达对以后顾客的到达无影响。5）输入过程可以是平稳的（stationary），也可以是非平稳的。输入过程平稳的指顾客相继到达的间隔时间分布和参数（均值、方差）与时间无关；非平稳的则与时间相关。,3 . 1.1 输入过程,14,分为损失制、等待制、混合制三大类。(1)损失制 指如果顾客到达排队系统时，所有服务台都已被先来的顾客占用，那么他们就自动离开系统永不再来。典型例子是，如电话拔号后出现忙音，顾客不愿等待而自动挂断电话，如要再打，就需重新拔号。,3 . 1 . 2. 排队规则,15,(2)等待制 当顾客来到系统时，
6、所有服务台都不空，顾客加入排队行列等待服务。例如，排队等待售票，故障设备等待维修等。等待制中，服务台在选择顾客进行服务时，常有如下四种规则：先到先服务（FCFS ） 按顾客到达的先后顺序对顾客进行服务，这是最普遍的情形。此外还有后到先服务（LCFS），随机服务（RAND）和优先权服务（PR）三种情形。,16,(3)混合制这是等待制与损失制相结合的一种服务规则，一般是指允许排队，但又不允许队列无限长下去。具体说来，大致有三种： 队长有限。当排队等待服务顾客人数超过规定数量时，后来顾客就自动离去，另求他处服务。如水库的库容、旅馆的床位等都是有限的。另两种情况指等待时间和逗留时间限制的情形，略去。一
7、般的，损失制和等待制可认为是混合制的两种极端特殊情形。,17,3.1.3 服务机构,1）服务机构可以是单服务员和多服务员服务，这种服务形式与队列规则联合后形成了多种不同队列，不同形式的排队服务机构。如前图2-1到2-3：,2)服务方式分为单个顾客服务和成批顾客服务。 3)服务时间分为确定型和随机型。 4)服务时间的分布在这里我们假定是平稳的。,18,上述特征中最主要的、影响最大的是： 顾客相继到达的间隔时间分布 服务时间的分布 服务台数 D.G.Kendall在1953提出了分类法，称为Kendall记号(适用于并列服务台)即：X/Y/Z， 式中：X顾客相继到达间隔时间分布。 M负指数分布Ma
8、rkov， D确定型分布Deterministic， EkK阶爱尔朗分布Erlang，,3.2 排队系统的描述符号与模型分类,19,GI 一般相互独立随机分布(General Independent)， G 一般随机分布， Y填写服务时间分布（与上同）， Z填写并列的服务台数。 如 M/M/1即为顾客到达为泊松过程，服务时间为负指数分布，单台的排队系统模型。在1971年的一次国际会议上，将Kendall记号扩充为 ： X/Y/Z/A/B/C。其中前三项意义不变，后三项为 A排队系统的最大容量 B顾客源数量 C排队规则 并约定，如略去后三项，即指 X/Y/Z/FCFS。M/M/1/FCFS，可简
9、写为M/M/1，指顾客到达为泊松过程，服务时间为负指数分布，单台，无限容量，无限源，先到先服务的排队系统模型。,20,3.3 排队论研究的基本问题1.排队系统的统计推断:即判断一个给定的排队系统符合于哪种模型，以便根据排队理论进行研究。 2.系统性态问题:即研究各种排队系统的概率规律性，主要研究队长分布、等待时间分布和忙期分布等统计指标,包括了瞬态和稳态两种情形。3.最优化问题：即包括最优设计(静态优化)，最优运营（动态优化）。,21,求解一般排队系统问题的目的主要是通过研究排队系统运行的效率指标，估计服务质量，确定系统的合理结构和系统参数的合理值，以便实现对现有系统合理改进和对新建系统的最优
10、设计等。排队问题的一般步骤：1. 确定或拟合排队系统顾客到达的时间间隔分布和服务时间分布。2. 研究分析排队系统理论分布的概率特征。3. 研究系统状态的概率。系统状态是指系统中顾客数。状态概率用Pn(t)表示,即在t时刻系统中有n个顾客的概率，也称瞬态概率。,22,求解状态概率Pn(t)方法是建立含Pn(t)的微分差分方程，通过求解微分差分方程得到系统瞬态解，由于瞬态解一般求出确定值比较困难，即便求得一般也很难使用。因此常常使用它的极限(如果存在的话)：,稳态的物理意义图，系统的稳态一般很快都能达到，但实际中达不到稳态的现象也存在。要注意的是求稳态概率Pn并不一定求t的极限,只需求Pn(t)=
11、0 。,称为稳态(steady state)解，或称统计平衡状态 (Statistical Equilibrium State)的解。,23,4.根据排队系统对应的理论模型求用以判断系统运行优劣的基本数量指标的概率分布或特征数。 数量指标主要包括:(1)平均队长（Ls）:系统中的顾客数。平均队列长（Lg）:系统中排队等待服务的顾客数。(2)平均逗留时间(Ws):指一个顾客在系统中的停留时间。平均等待时间(Wg):一个顾客在系统中排队等待的时间。(3)忙期：指从顾客到达空闲服务机构起到服务机构再次为空闲这段时间长度。（忙期和一个忙期中平均完成服务顾客数都是衡量服务机构效率的指标，忙期关系到工作强
12、度）,24,3.4 理论分布,1.泊松分布在概率论中，我们曾学过泊松分布，设随机变量为X，则有：,n=0,1,2, （1）,与时间有关的随机变量的概率，是一个随机过程，即泊松过程。,t0，n=0,1,2, （2）,25,（t2t1，n0）,若设N(t)表示在时间区间0,t)内到达的顾客数(t0),Pn(t1,t2)表示在时间区间t1,t2)(t2t1)内有n(0)个顾客到达的概率。即：,当Pn(t1,t2)符合下述三个条件时，顾客到达过程就是泊松过程(顾客到达形成普阿松流)。,26,平稳性：即对于足够小的t，有：,普阿松流具有如下特性：,在t,t+t内有一个顾客到达的概率与t无关,而与t成正比
13、。,27, 普通性：对充分小的t,在时间区间（t,t+t）内有2个或2个以上顾客到达的概率是一高阶无穷小.,由此知，在(t，t+t)区间内没有顾客到达的概率为：,令t1=0,t2=t,则P(t1,t2)=Pn(0,t)=Pn(t),0 是常数，它表示单位时间到达的顾客数，称为概率强度。,即,P0+P1+P2=1,28,2.负指数分布当输入过程是泊松流时，我们研究两顾客相继到达的时间间隔的概率分布。设T为时间间隔，分布函数为FT（t），则： FT（t）=PTt此概率等价于在0，t）区间内至少有1个顾客到达的概率。,没有顾客到达的概率为： (由（5）式而来),间隔：,间隔：,间隔,对分布函数微分,
14、29,表示单位时间内顾客平均到达数。1/表示顾客到达的平均间隔时间。可以证明，间隔时间T独立且服从负指数分布与顾客到达形成泊松流是等价的。,对顾客的服务时间：系统处于忙期时两顾客相继离开系统的时间间隔，一般地也服从负指数分布，设：,即T服从负指数分布，它的期望及方差为：,接受服务，然后离开,服务时间的分布：,30,其中：表示单位时间内能被服务的顾客数，即平均服务率。 1/表示一个顾客的平均服务时间。,，则,令 ，则称为服务强度。,一般的，要设1（否则队列会无限， 永远达不到稳态）,31,对排队模型，在给定输入和服务条件下，主要研究系统的下述运行指标：(1)系统的平均队长Ls(期望值)和平均队列
15、长Lq(期望值)； (2)系统中顾客平均逗留时间Ws与队列中平均等待时间Wq；本节只研究M/M/1模型.,3.5 M/M/1 模型研究,32,标准的M/M/1模型（M/M/1/FCFS）,系统中有n个顾客,1.稳态概率Pn的计算,在任意时刻t，状态为n的概率Pn(t)（瞬态概率），它决定了系统的运行特征。,已知顾客到达服从参数为的泊松过程，服务时间服从参数为的负指数分布。现仍然通过研究区间t，t+t）的变化来求解。在时刻t+t，系统中有n个顾客不外乎有下列四种情况（ t，t+t）内到达或离开2个以上没列入）。,？,33,由于这四种情况是互不相容的，所以Pn(t+t)应是这四项之和，则有：,所有
16、的高阶无穷小和并,34,当n=0时，只有表中的（A）、（B）两种情况，因为在较小的t内不可能发生（D）（到达后即离去），若发生可将t取小即可。,35,这种系统状态(n)随时间变化的过程就是生灭过程（Birth and Death Process）。方程(1)、(2) 解是瞬态解,无法应用。,它对时间的导数为0，所以由(1)、(2)两式得：,稳态时，,Pn(t)与时间无关,可以写成Pn,36,由此可得该排队系统的状态转移图：,将其代入（3）式并令n=1,2,(也可从状态转移图中看出状态平衡方程)得：,关于Pn的差分方程,37,n=1,n=2,38,以此类推，当n=n时，,(5),及概率性质知：,
17、否则排队无限远,系统稳态概率,39,2. 系统的运行指标计算(1) 系统中的队长Ls（平均队长）,(01),期望,40,(2) 队列中等待的平均顾客数Lq,(8),(3) 顾客在系统中的平均逗留时间Ws顾客在系统中的逗留时间是随机变量，可以证明，它服从参数为-的负指数分布，分布函数,41,(4)顾客在队列中的平均逗留时间Wq,等待时间,顾客在队列中的平均逗留时间应为Ws减去平均服务时间。,考虑LS与WS的关系,42,四个指标的关系为(Little 公式)：,43,知识点评：,从上述分析可以看出：排队论是有用的，同时又非常难。单是一个M/M/1模型就够复杂的了，何况其他模型。此外，对于生活中许多
18、一般化的排队问题，是求不出其概率分布，从而得不到其平均队长、等待时间等标志性参数的，怎么办？这就需要仿真或模拟。模拟的优点是不需要研究、知道模型的内在机理，只是通过多次仿真实验利用概率论知识求得相关参数的近似值（这就够了）以达到问题的解决。下面介绍著名的蒙特卡罗仿真方法。,44,四、蒙特卡罗方法简介,蒙特卡罗(Monte Carlo，美国一个赌城的名称)方法，也称（计算机）随机模拟方法，是一种基于计算机产生的“随机数”的数值计算方法，主要用来解决与随机性有关的一些问题。蒙特卡罗方法的概率论依据：1 用A表示任一随机事件，P(A)，fn(A)分别表示A发生的概率和频率，则当n很大时有P(A)fn
19、(A)。2 设X是一个随机变量（总体）, E(X)=,是X 的期望，X1, X2,Xn, 是来自X的一个样本，则 是 的无偏（有效）估计。,45,随机数的产生,46,仿真与模拟的目的和原理,仿真和模拟可以说是针对同一内容的不同角度的看法描述，当需要对某一问题观察研究而相应的观察和实验时间和成本花费太高时，可以考虑用一个模型代替原型，用模型的研究达到原型的研究的目的（以节约时间和成本），这就是仿真，其在计算机上等的实现过程也称为模拟。,47,在仿真中模型和原型重在本质结构上相同或相似（否则没有意义），不必和不能追求它们外型的相似或相同。理解了这一点，容易掌握蒙特卡罗方法。下面通过几个例子介绍蒙特
20、卡罗方法的特点和应用原理。,48,例1：掷硬币仿真,大家知道，掷硬币有两个结果：正面朝上（用A表示）和正面朝下（用表示），其概率都是0.5。掷硬币实验就是多次将硬币抛掷以观察正面与反面出现的次数。仿真不用真的抛掷硬币，而是构造一个随机变量X，其在某范围内取值的概率是0.5 （作为事件A，对立事件是 ），用计算机产生随机数，根据随机数判断A是否发生（认为硬币正面是否出现），从而达到掷硬币实验的目的。,五：仿真例子与分析,49,最简单的，令XU(0,1)，则P(X1/2)=0.5,（X1/2表示正面出现），用计算机产生服从U(0,1) 分布的一些随机数Xi,若Xi1/2, 认为X1/2发生（硬币正
21、面出现）了，将结果记录与累加就成了掷硬币仿真，程序如下,clc n=1000 %test number m=0; % face frequency for i=1:nr=rand;if r0.5m=m+1;end end m rate=m/n %face rate,50,例2：蒙特卡罗法求的（近似）值,求的值已有多种方法，而且要多精确都可以。蒙特卡罗方法求的值效果并不好，这里主要介绍方法。蒙特卡罗方法可以解决随机性问题，也可解决非随机性问题。当解决随机性问题时，应使模拟与原问题一致。当求解非随机性问题时，通常是构造一个概率模型，使所要求解的数值是该模型的一个参数，而后进行模拟，根据模拟结果得参
22、数的近似值，将问题解决。下面通过本例说明。,51,模型构造：设XU（-1，1）， YU（-1，1），是两个随机变量，且相互独立，则（X，Y）在如右图所示的正方形内服从均匀分布。今考虑事件X2+Y21，该事件就相当于随机变量（X，Y）落在圆周内。利用均匀分布的特征容易得到,从而 4fn(A)，程序如下,52,n=10000; % test number m=0; for i=1:nx=unifrnd(-1,1); %generate the rand number of xy=unifrnd(-1,1); %generate the rand number of yif x2+y21 % jud
23、ge if (X,Y) is in the circlem=m+1; % add the frequencyend end fn=m/n % rate yuanzhoulv=4*fn %result,53,例3：蒙特卡罗法模拟中子穿过原子反应堆屏障问题,原子反应堆外的铅屏障是用来屏障原子反应堆中逸出的中子的，以免给人类造成危害。经试验和分析，可做以下简单假设：假设每一个进入屏障的中子在撞到铅原子前行进的距离为D，然后这个中子以随机方向弹回来，并且在它的下一次撞击中又行进距离D。假设屏障厚度为3D，每一个中子能经受10次撞击，问进入的中子能以多大的比例穿透铅屏障？,54,分析：该问题显然难以用定
24、性理论（概率论）解决（求概率），用蒙特卡罗方法却很容易（多次模拟,用频率代替概率）。为方便不妨做进一步假设：1 中子反弹回反应堆认为不能穿过屏障。2 与铅原子相撞后任意方向等可能反弹。画图如右，模拟流程图如下,55,中子撞击铅屏模拟流程图,初始化系统状态,产生一个新中子的初试方向和运行终点,中子回到 反应堆了吗,求频率,结束,Y,N,碰撞，产生新方向和运行终点,模拟次数 到了吗,N,Y,中子出了 铅屏了吗,碰撞次数到了吗,N,频数增加,Y,Y,N,56,程序为,n=10000; % test number m=0; %frequency for i=1:ntheta=unifrnd(-pi,p
25、i); % initinal anglex=cos(theta); % only x needed for j=1:10theta=2*pi*rand;%new angle hittedx=x+cos(theta);if x3m=m+1;break;endend end fn=m/n % rate,57,例4： M/M/1模型及其模拟流程图,M/M/1排队模型即是顾客到来的时间间隔和服务时间都服从负指数分布，只有一个服务台（员），队长无限（来了就排队）。我们给出适当参数计算顾客的平均等待时间和平均逗留时间，模拟流程图如下,58,M/M/1模型流程图,初始化系统状态,第一个人接受服务,有人等吗,
26、计算相关结果 结束,时间推到该人服务结束 同时他离开系统,Y,时间推到下一个顾客到来,服务人数 到了吗,N,服务中会有 人来吗,产生新到达顾客,N,N,Y,Y,59,几点说明,1 一般模拟均值，实验次数都上万，但显示窗口只显示后面部分，为方便把模拟次数减少。 2 为了直观，我们给出了部分顾客逗留时间的散点图，它们可以没有。 3 当模拟 100次得 Wq= 9.0459, Ws= 12.5612,而当模拟 10000次得 Wq= 10.1294, Ws= 13.1731,而从理论上应该是Ws=3/4/(1/3-1/4)=9, Wq=1/(1/3-1/4)=12，它们还是接近的。,60,作业,1
27、两人约定下午1-3点在某处见面，设两人的到达分别服从分布U(1,3)和N(2,0.5)，且一人最多等另一人15分钟，求两人见面的概率。 2 编程模拟顾客到达服从负指数分布（1/=10），服务时间分别服从负指数分布（=8）和正态分布（均值为9，方差为4）的排队过程，观察模拟排队过程并求顾客平均等待时间和逗留时间。 3（选做）一列火车大约在下午1点离开A站，规律如下：离站时间 13：00 13：05 13：10概率 0.7 0.2 0.1 火车从A到B途中所需的平均时间为30分，有2分钟的标准差。如果你要赶的是这趟火车的下一站B，而你到达的时间分布为时间 13：28 13：30 13：32 13：34概率 0.3 0.4 0.2 0.1 问你能赶上这列火车的概率是多少？,
