系统工程-第四章

1、第四章 系统仿真w 4.1 基本概念w 4.2 蒙特卡罗仿真w 4.3 均匀随机数的产生方法w 4.4 随机变量的产生方法w 4.5 离散事件随机系统的仿真方法w 4.6 4.6 随机分布的确定随机分布的确定w 4.7 系统动力学仿真w 4.8 4.8 仿真语言仿真语言w 4.9 4.9 可视化仿真技术可视化仿真技术w 4.10 4.10 虚拟现实技术虚拟现实技术w 4.11 4.11 仿真技术在矿业中的应用仿真技术在矿业中的应用w 作业4.1 4.1 基本概念基本概念 在下图中，已知采区来煤量q1是一个随机变量，根据统计结果，它服从正态分布，已知其均值=1.45t/ min，其标准差=0.0

2、95t/min。大巷是胶带输送机，其运输量是常数：q2=1.33t/min。求8h后，采区最大煤仓容量。 q2q1煤仓投掷硬币投掷硬币随机模拟实验随机模拟实验 建立一个计算机能够运算的模拟模型,问题的关键是用随机数去模拟投掷硬币的实验。 所谓随机数就是一组随机出现的数列，通常用的是0，1区间的随机数。这时要求0，1区间的随机数必须具备下列性质：一是它的均匀性，即这些随机数落在0，1区间内任一位置的概率是相同的；另一个是它的随机性，即这些数在 0，1区间内任一位置出现与否是随机的、独立的。如以R表示随机数，则出现0R0.5和 0.5R1的概率相等，都为0.5.这样如以随机数R0.5模拟出现正面，

3、以R0.5模拟出现反面，则模拟模型就建立起来了。 设设 I I 是投掷次数是投掷次数 N N 为正面出现次数为正面出现次数 R R 为随机数为随机数 P P是正面出现频率是正面出现频率 I0.50.5实验开始 I=I+1产生随机数R记录出现正面次数N=N+1计算并输出结果N，P模拟结束初始状态I=0，N=0YNYN图3-2 投掷硬币实验模拟框图返回4.2 4.2 蒙特卡罗仿真蒙特卡罗仿真 设f(x)是区间（a,b）内的连续函数，且已知在该区间内f(x)=0,要求计算在（a，b）和（0，c）区间中分别产生两个均匀分布的随机数xi和yi，由此构成的坐标点 Pi(xi，yi)必定在矩形abde内。若

