第八章运输系统仿真

第八章

运输系统仿真

什么是系统仿真？------概念!

为什么要系统仿真？------作用!

如何进行系统仿真？------方法!本章

问题Contents系统仿真概念系统仿真概述系统仿真应用领域、分类蒙特卡洛模拟MC系统仿真随机数离散事件系统仿真§8-1系统仿真概述

1.系统仿真的概念

系统仿真（亦称系统模拟，SystemSimulation）是利用模型对实际系统进行试验研究的过程。就是根据系统分析的目的，在分析系统各要素性质及其相互关系的基础上，建立能描述系统结构或行为过程的、且具有一定逻辑关系或数量关系的仿真模型，据此进行试验或定量分析，以获得正确决策所需的各种信息。通过模型对系统时——空变化作动态的描述。这种模型称为模拟模型，它是某一系统的某些重要因素的外貌或形式的抽象，而不是实际存在事物的复制。模拟法涉及大量逻辑的、数字的和薄记的工作，比较复杂的作业非借助电子计算机不可。系统仿真是20世纪40年代末以来伴随着计算机技术的发展而逐步形成的一门新兴学科。仿真（Simulation）就是通过建立实际系统模型并利用所建模型对实际系统进行实验研究的过程。最初，仿真技术主要用于航空、航天、原子反应堆等价格昂贵、周期长、危险性大、实际系统试验难以实现的少数领域，后来逐步发展到电力、石油、化工、冶金、机械等一些主要工业部门，并进一步扩大到社会系统、经济系统、交通运输系统、生态系统等一些非工程系统领域。2.系统仿真应用领域不易控制实验的动态系统模拟模型的最大优点就是它的可实验性，它比实际系统更容易改变条件。凡无法用数学分析方法求解的复杂问题，或不易控制实验的动态系统，其输入大多是不确定的、可能含有许多随机变量和参数，且其转换过程又往往含有许多函数关系，而此种动态函数关系又非数学方程式所能表达者，则只能用模拟法试验求解。现实不允许进行有效试验的系统即使分析人员有信心和能力对一个大系统的行为，从理论上进行推测，也许现实又不允许进行有效的试验，例如，人类决不允许对全球战略性战争的结论进行试验。在这种情况下，唯一的办法是进行某种形式的模拟。系统太复杂有时候由于现象或系统太复杂，不可能用一个简单的数学模型来表达。例如，不可能把大规模的政治、经济或社会活动简化为一个数学方程式。而且，并不是所有的数学方程式都可以用分析方法求解的。2.系统仿真应用领域3.系统仿真的分类

实施一项系统仿真的研究工作，包括三个基本要素，即系统对象、系统模型以及计算机。系统仿真一般可以分为连续系统模拟和离散事件系统模拟。连续系统是指系统状态随时间连续变化的系统，系统行为通常是一些连续变化的过程。3.系统仿真的分类离散事件系统模拟指表征系统性能的状态只在随机的时间点上发生跃变，且这种变化是由随机事件驱动的，在两个时间点之间，系统状态不发生任何变化。例如，医院门诊病人数量、路口车辆通过数量、公共汽车上乘车人数的变化。还有一些用于系统(特别是社会经济和管理系统)仿真的特殊而有效的方法，如系统动力学方法、蒙特卡洛法等。系统动力学方法通过建立系统动力学模型(流图等)、利用DYNAMO仿真语言在计算机上实现对真实系统的仿真实验，从而研究系统结构、功能和行为之间的动态关系。4.系统仿真的功能评估系统各部分或子系统彼此之间的影响和对系统整体性能的影响比较各种设计方案，获得最优的设计在系统发生故障后，使之重演，以便研究故障原因进行假设检验培训系统操作人员问题描述与定义建立系统模拟模型数据采集模拟模型的验证对模拟结果的评价5.系统仿真的基本步骤

6.MC系统仿真随机数

随机数是专门的随机试验的结果。在统计学的不同技术中需要使用随机数，比如在从统计总体中抽取有代表性的样本的时候，或者在将实验动物分配到不同的试验组的过程中，或者在进行蒙特卡罗模拟法计算的时候等等。

随机数定义：(1).随机数生成器产生(2).计算机程序(如用C语言的函数rand()产生随机数)(3).随机数表

