计算机仿真.ppt

《计算机仿真.ppt》由会员分享，可在线阅读，更多相关《计算机仿真.ppt（58页珍藏版）》请在一课资料网上搜索。

1、 计算机仿真、随机模拟计算机仿真、随机模拟基本内容：基本内容：1 1、计算机仿真的基本概念、计算机仿真的基本概念2 2、随机数的产生、随机数的产生3 3、时间步长法、时间步长法4 4、事件步长法、事件步长法6 6、实验作业、实验作业5 5、Monte Carlo方法方法一、计算机仿真的基本概念一、计算机仿真的基本概念仿真：仿真：就是将所研究的对象用其它手段加以模仿就是将所研究的对象用其它手段加以模仿的一种活动。的一种活动。实物仿真实物仿真非实物仿真非实物仿真如曹冲称象、军事演习、如曹冲称象、军事演习、风洞试验等风洞试验等计算机仿真：计算机仿真： 是一种非实物仿真方法，是一种非实物仿真方法，通过

2、建立数学模型、通过建立数学模型、编制计算机程序实现对真实系统的模拟。从而了解编制计算机程序实现对真实系统的模拟。从而了解系统随时间变化的行为或特性系统随时间变化的行为或特性, , 以评价或预测一个以评价或预测一个系统的行为效果系统的行为效果, ,为决策提供信息的一种方法为决策提供信息的一种方法. .它是它是解决较复杂的实际问题的一条有效途径解决较复杂的实际问题的一条有效途径。 鼠疫的检测和预报，三峡的安全、生态，公鼠疫的检测和预报，三峡的安全、生态，公交车的调度交车的调度, ,航空管理航空管理, , 经营投资，道路的修建经营投资，道路的修建, ,通信网络服务，通信网络服务，电梯系统服务电梯系统

3、服务等等。等等。 计算机仿真举例：计算机仿真举例： 难以用数学公式表示的系统，或者没有建立和求解难以用数学公式表示的系统，或者没有建立和求解数学模型的有效方法数学模型的有效方法 虽然可以用解析的方法解决问题，但数学的分析与虽然可以用解析的方法解决问题，但数学的分析与计算过于复杂，这时计算机仿真可能提供简单可行计算过于复杂，这时计算机仿真可能提供简单可行的求解方法的求解方法 希望能在较短的时间内观察到系统发展的全过程，希望能在较短的时间内观察到系统发展的全过程，以估计某些参数对系统行为的影响以估计某些参数对系统行为的影响 难以在实际环境中进行实验和观察时，计算机仿真难以在实际环境中进行实验和观察

4、时，计算机仿真是唯一可行的方法，例如太空飞行的研究是唯一可行的方法，例如太空飞行的研究 需要对系统或过程进行长期运行比较，从大量方案需要对系统或过程进行长期运行比较，从大量方案中寻找最优方案中寻找最优方案适用计算机仿真解决的问题：适用计算机仿真解决的问题： 便于重复进行试验，便于控制参数，时间短，便于重复进行试验，便于控制参数，时间短， 代价小；代价小；计算机仿真的优点：计算机仿真的优点： 可以在真实系统建立起来之前，预测其行为效可以在真实系统建立起来之前，预测其行为效果，从而可以从不同结构或不同参数的模型的果，从而可以从不同结构或不同参数的模型的结果比较之中，选择最佳模型；结果比较之中，选择

5、最佳模型； 对于缺少解析表示的系统，或虽有解析表示但对于缺少解析表示的系统，或虽有解析表示但无法精确求解的系统，可以通过仿真获得系统无法精确求解的系统，可以通过仿真获得系统运行的数值结果；运行的数值结果； 对于随机性系统，可以通过大量的重复试验，对于随机性系统，可以通过大量的重复试验，获得其平均意义上的特性指标。获得其平均意义上的特性指标。仿真常用术语仿真常用术语:系统：系统： 一些具有特定功能相互之间以一定的规律联一些具有特定功能相互之间以一定的规律联系着的物体所组成的总体系着的物体所组成的总体系统的对象、系统的组成元素都可以称为实体系统的对象、系统的组成元素都可以称为实体.实体：实体： 属

6、性是对实体特征的描述，可以是文字型、属性是对实体特征的描述，可以是文字型、数字型或逻辑型数字型或逻辑型属性属性: 实体在一段时间内持续进行的操作或活动实体在一段时间内持续进行的操作或活动活动活动:例如：例如：银行系统，销存系统银行系统，销存系统 改变系统状态的瞬间变化的事情改变系统状态的瞬间变化的事情.事件事件: 为了限制所研究问题涉及的范围为了限制所研究问题涉及的范围, , 用系用系统边界把所研究的系统与影响系统的环境区分开来统边界把所研究的系统与影响系统的环境区分开来. .系统边界：系统边界： 事件表一般是一个有序的记录列，每个记事件表一般是一个有序的记录列，每个记录包括事件发生时间、事件