4、yif(xi)，则接受该点并作计数统计；否则舍弃该点不作统计. badxxf)(拒绝接受dxxbxjxia0cef(xi)f(xj) 设共产生N个随机仿真点，如果共接受个m点，（即恰在曲线上或在曲线以下的点），当N充分大时有 曲线下面积 badxxfcabNm)()(返回4.3 4.3 均匀随机数的产生方法均匀随机数的产生方法 4.3.1 4.3.1 均匀随机数均匀随机数 0，1区间上的均匀随机数，它的密度函数f(x)和分布F(x)分别为 其他0101)(xxf111000)(xxxxxF 4.3.2 4.3.2 均匀随机数的产生方法均匀随机数的产生方法 用数学公式产生的均匀随机数，称为伪随机

5、数,常用的公式方法是线性同余法 . 若有a、b两个整数被某个整数M除后，它们的余数相同，则称a、b两数对模M同余，记作 a=b(ModM) 求余数的公式为例如, MMyyc1510505151010515515c产生随机数的线性同余法是一个递推公式： 即 通过下式得到落在0，1区间内的随机数： 两个较好的计算公式如下 niMcaZZii,2, 1)(mod(1MMcaZcaZZiii11Mxrn1)2)(mod15(35115iiZZ)2)(mod453806246314159269(311iiZZ4.3.3 4.3.3 均匀随机数检验方法均匀随机数检验方法 设N个随机数从小到大排列起来并分成

6、互不重叠的等区间的k段.从理论上讲，这些随机数应该均匀分布在每一段上，在每一段上理论频数为如果统计出随机数实际落在每段的频数为ni，则令 x2值越大表明随机数列越不均匀。 kikNmi,2,1,kiiiimmnx122)( 给定判定标准称为显著度,根据自由度r=k-1，从x2表中查得 值，如果计算所得的x250），若取零假设H0：=0，并取统计量 (近似服从N（0，1）标准正态分布)，同时选定显著度=0.05，则根据概率论，当|u|1.96时为差异显著，拒绝H0假设；当|u|1.96时接受假设，即检验的随机数符合独立性要求。4.3.4 4.3.4 计算机函数计算机函数 计算机高级语言有产生随机

7、数的标准函数，可直接调用。 202/)(1SZZZNjNijiijjNuj返回4.4 4.4 随机变量的产生方法随机变量的产生方法 4.4.1 4.4.1 逆变换法逆变换法 已知某随机变量的分布函数F（x） ，并知其定义域为0，1区间，这时可用0，1区间的均匀分布随机数R来表示F（x），即F（x）=R，由此可得到随机变量的样本值x=F-1(R)。F-1为逆变换的表达式。xx0R1F(x)例如，已知负指数分布的密度函数为：例如，已知负指数分布的密度函数为：其分布函数为 令上式两边取对数后有 000)(xxexfx01)(0 xedtexFxxtuReexFRxx11)(uxuxlnln4.4.2

8、 4.4.2 近似法近似法 正态分布的密度函数为 分布函数为 当=0，=1时称标准正态分布: xeaxfx, 021)(222)(xtdtexF222)(21)(2221)(xexfxtdtexF2221)( 任何非标准的正态分布随机变量Y，通过下式都可以化为标准正态分布的随机变量李雅普诺夫的中心极限定理 设 是相互独立的随机变量，它们的数学期望和方差为令 当满足 yxn21,nkDEkkkk2 , 1,)(,)(2nknkkknkk11221,0)(lim212nkkkEE则随机变量 的分布近似服从标准正态分布 如为 为0，1区间均匀随机数，其期望 ， 方差 xnkknknkkkx1211n

9、21,2/ 12/ )()(abEk)1,0(12112)()(2baabD有 若已知一非标准正态，其期望值为 ，均方差为 ，按前式，则非标准的正态分布随机变量为将标准正态分布的随机变量代入后为 当n=5时 )2(121221212111111nkknkknknknkknnnnxxy)2(121nkknny) 3(261kky返回4.5 4.5 离散事件随机系统的仿真方法离散事件随机系统的仿真方法4.5.1 4.5.1 仿真时钟仿真时钟 仿真时钟是随着仿真的进程而不断更新的时间机构。通常，在仿真开始时将仿真时钟置为零。随后，仿真时钟不断地给出仿真时间的当前值。仿真时间是仿真模型中的时间指示，代

10、表模型运行的真实时间。4.5.2 4.5.2 时间步长仿真法时间步长仿真法 一个系统的各种变量和参数往往随着时间变化而变化，因此可把时间分成段（如1天、1小时、1分钟或1秒钟），仿真时钟按时间段等距推进，每次推进需要扫描所有的活动，模拟各段时间内系统各种变量和参数的变化，这种仿真方法称为时间步长法。 YN初始条件和数据增加一个时间步长在当前的时间步长内考察、分析、计算和记录系统的活动。时钟是否到达模 拟总时间?对模拟数据进行分析，输出结果。图3-5 时间步长仿真法 以井底煤仓系统为例构造仿真模型，分析煤仓容量变化情况. 煤煤仓箕斗 主井煤仓系统示意图 根据实测数据统计分析，单位时间内（根据实测

11、数据统计分析，单位时间内（10min10min） 煤仓来煤量为符合正态分布的随机变量，其均值为1=125t/10min，均方差为1=12.5t/10min；同时单位时间内箕斗提升次数也是随机的，统计结果见表3-1，符合正态分布，每箕斗的提升量为10t。根据表3-1中的数据可求出其均值为2=46.5t/10min，均方差为2=9.63t/10min。现在要求模拟300天，求出： 每天最大煤仓容量和300天的平均值； 300天的平均日产量。 设10min为仿真时钟的时间推进步长，每隔10min考察一次来煤量、煤仓容量提升量，构造仿真模型，其结构见图3-7。 产生来煤量随机数R1，产生箕斗提升随机数

12、R2， 模拟天数D=0煤仓容量V=0，最大煤仓容量VD=0，日产量QD=0，仿真时钟T=0。时钟推进：T=T+1令 V=V+R1-10R2，Rd= Rd+ R1。判断VD是否小于V，若小于V则令VD=V D=D+1T300d ?输出煤仓容量每天最大值和300d平均值以及日产量平均值。NNY图3-7 主井煤仓仿真结构4.5.3 4.5.3 事件步长仿真法事件步长仿真法 事件步长法是以事件的发生来推进模拟时钟的增长。 初始状态 选取最短时间事件 考察此事件所处的状态 改变事件的状态，并计算有关的数据。 累加事件的当前时间是否模拟完毕 ？ 输出结果NY 图3-8 事件步长仿真法汽车运输系统汽车运输系

13、统 已知一汽车队共有N辆汽车，每辆车的载重量为Qt，每天工作7h。图3-9为运输系统示意图。车库和仓库在一起，工作开始时汽车在仓库前排队装车，每次一辆；装完后汽车重载运行到火车站卸货；卸车后空载运行回到仓库再装车。经统计知道，汽车运输的装车、重载运行、卸车、空载运行时间都是随机的，符合正态分布，其均值和均方差见表3-2。因汽车运行一个周期要花一个多小时，因此汽车到达仓库时间等于或超过6点半时就不再发车。 现在需要通过建立仿真模型，了解以下情况： 在汽车数N一定的情况下，每班可运输多少吨货物（QA）？ 汽车排队总等待时间（F）； 汽车的利用率（FS）； 仓库的空闲时间（A）及仓库装车的利用率（A

14、S）； 每辆车平均运行时间（Ti） 仓库 3 N 2 1 汽车空载运行重载运行火车站卸车图3-9 汽车运输系统 设置初始状态计算第一轮装车及有关参数求最短时间事件（即最先回到仓库的汽车当前时间及号码）求装车、卸车和空、重载运行的正态分布运行时间。汽车到达时间仓库当前时间？计算汽车等待时间及有关参数计算卸车吨位数是否到收车时间？计算收车数目及运行时间汽车是否全部回到车库？计算有关参数及输出每天数据模拟天数是否完毕？输出分析结果计算仓库空闲时间及有关参数YNNNYYNY4.5.4 4.5.4 仿真结果的分析仿真结果的分析 仿真结果的平均值 4.5.5 误差分析误差分析 设用 作为未知参数的估计量，

15、两者误差不超过某一正数的概率为1-，即 区间 称为置信区间 , 1-称为置信概率，称为显著性水平（或置信度） nxxxxn211)(p),(由中心极限定理可知，当试验次数足够大时服从正态分布，由于2未知，需用样本方差S2来表示，它服从t分布，即：对给定置信度，有nxU/) 1()(ntSnx1)1()(2/ntSnxp即x在置信度下的置信区间是例如，5次仿真结果：12500，12650，12600，12750。试求置信概率为0.99的置信区间。 ) 1(,) 1(2/2/nSntxnSntx125905151iixx45700)(11122niixxnS试验次数n=5，n-1=41-=0.99

16、, /2=0.005 因此有：置信区间为（12344.19，12835.81） 6041. 4) 4() 1(005. 02/tnt81.128355457006041. 412590) 1(19.123445457006041. 412590) 1(2/2/nSntxnSntx 定置信区间（即允许波动范围）及置信度，欲求重复仿真次数n值，即寻找最小的n，满足 可利用试凑法，即令n=0，n=1，n=2，直到满足上式为止。 例如，已知=0.1，=1，S=1.78，分别将n=8、9、10、11、12代入上式，可发现n=11为符合要求的最少仿真次数。 nSnt) 1(2/返回4.6 随机分布的确定随

17、机分布的确定4.6.1 分布函数的假设分布函数的假设样本直方图 有一组样本值x1，x2，xn。样本值中的最大值xmax和最小值xmin。选b略小于xmax，选a略大于xmin，将区间a，b等分为m个互不相交的子区间，即有下列分割点ti： a = t0 t1 t2 ti tm = b其中 mabttii1 样本值落在区间ti-1，ti）中的个数，记作Vi ,则样本落入本区间的频率 以组距（ti-1，ti）为底边，以 为高度绘制直方矩形,构成直方图.画一条曲线大致经过各个直方矩形的上边，便得出近似的概率密度的函数。 mVfii1iiiittfyx0yti-1 ti4.6.2 参数估计参数估计 设总

18、体X的概率密度的形式f（x，）为已知，它只含一个未知参数。若x1，x2，xn为X的一组观察值，则它们的联合概率密度等于记作若存在的一个值 ，使得似然函数在时= 时则称 是的极大似然估计。 ),(),(121niinxfxxxf),(),(121niinxfxxxLLmax),(21nxxxL 求极大: 令或例如，设X服从指数分布 0Ldd0lnLdd000)(xxexfxniniinxnxexxxLi1121)exp(),()(ln)(lnln1xnxnLniiniixnx11令 得0lnLdxdniixnx114.6.3 X2检验检验 用n个观测值来检验总体分布的假设H0。 H0：总体X的分

19、布函数为F（x）。 将x轴分为k个不相交的区间（ai，ai+1），i =1，2，k,其中ai、ak+1可分别取-及+ ,这时有 pi = F（ai+1）- F（ai）= p（ai50）时，则不论总体属什么分布，皮尔逊统计量总是近似服从自由度为k-1的分布，其中k是所分区间个数，是估计参数的个数。 于是，若在假设下算得的 ，则在水平下拒绝H0，否则就接受H0。) 1(22rkxx返回4.10 系统动力学系统动力学4.10.1 基本概念及特点 系统动力学（System Dynamics，简称SD）是一种研究复杂系统的方法1-4。该方法基于系统论，并融会了控制论、信息论、组织理论等多种理论，通过对系

20、统结构因果关系的分析，设计出反映系统行为的反馈回路，建立相应的方程式，最终建立系统动态模型，以此为基础对系统进行仿真实验。 特点:w 系统动力学问题是动态的问题，这些问题通常是用随时间连续变化的量来表示的，如人口变化、资源变化等。w 应用系统动力学研究复杂系统，能够容纳大量变量，一般可达数千个以上。w 系统动力学模型，既有描述系统各要素之间因果关系的结构模型，以此来认识和把握系统结构，又有专门形式表现的数学模型，据此进行仿真实验和计算，以掌握系统的未来动态行为。因此，系统动力学是一种定性分析和定量分析相结合的仿真技术。 w 在系统动力学模型中，能够设定各种控制因素，当改变输入的控制因素（例如不

21、同的组织状态、经济参数或不同的政策因素）时，可观察系统的行为和变化的趋势，对系统进行动态模拟。w 系统动力学通过模型进行仿真计算的结果，都用预测未来一定时期内各种变量随时间而变化的曲线来表示，也就是说，系统动力学能处理高阶次、非线性、多重反馈的复杂时变的复杂系统的有关问题。4.10.2 回路w 复杂系统总是由许多因果反馈回路偶合而成的。系统动力学就是使用反馈来揭示原因和寻找解决问题的办法，因果关系是构成系统动力学模型的基础。所谓反馈回路是指由两个以上具有因果关系的变量,以因果关系彼此连接而成的闭合结构.为了建立系统动力学方程，通常用反馈决策回路表示系统各元素间的关系. 回路分正回路和负回路:

22、生产增加商品减少收入增加+-订货速度库存差额库存量-+-a-正反馈回路b-负反馈回路图 2-15 正、负反馈回路+-+人口密度平均寿命死亡率人口生育率出生率+-+- 图 2-16 复杂因果关系回路决策回路： w因果反馈回路表达了系统发生变化的原因，但这种定性描述还不能确定使回路中的变量发生变化的机制。为了建立系统动力学方程，可进一步用反馈决策回路表示系统各元素间的关系，其基本结构可用下图表示。 信息水平(系统的表现水平)水平(系统状态)决策源决策反馈回路基本结构4.10.3 系统动态行为 w系统动态行为是指系统特征变量随时间变化的情况，典型情况：w曲线A表示的是简单的反馈系统，变量随时间增加而

23、增加，但增加速率随变量值增加而逐渐减小，最后趋于某数值，属负反馈。曲线B则为复杂系统的动态行为，随时间而振荡，然后趋于某一目标值，属负反馈。曲线C表示系统变量随时间作指数增长，属正反馈。曲线D表示的动态行为开始一段时间如同曲线C，而后同曲线B呈振荡状，属组合型系统动态行为。 4012312345BADC 动态行为曲线4.10.4 模型方程模型方程 w 系统是由许多反馈回路组合而成的，而每一个反馈回路模型则是通过一组微分方程来表示的，这就是模型方程。w 模型方程能反映系统行为随时间改变的动态行为，其中最重要的两个概念是：水平变量（也称状态变量）和速率变量（也称决策变量）。所谓水平变量是指能表征系

24、统某种属性的变量，一般它应该是一个积累量，如矿产资源地质储量、库存量等；而速率变量是指状态变量变化的速度，在系统中描述的是物质的实际流动，如矿产资源的增长与消耗、库存的入库与出库等。下下图反映了系统动态行为过程: 水平方程。水平方程标识符为L，如 L LL.K=LL.J+（DT）（RA.JK-RS.JK） 速率方程。速率表示单位时间内水平的变化，标号为R。以订货速率为例： R OR.KL=（1/AT）（DI-IK） w 辅助方程。 辅助方程的标号为A。上式的速率方程中，DI要求的库存量假如为一变量（等于要求库存量平均销售的周数WID与平均销售率AS.RK的乘积），通过辅助方程独立列出：w R

25、OR.KL=（1/AT）（DI.K-IK）w A DI.K=（WID）（AS.RK） J DT K DT LL1.LL2.LL1.KL2.KL1.jL2.jR1.JKR2.JK时间 图2-19 系统动态行为过程描述返回系统动力学建模步骤 w 首先，确定系统的边界，画出因果图；w 选择模型的基本变量水平变量；w 以水平变量为中心构造各自的子系统；w 然后，根据因果图，连接各子系统；w 根据以上的描述，写出方程式；w 最后，进行仿真运算，进行仿真结果及政策分析。 矿区规划系统动力学模型 从矿区煤炭生产系统的整体性出发，研究矿区发展的规模、矿区产量、投资、效益等彼此相互影响和制约的关系。通过对矿区生

26、产和自然资源内部结构和变量之间的因果关系、社会行为和发展过程的仿真实验，对矿区未来的发展趋势进行多种可行方案的预测。 矿区规划的总模型 技术进步职工培训固定资产生产效率投资工业产值职工总人数矿区生产能力施工工期自投资矿区人员矿井衰老物资消耗利润留成矿区环境矿区总产量销售利润销售总收入+-共有5个正反馈回路和2个负反馈回路。 正回路： 投资-矿区生产能力-产量-收入利润-利润留成-自投资-投资。 投资-固定资产-产值在企业收入中的贡献-收入-利润-利润留成- 自投资-投资。 投资-职工培训-生产效率-职工数减少-物质消耗-利润-利润留成- 自投资-投资。 投资-技术进步-生产效率-工期缩短-物质

27、消耗-利润-利润留成- 自投资-投资。 投资-技术进步-矿区环境-收入-利润-利润留成-自投资-投资。负回路： 矿区生产能力-产量-矿井衰老与报废-矿区生产能力。 矿区生产能力-产量-矿区环境-收入-利润-利润留成-自投资-投资 -矿区生产能力。矿区总模型划分为矿区总模型划分为5个子模型个子模型:w矿区生产能力与产量子模型：从已有的矿井生产能力增加或减少、新建矿井的增加两个方面考虑。 w矿区经济子模型:矿区经济子模型以固定资产和洗精煤产量的变化作为变量，主要反映矿区收入、利润和积累的状况，综合体现矿区不同规划方案的经济效果。w矿区职工人数和效率子模型:矿区职工人数子模型通过工人数和效率反映矿区

28、生产发展、工资增长和需要进行生活福利设施建设等方面的情况。w矿区人口和住宅子模型:矿区人口和住宅子模型主要反映矿区人口发展状况及对住房的要求。w矿区环境子模型:矿区环境子模型主要考虑矿区规划中社会效益方面的内容，按照煤矿开采的特点，从矸石、废水、废气、地表塌陷破坏等方面对矿区的环境状况进行评价。矿区生产能力与产量子模型： DELAY3DELAY3量增加老井储开采量老井老井可采储量新井可采储量量增加新井储开采量新井计划产量矿区总产量老井报废需求差额改扩建投资新井投资矿区总生产能力老井能力增加老井能力减少改扩建开工率新井开工率新井能力老井能力图2-16 矿区生产能力与产量模型流程图系统模拟：w方案

29、：对8个生产矿井进行改扩建，并陆续新建3对矿井。w方案：对5个生产矿井进行改扩建，并陆续新建3对矿井。w方案：只对8个生产矿井进行改扩建，不建新矿井。w方案：维持现有状况，既不改扩建也不新建矿井。w利用系统动力学模型对矿区30年内的生产能力、产量、投资、经济效益、生产效率、职工人数、人口、生态环境等各项参数的动态发展和变化进行模拟计算，得到大量数据结果供决策层分析决策。 煤产量模拟结果比较煤产量模拟结果比较 方案矿区产量，万吨各方案产量比较1990199520002005201020151990199520002005201020151882.42069.52134.72345.72476.4

30、2256.90000001822.41936.01925.22119.42257.22048.30-133.5-209.5-226.3-219.2-208.61822.42069.52016.42018.11868.81479.300-118.3-327.6-327.6-777.61822.41936.01806.91791.81649.61270.60-133.5-328.1-553.9-826.8-986.34.8 4.8 仿真语言仿真语言 开发仿真系统的专用语言,能够提供离散事件仿真系统需要的各种功能：产生随机数和具有特定分布性质的随机变量；推进仿真时间；根据事件表确定下次事件，并把控制

31、转移到相应的程序模块；往事件表中增加记录或从表中删除记录；收集和分析数据；生成仿真结果报表；检查各种错误等。 具有较长发展历史的仿真语言有GPSS（GPSSII、III、360、V）、GASP、和NGPSS以及SIMSCRIPTII，II-PL-US和SIMULA67等。 发展趋向：连续、离散混合仿真; 多种程序结构; 综合化、一体化。 返回4.9 4.9 可视化仿真技术可视化仿真技术 可视化仿真的实质是采用图形或图象处理方式对仿真计算过程的跟踪、控制和结果的后处理，同时实现仿真软件界面的可视化。 OpenGL是包括Microsoft、SGI、Intergragh、IBM以及Intel等在内的ARB机构推出的3D图形应用程序接口的工业标准。Ope