产生随机数有多种不同的方法。这些方法被称为随机数发生器。随机数最重要的特性是：它所产生的后面的那个数与前面的那个数毫无关系。真正的随机数是使用物理现象产生的：比如掷钱币、骰子、转轮、使用电子元件的噪音、核裂变等等。这样的随机数发生器叫做物理性随机数发生器，它们的缺点是技术要求比较高。在实际应用中往往使用伪随机数就足够了。这些数列是“似乎”随机的数，实际上它们是通过一个固定的、可以重复的计算方法产生的。计算机或计算器产生的随机数有很长的周期性。它们不真正地随机，因为它们实际上是可以计算出来的，但是它们具有类似于随机数的统计特征。这样的发生器叫做伪随机数发生器。在真正关键性的应用中，比如在密码学中，人们一般使用真正的随机数。

随机数生成

随机数表例1:在50名同学（学号20—69）中随机抽出五位同学参加某项活动。用“随机数”模拟抽样的方法在excel中输入=RAND()*49+20

例2:设某商店每天到达的顾客数是随机的，一个月的统计结果表明，一天最多有60人，最少有25人。现模拟今后一个星期内平均每天有多少顾客。假定每天到达的顾客多少的机会相同，即不存在月初、月中或月尾及休息天的差别。在excel中输入=RAND()*35+25

求出模拟事件出现的概率计算累计频率将累计频率换算为随机概率抽取随机数将随机数作为事件出现的概率进行模拟某汽车修理车间修理每辆汽车所需要的时间进行统计，共统计了100次，各种时间出现的次数如表所示，试对汽车修理车间修理每辆汽车的服务时间进行模拟。例题3服务时间（h)发生次数分布频率累计频率随机数6100.107250.258350.359200.2010100.10服务时间（h)发生次数分布频率累计频率随机数6100.100.100.00-0.097250.250.350.10-0.348350.350.700.35-0.699200.200.900.70-0.8910100.101.000.90-0.99在excel中输入=RAND()*0.99+0.000.010.640.440.530.090.850.550.610.020.750.260.350.350.550.530.490.580.660.556888698869788888888随机数修理时间§8-2蒙特卡罗模拟蒙特·卡罗方法（MonteCarlomethod），也称统计模拟方法、随机模拟方法，是二十世纪四十年代中期由于科学技术的发展和电子计算机的发明，而被提出的一种以概率统计理论为指导的一类非常重要的数值计算方法。是指使用随机数（或更常见的伪随机数）来解决很多计算问题的方法。蒙特·卡罗方法的名字来源于摩纳哥的一个城市蒙地卡罗，该城市以赌博业闻名，而蒙特·卡罗方法正是以概率为基础的方法。原理：当所求解问题是某种随机事件出现的概率，或者是某个随机变量的期望值时，通过某种“实验”的方法，以这种事件出现的频率估计这一随机事件的概率，或者得到这个随机变量的某些数字特征，并将其作为问题的解。蒙特卡罗方法研究的问题大致可分为两种类型，一种是问题本身就是随机的，另一种本身属于确定性问题，但可以建立它的解与特定随机变量或随机过程的数字特征或分布函数之间的联系，因而也可用随机模拟方法解决。如计算多重积分，求解积分方程、微分方程、非线性方程组，求矩阵的逆等。例4:假设我们要计算一个不规则图形的面积，那么图形的不规则程度和分析性计算（比如:积分）的复杂程度是成正比的。蒙特卡罗方法是怎么计算的呢？首先作一个包含复杂图形的矩形D，假想你有一袋豆子，把豆子均匀地朝这个矩形上撒，然后数矩形上和复杂图形之中有多少颗豆子，复杂图形之中豆子的数目比矩形上豆子的数目再乘以矩形D的面积就是复杂图形的面积。当你的豆子越小，撒的越多的时候，结果就越精确。在这里我们要假定豆子都在一个平面上，相互之间没有重叠。例5.计算π的近似值原理：

一个边长为2的正方形，面积也就是4，中间有一个半径为1的圆，按照公式，圆的面积是π

，那么，模拟正方形内N个均匀分布的随机点，如果落在圆内就编码为1，如果落在圆外的正方形内，就编码为0，然后所有编码值加起来除以N，就是圆的面积占正方形的面积的近似比例，也等于π