7、类型等一些内容录包括事件发生时间、事件类型等一些内容事件表事件表: 系统的状态是指在某一时刻实体及其属性系统的状态是指在某一时刻实体及其属性值的集合值的集合状态状态: 表示仿真时间的变量表示仿真时间的变量仿真时钟仿真时钟:仿仿真真研研究究步步骤骤问题的阐述设置目标及完整的项目研究计划建立模型收集数据编程序程序验证模型确认试验设计运行与分析进一步运行仿真结束输出结果是是是是否否否否系统分析系统分析模型构造模型构造模型运行模型运行输出结果输出结果明确问题和提出总体方案明确问题和提出总体方案: :把被仿真系统的内容表达清楚；把被仿真系统的内容表达清楚；弄清仿真的目的、系统的边界；弄清仿真的目的、系统

8、的边界；确定问题的目标函数和可控变量；确定问题的目标函数和可控变量；找出系统的实体、属性和活动等。找出系统的实体、属性和活动等。系统分析建立模型；建立模型； 选择合适的仿真方法（如时间步长法、事选择合适的仿真方法（如时间步长法、事件表法等）；确定系统的初始状态；设计整个件表法等）；确定系统的初始状态；设计整个系统的仿真流程图系统的仿真流程图。收集数据；收集数据；编写程序、程序验证；编写程序、程序验证；模型确认模型确认。模型构造运行：运行：确定具体的运行方案确定具体的运行方案, ,如初始条件、如初始条件、参数、步长、重复次数等参数、步长、重复次数等, ,然后输入数据然后输入数据, ,运运行程序。

9、行程序。改进：改进：将得出的仿真结果与实际系统比较将得出的仿真结果与实际系统比较, ,进一步分析和改进模型进一步分析和改进模型, ,直到符合实际系统直到符合实际系统的要求及精度为止的要求及精度为止。模型的运行与改进 设计出结构清晰的仿真结果输出。包括设计出结构清晰的仿真结果输出。包括提供文件的清单提供文件的清单, ,记录重要的中间结果等。记录重要的中间结果等。 输出格式要有利于用户了解整个仿真过输出格式要有利于用户了解整个仿真过程程 , ,分析和使用仿真结果分析和使用仿真结果. .设计格式输出仿真结果 (库存问题库存问题) 某电动车行的仓库管理人员采取一种简单的订某电动车行的仓库管理人员采取一

10、种简单的订货策略，当库存量降低到货策略，当库存量降低到P辆电动车时就向厂家订辆电动车时就向厂家订货，每次订货货，每次订货Q辆，如果某一天的需求量超过了库辆，如果某一天的需求量超过了库存量，商店就有销售损失和信誉损失，但如果库存存量，商店就有销售损失和信誉损失，但如果库存量过多，会导致资金积压和保管费增加。若现在已量过多，会导致资金积压和保管费增加。若现在已有如下表所示的两种库存策略，试比较选择一种策有如下表所示的两种库存策略，试比较选择一种策略以使总费用最少。略以使总费用最少。 计算机仿真举例：计算机仿真举例： 这个问题的已知条件是：这个问题的已知条件是：(1)从发出订货到收到货物需隔从发出订

11、货到收到货物需隔3天。天。(2)每辆电动车保管费为每辆电动车保管费为0.50元元/天，每辆电动车的缺货天，每辆电动车的缺货损失为损失为1.60元元/天，每次的订货费为天，每次的订货费为75元。元。(3)每天电动车需求量是每天电动车需求量是0到到99之间均匀分布的随机数。之间均匀分布的随机数。(4)原始库存为原始库存为110辆，并假设第一天没有发出订货。辆，并假设第一天没有发出订货。分析：分析：这一问题用解析法讨论比较麻烦，但用计算这一问题用解析法讨论比较麻烦，但用计算机按天仿真仓库货物的变动情况却很方便。我们以机按天仿真仓库货物的变动情况却很方便。我们以30天为例，依次对这两种方案进行仿真，最

12、后比较天为例，依次对这两种方案进行仿真，最后比较各方案的总费用，从而就可以做出决策。各方案的总费用，从而就可以做出决策。 计算机仿真时的工作流程是早上到货、全天计算机仿真时的工作流程是早上到货、全天销售、晚上订货，以一天为时间步长进行仿真。销售、晚上订货，以一天为时间步长进行仿真。 首先检查这一天是否为预定到货日期，如果是，首先检查这一天是否为预定到货日期，如果是，则原有库存量加则原有库存量加Q，并把预定到货量清为零；如果，并把预定到货量清为零；如果不是，则库存量不变。不是，则库存量不变。（程序附后）（程序附后） 接着仿真随机需求量，这可用计算机语言中的接着仿真随机需求量，这可用计算机语言中的