/4，这样，N越大，算出来的π值越精确，模拟了1000个点，刚好算出来π值等于3.14，因为将这1000个随机点距离圆心的距离排序之后，第785个值开始大于1，意味着786之后的点都落在圆外，785×4/1000=3.14。如果模拟更多的点，就可以得到更加精确的π值。例6．某理发店有3个座位(1号位、2号位、3号位)，每5分钟进来的顾客数为0，1或2位，且进来的概率相同（都是1/3）。顾客在每个理发椅上的服务时间是15或20分钟，且概率也相同（都是1/2），现模拟1小时内该店顾客排队或理发师空闲情况。蒙特卡罗(MC)仿真方法应用

仿真：依所要求的概率分布产生的随机数来模拟可能出现的随机现象。

第二步产生仿真变量的随机数得到仿真量：第一步确定仿真变量的概率分布：第三步仿真（模拟）座椅被占用的情况：第四步理发店营业情况分析。例7．某工厂从外地采购原料，到货天数是一个随机变量（设为

X)。根据过去的资料，在100次到货中，到货天数与次数的关系如表1到货天数X2357812次数204082552现模拟今后10批货物到达的平均天数解：①根据已知条件，到货天数X的概率见表到货天数X2357812概率P0.200.400.080.250.050.02②变换：③产生均匀分布的随机数：从随机数表中按行读出10个随机数：

68、34、30、13、70、55、74、30、77、40

④10批货物平均到货天数:(7+3+3+2+7+3+7+3+7+3)/7到货天数X2357812概率P0.200.400.080.250.050.02对应随机数00～1920～5960～6768～9293～9798～99例8：12132解：①根据已知条件，每天销售量X与到货天数T的概率见表3。每天销售量X概率P对应的随机数每天销售量X概率P对应的随机数700.0400～03950.1440～53750.0404～071000.1954～72800.0908～161050.1473～86850.0917～251100.0987～95900.1426～391200.0496～99

到货天数T2346812概率P0.170.250.330.170.040.04对应随机数00～1617～4142～7475～9192～9596～99产生X的均匀分布随机数：从随机数表第21行起读取30个随机数：②变换：对应的随机数对应的销售量：100、90、90、80、100、100、105、90、…对应的到货天数：4、3、3、2、4、4、4、…

③仿真：④分析：

日期123456库存量30020011020240140订货量300300到货量300300销售量100909080100100日期123456库存费用300×0.2=60200×0.2=40110×0.2=2220×0.2=4240×0.2=48140×0.2=28销售利润100×20=200090×20=180090×20=180080×20=1600100×20=2000100×20=2000净利润2000-60=19401800-40=17601800-22=17781600-4=15962000-48=19522000-28=1972

注意：

应继续模拟，比如：

定货量300台/次改成350台/次或250台/次看其效果;改变月初存货量看其效果;

选择一个较好的可能方案。一般排队系统的模拟考虑发展因素排队问题用随机概率模拟排队论问题一般排队系统的模拟计算排队系统的随机概率排队系统模拟排队系统运行指标计算排队系统分析某运输公司仓库有一个管理员负责发料工作，根据过去的记录，得知领料人的达到时间间隔和管理员发料时间长度的频率如下表所示。例8：领料人达到时间间隔到达时间间隔345678频率0.050.200.350.250.100.05发料时间长度34567频率0.100.200.400.200.10试用蒙特卡洛模拟求：领料人的平均等待时间等待行列的平均顾客人数发料人的平均服务时间领料人平均达到时间领料人在仓库的平均消耗时间领料人到达时间间隔累计频率及随机概率到达间隔频率累计频率随机数30.050.050.00-0.0440.200.250.05-0.2450.350.600.25-0.5960.250.850.60-0.8470.100.950.85-0.9480.051.000.94-0.99发料时间长度频率累计频率随机数30.100.100.00-0.0940.200.300.10-0.2950.400.700.30-0.6960.200.900.70-0.8970.101.000.90-0.99管理员服务时间长度累计频率及随机概率在excel中输入=RAND()*0.99+0.00随机数到达间隔随机数发料时间长度队列长度0.4850.80610.1740.45510.3350.50510.6660.12400.3850.03300.9270.6050302520151051离散事件系统仿真