13、随机函数得到。若库存量大于需求量，则新的库存随机函数得到。若库存量大于需求量，则新的库存量减去需求量；反之，则新库存量变为零，并且要量减去需求量；反之，则新库存量变为零，并且要在总费用上加缺货损失。在总费用上加缺货损失。 然后检查实际库存量加上预定到货量是否小于然后检查实际库存量加上预定到货量是否小于重新订货点重新订货点P，如果是，则需要重新订货，这时就加，如果是，则需要重新订货，这时就加一次订货费。如此重复运行一次订货费。如此重复运行30天，即可得所需费用天，即可得所需费用总值。总值。 由此比较这两种方案的总费用，可以得到最好由此比较这两种方案的总费用，可以得到最好的方案。的方案。1 1均匀

14、分布的随机数及其产生均匀分布的随机数及其产生 对随机现象进行模拟，实质上是要给出随机变对随机现象进行模拟，实质上是要给出随机变量的模拟也就是说利用计算机随机地产生一系列量的模拟也就是说利用计算机随机地产生一系列数值，它们的出现服从一定的概率分布数值，它们的出现服从一定的概率分布, ,则称这些数则称这些数值为值为随机数。随机数。 最常用的是在最常用的是在(0,1)(0,1)区间内均匀分布的随机数，区间内均匀分布的随机数，也就是我们得到的这组数值可以看作是也就是我们得到的这组数值可以看作是(0,l)(0,l)区间内区间内均匀分布的随机变量的一组独立的样本值均匀分布的随机变量的一组独立的样本值其它分

15、布的随机数可利用均匀分布的随机数产生其它分布的随机数可利用均匀分布的随机数产生二、随机数的产生二、随机数的产生产生模拟随机数的计算机命令产生模拟随机数的计算机命令 在在Matlab软件中，可以直接产生满足各种分布的软件中，可以直接产生满足各种分布的随机数，命令如下：随机数，命令如下： 1、产生、产生m*n阶阶(0,1)均匀分布的随机数矩阵：均匀分布的随机数矩阵：rand (m, n) 产生一个产生一个(0,1)均匀分布的随机数：均匀分布的随机数：rand2、产生、产生m*n阶阶(a, b)均匀分布均匀分布U(a, b)的随机数矩阵：的随机数矩阵： unifrnd (a, b, m, n) 产生

16、一个产生一个(a, b)均匀分布的随机数：均匀分布的随机数：unifrnd(a, b)randn(m,n)： m*n阶阶N(0,1)标准正态分布随机数矩阵标准正态分布随机数矩阵3、正态分布、正态分布 )n,m,(normrnd： To Matlab（rnd） ),(normrnd： 若连续型随机变量若连续型随机变量X 的概率密度函数为的概率密度函数为 其中其中 0为常数，则称为常数，则称X 服从服从参数参数为为 的指数分布。的指数分布。 000 xxe)x(fx 指数分布的期望值为指数分布的期望值为 1电子元件的寿命，电话的通话时间，微生物的寿命，随电子元件的寿命，电话的通话时间，微生物的寿命

17、随机服务系统中的服务时间都服从指数分布。机服务系统中的服务时间都服从指数分布。指数分布在排队论、可靠性分析中有广泛应用。指数分布在排队论、可靠性分析中有广泛应用。 )n,m,(exprnd： 注意：注意：Matlab中，产生中，产生参数参数为为 的指数分布的的指数分布的 命令为命令为exprnd( ).exprnd( ). 1例例 顾客到达某商店的间隔时间服从参数为顾客到达某商店的间隔时间服从参数为0.10.1的指的指数分布数分布 指数分布的均值为指数分布的均值为1/0.1=101/0.1=10。 指两个顾客到达商店的平均间隔时间是指两个顾客到达商店的平均间隔时间是1010个个单位时间单位时

18、间. .即平均即平均1010个单位时间到达个单位时间到达1 1个顾客个顾客. . 顾顾客到达的间隔时间可用客到达的间隔时间可用exprnd(10)exprnd(10)模拟。模拟。设离散型随机变量设离散型随机变量X X的所有可能取值为的所有可能取值为0,1,2,0,1,2,且取各个值且取各个值的概率为的概率为其中其中 0 0为常数，则称为常数，则称X服从参数为服从参数为 的的泊泊松分布。松分布。, 2 , 1 , 0k,!ke)kX(Pk 容器内的细菌数，十字路口的交通事故，寻呼台的寻呼容器内的细菌数，十字路口的交通事故，寻呼台的寻呼次数及每天到商店购买商品的顾客数等都服从泊松分布。次数及每天到

19、商店购买商品的顾客数等都服从泊松分布。泊松分布的期望值为泊松分布的期望值为 ),( poissrnd：nm 如相继两个事件出现的间隔时间服从参数为如相继两个事件出现的间隔时间服从参数为 的指数分布，的指数分布，则在单位时间间隔内事件出现的次数服从参数为则在单位时间间隔内事件出现的次数服从参数为 的泊松分的泊松分布即单位时间内该事件出现布即单位时间内该事件出现k k次的概率为：次的概率为：, 2 , 1 , 0k,!ke)kX(Pk 反之亦然。反之亦然。指数分布与泊松分布的关系：指数分布与泊松分布的关系： (1)(1)指两个顾客到达商店的平均间隔时间是指两个顾客到达商店的平均间隔时间是1010个

20、单位时个单位时间间. .即平均即平均1010个单位时间到达个单位时间到达1 1个顾客个顾客. . (2) (2)指一个单位时间内平均到达指一个单位时间内平均到达0.10.1个顾客个顾客例例 (1)(1)顾客到达某商店的间隔时间服从参数为顾客到达某商店的间隔时间服从参数为0.10.1的指数分的指数分布布 (2)(2)该商店在单位时间内到达的顾客数服从参数为该商店在单位时间内到达的顾客数服从参数为0.10.1的泊松分布的泊松分布 例例2 2敌坦克分队对我方阵地实施突袭，其到达规敌坦克分队对我方阵地实施突袭，其到达规律服从泊松分布，平均每分钟到达辆（律服从泊松分布，平均每分钟到达辆（1 1）模）模拟

21、敌坦克在分钟内到达目标区的数量，以及在拟敌坦克在分钟内到达目标区的数量，以及在第、分钟内各到达几辆坦克（第、分钟内各到达几辆坦克（2 2）模拟）模拟在在3 3分钟内每辆敌坦克的到达时刻。分钟内每辆敌坦克的到达时刻。 （1）用）用poissrnd(4)进行模拟进行模拟 To Matlab（poiss）（2）坦克到达的间隔时间应服从参数为）坦克到达的间隔时间应服从参数为4的的 指数分布，用指数分布，用exprnd（1/4）模拟。）模拟。 To Matlab（time）2 2随机变量的模拟随机变量的模拟 利用均匀分布的随机数可以产生具有任意分布利用均匀分布的随机数可以产生具有任意分布的随机变量的样本

22、从而可以对随机变量的取值情的随机变量的样本，从而可以对随机变量的取值情况进行模拟况进行模拟(1)离散型随机变量的模拟离散型随机变量的模拟 设随机变量设随机变量X的分布律为的分布律为P (X=xi)=pi,i=1,2,令令p(0)=0, p(n)= pi, , n=1,2,将将p(n)作为分点，把区间作为分点，把区间(0 , 1)分为一系列小区间分为一系列小区间(p(n-1), p(n)对于均匀的随机变量对于均匀的随机变量RU(0,l)，则有，则有P( p(n-1)R p(n)= p(n)-p(n-1)= pn ,n=1,2, 由此可知，事件由此可知，事件(p(n-1)R p(n)和事件和事件

23、X=xn) 有相同的发生的概率因此我们可以用随机变量有相同的发生的概率因此我们可以用随机变量R落在小区间内的情况来模拟离散的随机变量落在小区间内的情况来模拟离散的随机变量X的的取值情况取值情况 具体执行的过程是：每产生一个具体执行的过程是：每产生一个(0,1)上均匀分上均匀分布的随机数布的随机数r，若，若p(n-1)r p(n)则理解为发生事件则理解为发生事件“X=xn”例例3 3 随机变量随机变量 x = 0，1，2表示每分钟到达超市表示每分钟到达超市收款台的人数，有分布列收款台的人数，有分布列 xk 0 1 2 pk 0.4 0.3 0.3模拟十分钟内顾客到达收款台的状况模拟十分钟内顾客

24、到达收款台的状况. . To Matlab（paidui） 仿真分为静态仿真和动态仿真。动态仿真可分仿真分为静态仿真和动态仿真。动态仿真可分为连续系统仿真和离散系统仿真。离散系统是指状为连续系统仿真和离散系统仿真。离散系统是指状态变量只在某个离散时间点集合上发生变化的系统。态变量只在某个离散时间点集合上发生变化的系统。例如：电梯系统服务，排队系统，通信网络服务的例如：电梯系统服务，排队系统，通信网络服务的仿真等。连续系统是指状态变量随时间连续改变的仿真等。连续系统是指状态变量随时间连续改变的系统。例如：传染病的检测和预报等。系统。例如：传染病的检测和预报等。 为了仿真系统必须设置一个仿真时钟将

25、时间从为了仿真系统必须设置一个仿真时钟将时间从一个时刻向另一个时刻推进，并且可随时反映系统一个时刻向另一个时刻推进，并且可随时反映系统时间的当前值。模拟时间推进方式有两种：时间步时间的当前值。模拟时间推进方式有两种：时间步长法和事件步长法。模拟离散系统常用事件步长法。长法和事件步长法。模拟离散系统常用事件步长法。连续系统常用时间步长法（也称固定增量推进法或连续系统常用时间步长法（也称固定增量推进法或步进式推进）。步进式推进）。三、时间步长法三、时间步长法时间步长法时间步长法: : 按照时间流逝的顺序，一步一步地对按照时间流逝的顺序，一步一步地对系统的活动进行仿真。在整个仿真过程中，时间步系统的

26、活动进行仿真。在整个仿真过程中，时间步长固定不变。长固定不变。基本思想基本思想: : 在进行系统仿真的过程中，可以把整个过程在进行系统仿真的过程中，可以把整个过程分成许多相等的时间间隔，时间步长的长度可以分成许多相等的时间间隔，时间步长的长度可以根据实际问题分别取作秒根据实际问题分别取作秒, ,分分, ,小时，天等。程序小时，天等。程序中按照这个步长前进的时钟就是仿真的时钟。中按照这个步长前进的时钟就是仿真的时钟。基本步骤基本步骤: : 首先选取对象系统的一个初始起点作为仿真首先选取对象系统的一个初始起点作为仿真时钟的零点，然后根据实际问题的需要，选定一时钟的零点，然后根据实际问题的需要，选定

27、一个时间步长。于是从仿真时钟的零点开始，每推个时间步长。于是从仿真时钟的零点开始，每推进一个时间步长就对系统的活动和状态按照预定进一个时间步长就对系统的活动和状态按照预定的规则和目的进行考察、分析、计算和记录，直的规则和目的进行考察、分析、计算和记录，直到预定仿真结束时刻为止。到预定仿真结束时刻为止。时时间间步步长长法法初始状态时间步长加1在当前步长内，考察分析，计算和记录系统的活动仿真时间到否？结束输出结果是否应用举例应用举例-池水含盐量问题池水含盐量问题例例6:6:池水含盐量问题池水含盐量问题 某水池有某水池有2000m2000m3 3水，其中含盐水，其中含盐2kg2kg，以每分钟，以每分

28、钟6m6m3 3的速率向水池内注入含盐率为的速率向水池内注入含盐率为0.5kg/m0.5kg/m3 3的盐水的盐水，同时又以每分钟，同时又以每分钟4m4m3 3的速率从水池流出搅拌均匀的速率从水池流出搅拌均匀的盐水试用计算机仿真该水池内盐水的变化过的盐水试用计算机仿真该水池内盐水的变化过程，并每隔程，并每隔10min10min计算水池中水的体积、含盐量和计算水池中水的体积、含盐量和含盐率欲使池中盐水的含盐率达到含盐率欲使池中盐水的含盐率达到0.2kg/m0.2kg/m3 3，需，需经过多少时间？经过多少时间？池水含盐量随时间变化的示意图池水含盐量随时间变化的示意图 分析分析: :系统中，实体是

29、水，属性是水的体积、含盐系统中，实体是水，属性是水的体积、含盐量和含盐率，活动是水的注入与流出，由于注入和流量和含盐率，活动是水的注入与流出，由于注入和流出活动的作用，池中水的体积与含盐量、含盐率均随出活动的作用，池中水的体积与含盐量、含盐率均随时间变化，初始时刻含盐率为时间变化，初始时刻含盐率为0.001kg/m0.001kg/m3 3, ,以后每分钟以后每分钟注入含盐率为注入含盐率为0.5 kg/m0.5 kg/m3 3的水的水6m6m3 3，流出混合均匀的盐水，流出混合均匀的盐水4m4m3 3, ,当池中水的含盐率达到当池中水的含盐率达到0.2kg/m0.2kg/m3 3时时, ,结束仿

30、真过程结束仿真过程. . 为了能定量地考察系统在任一时刻的属性，引入为了能定量地考察系统在任一时刻的属性，引入下列记号：下列记号： 注水速度：注水速度： WIWI6 m6 m3 3/min/min 排水速度：排水速度： WOWO4 m4 m3 3/min/min 注入水的含盐率：注入水的含盐率：SISI0.5 kg/m0.5 kg/m3 3 最终含盐率：最终含盐率：SFSF0.2 kg/m0.2 kg/m3 3 T T时刻水的体积：时刻水的体积：VTmVTm3 3 T T时刻水的含盐量：时刻水的含盐量：ST kgST kg T T时刻水的含盐率：时刻水的含盐率：SR=ST/VT kg/mSR=

31、ST/VT kg/m3 3对于这样一个连续系统仿真时，必须在一系列离散时对于这样一个连续系统仿真时，必须在一系列离散时间点间点t t0 0 t t1 1 t t2 2 t tn n上来进行考察，这些离散时间上来进行考察，这些离散时间点之间的间隔点之间的间隔T= ti- ti-1(i=1,2,n)就是时间步长若取就是时间步长若取步长为步长为1min时，就要每隔时，就要每隔1min考察一次系统的状态特考察一次系统的状态特性，并作出相应的记录与分析在本题中每隔性，并作出相应的记录与分析在本题中每隔1min池池水的动态变化过程是：水的动态变化过程是： 每分钟水的体积增加每分钟水的体积增加 6-42 m

32、3； 每分钟向池内注入盐每分钟向池内注入盐 60.5 3 kg； 每分钟向池外流出盐每分钟向池外流出盐 4SR kg； 每分钟池内增加盐每分钟池内增加盐 3-4SR kg取取T0作为系统仿真的初始时刻，取作为系统仿真的初始时刻，取T=0时系统的属时系统的属性为系统的初始状态取池中盐水的含盐率达到性为系统的初始状态取池中盐水的含盐率达到0.2kg/m3的时刻作为仿真的终止时刻的时刻作为仿真的终止时刻池池水水含含盐盐量量仿仿真真流流程程初始化按时间步长前进1min计算池水体积VT,含盐量ST,含盐率SR是否到达10min？结束输出结果是否含盐率达到SF?打印次数加1输出时间T,水体积VT含盐量ST

33、含盐率SR下一个10min开始，计时单元清零是否池池水水含含盐盐量量仿仿真真结结果果 在这个问题中，池水的体积和含盐量、含盐率都在这个问题中，池水的体积和含盐量、含盐率都是随时间连续变化的故这是一个连续系统的仿真问是随时间连续变化的故这是一个连续系统的仿真问题，我们是把这个连续过程离散化后用时间步长法进题，我们是把这个连续过程离散化后用时间步长法进行仿真的我们还可以用微分方程建立这个问题的数行仿真的我们还可以用微分方程建立这个问题的数学模型，并求出它的解析解，这个解析解就是问题的学模型，并求出它的解析解，这个解析解就是问题的精确解有兴趣可以按照这一思路求出该问题的精确精确解有兴趣可以按照这一

34、思路求出该问题的精确解，考察相应时刻精确解与仿真解的差异也可进一解，考察相应时刻精确解与仿真解的差异也可进一步调整前面仿真过程的时间步长，通过与精确解的比步调整前面仿真过程的时间步长，通过与精确解的比较来研究时间步长的大小对仿真精度的影响较来研究时间步长的大小对仿真精度的影响四、事件步长法四、事件步长法事件步长法事件步长法: 是以事件发生的时间为增量，按照事是以事件发生的时间为增量，按照事件发生的时间顺序，一步一步地对系统的行为进行件发生的时间顺序，一步一步地对系统的行为进行仿真，直到预定的时间结束为止。仿真，直到预定的时间结束为止。事事件件步步长长法法初始状态事件步长加1在当前步长内，考察分

35、析，计算和记录系统的活动仿真时间到否？仿真结束输出结果是否事件步长法与时间步长法的主要事件步长法与时间步长法的主要区别区别是：是：(1)事件步长法与时间步长法都是以时间为增量来考事件步长法与时间步长法都是以时间为增量来考察系统状态的变化但在时间步长法中，仿真时钟察系统状态的变化但在时间步长法中，仿真时钟以等步长前进，而在事件步长法中，仿真时钟的步以等步长前进，而在事件步长法中，仿真时钟的步长取决于事件之间的间隔长取决于事件之间的间隔(2)时间步长法在一个步长内，认为系统所处的状态时间步长法在一个步长内，认为系统所处的状态相同，因而所选步长的大小将影响仿真的精度而相同，因而所选步长的大小将影响仿

36、真的精度而在事件步长法中，每个事件的发生均有确切的时刻在事件步长法中，每个事件的发生均有确切的时刻，不需要人为地选取步长，步长的大小对仿真精度，不需要人为地选取步长，步长的大小对仿真精度影响较小影响较小(3) 时间步长法每前进一个步长就要对整个系统进时间步长法每前进一个步长就要对整个系统进行一次全面考察，即使状态没有发生变化时也要行一次全面考察，即使状态没有发生变化时也要扫描，而事件步长法只是在某一事件点上判断和扫描，而事件步长法只是在某一事件点上判断和比较事件是否出现比较事件是否出现 因此，一般地讲，当判断比较的数目较大时，因此，一般地讲，当判断比较的数目较大时，用时间步长法可以节省用机时间

37、而当相继两个用时间步长法可以节省用机时间，而当相继两个事件出现的平均间隔较长时，更适合于用事件步事件出现的平均间隔较长时，更适合于用事件步长法长法例例7 收款台前的排队过程的仿真收款台前的排队过程的仿真 假设：假设： 1、顾客的到达、顾客的到达收款台收款台是随机的，平均时间间隔是随机的，平均时间间隔 为为0.5分钟，即间隔时间服从分钟，即间隔时间服从 =2的指数分布。的指数分布。 2、对不同的顾客收款和装袋的时间服从正态分布、对不同的顾客收款和装袋的时间服从正态分布 N(1，1/9)。 现要求模拟现要求模拟20位顾客到收款台前的排队情况，位顾客到收款台前的排队情况，主要关心问题：每个顾客的平

38、均等待时间主要关心问题：每个顾客的平均等待时间,队长队长,服服务务员的工作效率员的工作效率.收款排队系统主控程序图收款排队系统主控程序图顾客到达子程序图顾客到达子程序图服务结束子程序图服务结束子程序图假设假设 t(i): 第第i位顾客到达时刻，位顾客到达时刻， t2(i):第第i位顾客受到的服务时间位顾客受到的服务时间(随机变量随机变量)，T(i): 第第i位顾客离去时刻位顾客离去时刻.将第将第i位顾客到达作为第位顾客到达作为第i件事发生件事发生; t(i+1)- t(i)= r(i) (随机变量随机变量)平衡关系：平衡关系：当当 t(i+1) T(i) 时时, T(i+1)=t(i+1)+t

39、2(i+1);否则否则, T(i+1)=T(i)+t2(i+1)Matlab(kefu)程序程序五、五、Monte Carlo方法方法基本思想：基本思想： Monte Carlo方法亦称为随机模拟方法亦称为随机模拟（Random simulation）方法，有时也称作随机抽样技术或统计方法，有时也称作随机抽样技术或统计试验方法。试验方法。此方法对研究的系统进行随机观察抽样，此方法对研究的系统进行随机观察抽样，通过对样本值的观察统计，求得所研究系统的某些通过对样本值的观察统计，求得所研究系统的某些参数参数 当试验次数充分多时，某一事件出现的频当试验次数充分多时，某一事件出现的频率近似于该事件发生

40、的概率率近似于该事件发生的概率. .例例8 8：用蒙卡模拟求圆周率用蒙卡模拟求圆周率的估计值。的估计值。11 设二维随机变量设二维随机变量(X,Y)在正方形服在正方形服从均匀分布，从均匀分布， 则则(X,Y)落在四分之一圆内的概率落在四分之一圆内的概率为为 PX2+Y2=1= /4 现产生现产生n对二维随机点对二维随机点(xi,yi), xi,yi是区间是区间0,1上上均匀分布的随机数，若其中有均匀分布的随机数，若其中有k对满足对满足 xi2+yi2=1，即相当于做即相当于做n次投点试验，其中有次投点试验，其中有k点落在点落在1/4圆内。圆内。计算落入计算落入1/4圆内的频率为圆内的频率为k/

41、n，由大数定律，事件，由大数定律，事件发生的频率依概率收敛于发生的概率，可得圆周率发生的频率依概率收敛于发生的概率，可得圆周率的估计值为的估计值为4k/n。次数越多，精度越大。次数越多，精度越大。Matlab(yuanzhoulv)例例9 9：用用Monte Carlo方法计算积分方法计算积分 badxxgJ)(分析：任取一列相互独立的、都具有分析：任取一列相互独立的、都具有a,b中均匀中均匀分布的随机变量分布的随机变量xi,则则g(xi)也是一列相互独立的也是一列相互独立的随机变量随机变量, ,而且：而且：abJbaabidxxgxEg)()(1所以所以)()(ixEgabJ只要求出只要求出

42、)(ixEg, ,便能得到便能得到J 的数值。的数值。为求为求Eg(xi), )()()()(21inxEgnxgxgxg 由大数定理知：由大数定理知：) ( n当当应用举例应用举例-报童的策略报童的策略例例10 10 报童每天清晨从报社购进报纸零售报童每天清晨从报社购进报纸零售, ,晚上将晚上将没有卖掉的报纸退回没有卖掉的报纸退回. .每份报纸的购进价为每份报纸的购进价为1.31.3元元, ,零售价为零售价为2 2元元, ,退回价为退回价为0.20.2元元. .报童售出一份报纸报童售出一份报纸赚赚0.70.7元元 , ,退回一份报纸赔退回一份报纸赔1.11.1元元. .报童每天如果购报童每天

43、如果购进的报纸太少进的报纸太少, ,不够卖时会少赚钱不够卖时会少赚钱, ,如果购得太多如果购得太多卖不完时要赔钱卖不完时要赔钱. . 试为报童筹划每天应如何确定试为报童筹划每天应如何确定购进的报纸数使得收益最大购进的报纸数使得收益最大. .报纸每捆报纸每捆1010张，只能张，只能整捆购买，报纸可以分为整捆购买，报纸可以分为3 3种类型的新闻日：好、种类型的新闻日：好、一般、差，它们的概率分别为一般、差，它们的概率分别为0.35,0.450.35,0.45和和0.2,0.2,在在这些新闻日中每天对报纸的需求分布的统计结果这些新闻日中每天对报纸的需求分布的统计结果如下表：如下表：试确定每天报童应该

44、订购的报纸数量。试确定每天报童应该订购的报纸数量。解：解：我们通过计算机仿真来解决此问题。最优策略我们通过计算机仿真来解决此问题。最优策略应该是每天的利润最大。应该是每天的利润最大。 利润利润= =销售收入销售收入- -报纸成本报纸成本+ +残值残值这是一个随机现象的计算机仿真问题，这是一个随机现象的计算机仿真问题，故先确定各种情况的随机数的对应关系。故先确定各种情况的随机数的对应关系。新闻日和需求量对应的随机数分别如下面两个表格所示新闻日和需求量对应的随机数分别如下面两个表格所示计算机仿真的流程：计算机仿真的流程：1 1）令每天的报纸订购数变化，）令每天的报纸订购数变化，40-10040-1

45、00；2 2）让时间从）让时间从1 1开始变化（循环）到开始变化（循环）到36365 5；3 3）产生新闻种类的随机数，确定当天的新闻类型；）产生新闻种类的随机数，确定当天的新闻类型；4 4）产生需求量随机数，确定当天的报纸需求量；）产生需求量随机数，确定当天的报纸需求量；5 5）计算当天的收入，计算累积利润，）计算当天的收入，计算累积利润，6 6）比较得出最优定货量。）比较得出最优定货量。需求量需求量好新闻的随机数区间好新闻的随机数区间一般新闻的随机数区间一般新闻的随机数区间差新闻的随机数区间差新闻的随机数区间 40 （0.00,0.03 （0.00,0.10 （0.00,0.44 50 0

46、03,0.08) 0.10,0.28) 0.44,0.66) 60 0.08,0.23) 0.28,0.68) 0.80,1.00) 70 0.23,0.43) 0.68,0.88) 0.80,1.00) 80 0.43,0.78) 0.88,0.96) 0.80,1.00) 90 0.78,0.93) 0.96,1.00) 100 0.93,1.00 Matlab编程编程(baotong程序程序)x1=rand(365,1);x2=rand(365,1);for n=4:10 paper=n*10; sb(n)=0; for i=1:365 if x1(i)0.35 if x2(i)0.0