离散事件系统仿真的基础

基本概念

1.实体（Entity）：一个系统边界内部的客观对象称为实体。•永久实体是指经常处于系统之内，其数量保持稳定的实体。如排队系统中的服务员，交通系统中的道路设施等。•临时实体（Transaction）又称主动实体、活动实体,是指先进入系统并经过相应的环节以后再离开系统，在系统中的数量经常变化的实体。2．属性（Attributes）是指系统实体的特性,它可以是文字型，数字型或逻辑型。3．活动（Activity）活动是占用一定时间和资源导致系统状态发生改变的一定过程。4.系统状态（Systemstate）是描述在某时间点系统的所有实体、属性和活动的变量集合。5.事件（Event）是引起系统状态发生瞬间变化的事实。•时间事件是依照系统的作业规则在预定时间发生的事件。•状态事件是当系统状态符合某种条件而发生的事件。6．时钟（Clock）研究系统一般是为认识其状态随时间变化的规律,所以需要一个仿真的时间变量，由于计算机仿真是一种数值的方法，它不提供解析解，所以对连续的系统进行仿真的时候，常在均匀的时间点上展现其状态值，这样仿真时钟的步进是一个常数，对离散的系统仿真的时候，只有事件发生的时候，系统的状态才发生变化，才有必要展现出系统的状态。7．统计计数器用来统计仿真次数。8．进程（Process）进程由若干个事件及若干活动组成，一个进程描述了它所包括的事件及活动间的相互逻辑关系及时序关系。9．系统环境（Environment）存在于系统周围的对象和过程（即实体和活动）称为系统的环境。•产生于系统内部的活动称为内生活动；•发生于系统环境并对系统产生影响的活动称为外生活动。离散事件系统仿真算法事件调度法（EventScheduling）基本思想离散事件系统中最基本的要素是事件，事件发生引起系统状态发生变化。用事件的观点来分析系统，通过定义事件及每个事件发生引起系统状态的变化，按时间顺序确定并执行每个事件发生时有关的逻辑关系。活动扫描法（ActivityScanning）基本思想•系统由实体组成，实体要产生活动，它标志系统状态的转移。•活动的发生必须满足某些条件，在仿真过程中，活动的发生时间也作为条件之一，而且比其他条件更有优先权。•在以活动扫描为基础的仿真模型中，要描述系统实体所进行的活动以及预定导致活动开始或结束的条件。活动开始或结束的事件并非由模型事先安排，而是由活动规定的条件所初始化，随着仿真时钟按一定步长推进，对每项活动开始或终止的条件进行扫描。如果预定的条件满足的话，那么相应的活动就得以进行。为了保证顾及每项活动，必须在每个时间点对整个活动集合进行扫描。进程交互法（ProcessInteraction）基本思想进程交互法采用进程描述系统，它将模型中主动实体历经系统时所发生的事件及活动按时间顺序进行组合，从而形成进程表，一个实体一旦进入进程，它将完成该进程的全部活动。例9:混凝土搅拌中心的位置

有K1(x1,y1),K2(x2,y2),....,Kn(xn,yn)共n个工地，各需混凝土Q1,Q2,...,Qn

吨，混凝土每吨公里的运费为C元。如何确定混凝土搅拌中心的位置K0(x0,y0)使得费用最少？解1：数学计算（解析解）记第K个工地的位置为（xk,yk）,中心的位置为(x0,y0)，则目标函数为关于求偏导，并令其等于零，计算方程组的解即为所求的解析解。解2：计算机的模拟。将平面区域划分成许多小方格，利用计算机循环计算搅拌中心设立在每个交叉点的费用，再经过比较就可以得到最佳位置。例10：捕捞队准备购置冷藏车运输所捕捞的鱼产品,有关资料如下，用模拟方法进行方案选择:资料统计,平均日捕鱼量为180吨,标准差为30吨,预计计划期捕鱼量可增加5%;每台冷藏车日运输量为20吨;捕获的鲜鱼当日不能运完,每吨损失2000元,冷藏车空驶,运输费用损失每吨100元;捕捞队初步提出购置9台、10台、11台冷藏车方案。将计划期发展因素加入到实际事件的平均数与标准差中:

σ=30*189/180=31.5选取均匀随机数:序号均匀随机数172,95,17,31,16,93,32,43,50,27289,87,19,20,15,37,00,49,32,85366,60,44,38,68,88,11,80,68,34430,13,70,55,74,30,77,40,44,22578,84,26,04,33,46,09,52,68,07将均匀随机数改造为常态随机数:序号常态随机数148,52,35,51,55,32,51,38,60,62259,48,53,53,60,55,48,48,40,42338,67,53,46,40,47,51,39,41,24450，49，42求常态随机数的平均数与标准差:3163.68160533∑556.96924.80155.52450.561075.84556.96184.9612.9640.96268.9623.613.6-3.66.416.42517554060526015121142535455565203040506070fm频数f组中值m分组求事件标准差与常态平均数标准差的比值e和模拟修正数Ee=31.5/9.79=3.22E=48.6*3.22-189=-32.5求模拟随机数m1=48*3.22-（-32.5）=187.06≈187m2=52*3.22-（-32.5）=199.94≈200序号模拟随机数1187，200，145，197，210，138，197，155，226，2322222，187，203，203，226，210，187，187，161，1683155，248，203，181，161，184，197，158，165，1104194，190，168模拟过程方案评价结论：10台购车方案总损失最小，为最优方案。报童问题－库存模型

报童从分销商处购入《金融时报》，每份购入价为$1.5，售出价为$2.5，若未售出，则分销商以每份$0.5回购。最佳的定购数量是多少？

目标：每天平均利润最大每天销售量在40-70份之间在40-70之间的频率大致相同市场需求状况候选策略：每天定购60份

报童问题模型（Excel）Profit=SalesRevenue-PurchasingCost+SalvageValueFrequencyChart在Preferences\Chart中选择显示均值线求输出在某一范围内的概率利润大于40的概率为68.40%

仿真精度控制的结果

仿真结果报童有95%的置信度相信，定购数量为60时，从长远看能获得的日均利润在45.65到47.57美元之间（46.610.96）。项目管理工期：44周P(T47)=0.84有什么问题？

项目管理问题的仿真模型（Excel）FrequencyChartPercentilesTable从百分比图可看出，有30%的试验次数超过了截止日期（47周）一周。现欲了解项目进度风险（在47周以内完工的概率）。FrequencyChart按期完工的概率为60.2%按PERT三点估计法所得到的概率为84%，故过于乐观SensitivityAnalysis时间步长法时间步长法就是按照时间流逝的顺序，一步一步的对系统的活动进行仿真。在整个仿真的过程中，时间步长的长度固定不变。它的基本思路是：在进行系统仿真的过程中，可以把整个过程分成许多相等的时间间隔，时间步长的长度可以根据实际问题分别取作秒、分、小时、天等。程序中按照这个步长前进的时钟就是仿真的时钟。选取系统的一个初始的状态作为仿真时钟的零点，仿真时钟每步进一次，就对系统的所有的实体和属性以及活动进行一次全部的扫描考察，按照预定的计划和目标进行分析，计算和记录系统状态的变化，这个过程一直进行到仿真的时钟结束为止。其流程图为：

某水池中有2000的水，其中含盐2kg。以每分钟6的速度向水池内注入含盐率为0.5kg/的盐水，同时又以每分钟4的速度从水池中流出搅拌均匀的盐水。使用计算机仿真该水池里盐水的变化过程，并每隔10分钟计算水池中水的体积、含盐量、含盐率，欲使池中盐水的含盐率到达0.2kg/,需经过多少时间？例11池水含盐量问题池水含盐量随时间变化的示意图

理论分析法

当系统开始运行后，池中的水和含盐量均随时间变化。记x(t)为t时刻池水中含盐量，考察[t,t+Δt]时间区间内盐的变化，这些变化是由于流入的盐水和流出的盐水所引起的，于是有平衡方程:x(t+Δt)-x(t)=[6×0.5－4×x(t)/(2000+2t)]Δt两边同时除以Δt,再令Δt趋于0得下面微分方程的初始值问题。dx(t)/dt=3-4x(t)/(2000+2t),x(0)=2.解和图形含