47、3 news=40; elseif x2(i)0.08 news=50; elseif x2(i)0.23 news=60; elseif x2(i)0.43 news=70; elseif x2(i)0.78 news=80; elseif x2(i)0.93 news=90; else news=100; endelseif x1(i)0.8 if x2(i)0.10 news=40; elseif x2(i)0.28 news=50; elseif x2(i)0.68 news=60; elseif x2(i)=news sale=news; remand=paper-news; els

48、e sale=paper;remand=0; end sb(n)=sb(n)+2*sale- 1.3*paper+0.2*remand; endendelseif x2(i)0.96 news=80; else news=90; end else if x2(i)0.44 news=40; elseif x2(i)0.66 news=50; elseif x2(i)0.82 news=60 elseif x2(i)=optmoney optnews=n*10; optmoney=sb(n); end n,sb(n)/365endoptnews,optmoney,optmoney/365 经过计

49、算机仿真后得到最优购货量是每天经过计算机仿真后得到最优购货量是每天6060份，份，平均每天利润平均每天利润34.434.4元。元。1、(射击命中率射击命中率) 在我方某前沿防守地域，敌人以在我方某前沿防守地域，敌人以一个炮排（含两门火炮）为单位对我方进行干扰和一个炮排（含两门火炮）为单位对我方进行干扰和破坏为躲避我方打击，敌方对其阵地进行了伪装破坏为躲避我方打击，敌方对其阵地进行了伪装并经常变换射击地点经过长期观察发现，我方指并经常变换射击地点经过长期观察发现，我方指挥所对敌方目标的指示有挥所对敌方目标的指示有50是准确的，而我方火是准确的，而我方火力单位，在指示正确时，有力单位，在指示正确时

50、有1/3的射击效果能毁伤敌的射击效果能毁伤敌人一门火炮，有人一门火炮，有1/6的射击效果能全部消灭敌人现的射击效果能全部消灭敌人现模拟我方将要对敌人实施的模拟我方将要对敌人实施的20次打击结果，并确定次打击结果，并确定有效射击的比率及毁伤敌方火炮的平均值。有效射击的比率及毁伤敌方火炮的平均值。实验作业实验作业2、已知零件、已知零件C由零件由零件A和零件和零件B连接而成，已知连接而成，已知A、B的长度均为随机变量，具体数值如下。试抽取的长度均为随机变量，具体数值如下。试抽取10个样本以计算个样本以计算C的平均长度。的平均长度。3、用、用Monte Carlo方法计算下列积分。方法计算下列积分。