计算机模拟（计算机仿真）

计 算 机 模 拟（计算机仿真）一、定义在一定的假设条件下，运用数学计算模拟系统的运行，称为数学模拟。现代的数学模拟都是在计算机上进行的，故又称为计算机模拟.二、基本步骤一 一 一形成数学模型，找出组成部件，并描述有关变量或参数；确定变量之间的函数关系程序验证,检验程序与数学模型之间的一致性，以及输入量的合理性进行模拟试验，对给定的输入在计算机上执行程序数据分析，收集和整理试验结果并作出解释检验由模型所得的结果与真实系统的性能数据的一致性程度两点说明：1.由于系统本身是随时间变化的，或者具有随机性，对真实系统数据和模型试验结果的比较常常需要采用时间序列分析方法或统计分析的方法1.1 时间序列分析方法1）定义：利用按时间顺序排列的一组数字序列，应用数理统计方法加以处理，以预测未来事物的发展的方法叫时间序列分析方法2）基本原理：应用过去数据，就能推测事物的发展趋势；考虑到事物发展的随机性，利用统计分析中加权平均法对历史数据进行处理。1.2 统计分析1）定义：统计分析常指对收集到的有关数据资料进行整理归类并进行解释的过程。是统计工作中统计设计、资料收集、整理汇总、统计分析、信息反馈五个阶段最关键的一步。2）统计分析的步骤：收集数据；整理数据；分 析 数 据（先用描述统计的方法计算出反映数据集中趋势、离散程度和相关强度的具有外在代表性的指标；再用推断统计的方法对数据进行处理，以样本信息推断总体情况,并分析和推测总体的特征和规律。2.按功能将统计分析方法分为描述统计和推断统计。2.1 描述统计:将研究中所得的数据加以整理、归类、简化或绘制成图表，以此描述和归纳数据的特征及变量之间的关系的一种最基本的统计方法。描述统计主要涉及数据的集中趋势、离散程度和相关强度，最常用的指标有平均数、标准差（a x）、相关系数（r）等；2.2 推断统计:指用概率形式来决断数据之间是否存在某种关系及用样本统计值来推测总体特征的一种重要的统计方法。推断统计包括总体参数估计假设检验，最常用的方法有Z检验、T检验、卡方检验等。三、模拟语言模拟语言一般是在别的通用编程语言的基础上建立的，它需要自己的编译程序进行预编译，把模拟语言程序转换成通用编程语言程序，再经过一次编译，转换成计算机上可执行程序。计算机模拟语言是一种高级的描述系统模型的编程语言，它提供表示系统模型中许多基本单元、部件和调度操作的模块。用户用它可较方便地确定模型的基本结构，只要再添加一些辅助的程序就可编出模拟程序。模拟语言可分为离散事件模拟语言（如GPSS及其各种改型、SIMSCRIPT.GASD、CSL、SIMULA等）和连续系统模拟语言（如DARE、ACSL、CSS、CSSL等）两大类型。四、基本方法：蒙 特 卡 洛（Monte Carlo）方法1.蒙 特 卡 洛（Monte Carlo）方法的定义：是一种应用随机数来进行计算机模拟的方法。此方法对研究的系统进行随机观察抽样，通过对样本值的观察统计，求得所研究系统的某些参数。2.用蒙特卡洛方法进行计算机模拟的步骤：设计一个能正确反映系统各部分运行时的逻辑关系的逻辑框图，即模拟模型通过具有各种概率分布的模拟随机数来模拟随机现象模拟随机数：在MATLAB软件中，可以直接产生满足各种分布的随机数的计算机命令，命令分别如下：1)产 生m x n阶 ab均 匀 分 布U(a,b)的 随 机 数 矩 阵：unifrnd(a,b,m,n)产 生 一 个 a,b均 匀 分 布 的 随 机 数：unifmd(a,b)当只知道一个随机变量取值在(a,b)内，但不知道(也没理由假设)它在何处取值的概率大,在何处取值的概率小，就只好用U(a,b)来模拟它。2)产 生m x n阶 0,1均匀分布的随机数矩阵:rand(m,n)产 生 一 个 0,1均 匀 分 布 的 随 机 数:rand3)产 生m*n阶 均 值 为 ，方 差 为cr的正态分布的随机数矩阵:normrnd(,cr,m,n)产 生 一个 均 值 为 从，方 差 为 b的 正 态 分 布 的 随 机 数：normrnd(,cr)当研究对象视为大量相互独立的随机变量之和，且其中每一种变量对总和的影响都很小时，可以认为该对象服从正态分布。机械加工得到的零件尺寸的偏差、射击命中点与目标的偏差、各种测量误差、人的身高、体重等，都可近似看成服从正态分布。4)产 生m x n阶 期 望 值 为M的 指 数 分 布 的 随 机 数 矩 阵：exprnd(,m,n)若 连 续 型 随 机 变 量X的概率密度函数为，,、工 e x 0,/(%)=0 x 0其中几0为 常 数，则 称X服 从 参 数 为4的 指 数 分 布。指 数 分 布 的 期 望 值 为%。在MATLAB中，产 生 参 数 为；I的 指 数 分 布 的 命 令 为e x p rn d()。A指数分布在排队论、可靠性分析中有广泛应用，排队服务系统中顾客到达率为常数时的到达间隔、故障率为常数时零件的寿命都服从指数分布。5)产 生m x n阶 参 数 为X的 泊 松 分 布 的 随 机 数 矩 阵：poissmd(2,m,n)设 离 散 型 随 机 变 量X的 所 有 可 能 取 值 为0,1,2,，且取各个值的概率为/CkeAP(X=k)=,k=0,1,2,k其 中/l0为 常 数，则 称X服 从 参 数 为X的泊松分 布。帕 松 分 布 的 期 望 值 为X。泊松分布在排队系统、产品检验、天文、物理等领域有广泛应用。指 数 分 布 与 泊 松 分 布 的 关 系：若 相 继 两 个 事 件 出 现 的 间 隔 时 间 服 从 参 数 为 的指 数 分 布，则 在 单 位 时 间 间 隔 内 事 件 出 现 的 次 数 服 从 参 数 为X的 泊 松 分 布.即 单 位时 间 内 该 事 件 出 现k次 的 概 率 为P(X=A)=%J =0,l,2,反 之 亦 然。k敏.顾 客 到 达 某 商 店 的 间 隔 时 间 服 从 参 数 为0.1的指数分布:指两个顾客到达商店 的 平 均 间 隔 时 间 是10个 单 位 时 间，即 平 均10个 单 位 时 间 到 达1个 顾 客，则指数 分 布 的 均 值 为1/0.1=10,顾 客 到 达 的 间 隔 时 间 可 用expmd(10)模拟该 商 店 在 单 位 时 间 内 到 达 的 顾 客 数 服 从 参 数 为0.1的 帕 松 分 布：指一个单位时间内平均到达0.1个顾客，即平均10个单位时间到达一个顾客。例2.敌坦克分队对我方阵地实施突袭，其到达规律服从泊松分布，平均每分钟到达4辆.(1)模拟敌坦克在3分钟内到达目标区的数量，以及在第1、2、3分钟内各到达几辆坦克；(2)模拟在3分钟内每辆敌坦克的到达时刻。(1)用poissmd(4)进行模拟；(2)坦克到达的间隔时间应服从参数为4的负指数分布，用exprnd(1/4)模拟3.通常蒙特卡罗方法可以粗略地分成两类一类是所求解的问题本身具有内在的随机性，借助计算机的运算能力可以直接模拟这种随机的过程。例如在核物理研究中，分析中子在反应堆中的传输过程：中子与原子核作用受到量子力学规律的制约，人们只能知道它们相互作用发生的概率,却无法准确获得中子与原子核作用时的位置以及裂变产生的新中子的行进速率和方向，科学家依据其概率进行随机抽样得到裂变位置、速度和方向，这样模拟大量中子的行为后，经过统计就能获得中子传输的范围，作为反应堆设计的依据；另一种类型是所求解问题可以转化为某种随机分布的特征数，比如随机事件出现的概率，或者随机变量的期望值。通过随机抽样的方法，以随机事件出现的频率估计其概率，或者以抽样的数字特征估算随机变量的数字特征，并将其作为问题的解。这种方法多用于求解复杂的多维积分问题。4.在解决实际问题的时候应用蒙特卡罗方法主要有两部分工作首先需要产生各种概率分布的随机变量；再用统计方法把模型的数字特征估计出来，从而得到实际问题的数值解。5.用蒙特卡洛方法进行计算机模拟的实例例3 在我方某前沿防守地域，敌人以一个炮排(含两门火炮)为单位对我方进行干扰和破坏。为躲避我方打击，敌方对其阵地进行了伪装并经常变换设计地点经过长期观察发现，我方指挥所对敌方目标的指示有50%是准确的，而我方火力单位，在指示正确时，有1/3的射击效果能毁伤敌人一门火炮，有1/6的射击效果能全部消灭敌人。现在希望能用某种方式把我方将要对敌人实施的20次打进结果显现出来，确定有效射击的比率及毁伤敌方火炮的平均值。分析：这是一个概率问题，可以通过理论计算得到相应的概率和期望值。但这样只能给出作战行动的最终静态结果，而显示不出作战行动的动态过程。为了能显示我方20次射击的过程，现采用模拟的方式问题分析需要模拟出以下两件事1)观察所对目标的指示正确与否模拟实验有两种结果，每一种结果出现的概率都是1/2.因此，可用投掷一枚硬币的方式予以确定，当硬币出现正面指示正确，反之不正确。2)当指示正确时，我方火力单位的射击结果情况模拟实验有三种结果：毁伤一门火炮的可能性是1/3（即 2/6）,毁伤两门的可能性是1/6,没能毁伤火炮的可能性是1/2（即3/6）。这时可以用投掷骰子的方法来确定：如果出现的是1、2、3 这三个点:则认为没能击中敌人；如果出现的是4、5 点：则认为毁伤敌人一门火炮；如果出现的是6 点：则认为毁伤敌人两门火炮。符号假设_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _i 要模拟的打击次数kl 没击中敌人火炮的射击总数k2 击中敌人一门火炮的设计总数k3 击中敌人两门火炮的射击总数E 有效射击比率E120次射击平均每次毁伤敌人的火炮数模拟框图模拟结果试验序号投硬币结 果指示正确指示不正掷骰子结消灭敌人火炮数确果0121正V4V2正V4V3反VV4正V1V5正V2V6反VV7正V3V8正V6V9反VV10反VV11正V2V12反VV13正V3V14反VV15正V6V16正V4V17正V2V18正V4V19反VV20正V6V1 3 4 3从以上模拟结果可以算出：E=7/2 0=0.3 5 ,E l=O x +l x +2 x =0.5.2 0 2 0 2 0理论计算假设1)j=0 x2 0在MATLAB软件包中编程，共需三个M 一文件：randlp.m,mmylp.m,lpconst.主程序为randlp.m%mylp.mfunction z=mylp(x)%目标函数z=2*x(1 )A2+x(2)A2-x(1 )*x(2)-8*x(1 )-3*x(2);%转化为求最小值问题%lpconst.mfunction lpc=lpconst(x)%约 束条件if 3*x(l)+x(2)-10=0.5%约 束条件lpc=l;elselpc=0end%randlp.mfunction sol,r 1 ,r2J=randlp(a,b,n)%随机模拟解非线性规划debug=l;a=0;b=10;%试验点下界和上界n=1000;%试验点个数r 1 =unifrnd(a,b,n,1);r2=unifrnd(a,b,n,1);%n&l 阶 的 a,b 均匀分布随机数矩阵sol=rl(l)r2(l);z0=inf;for i=l:nxl=rl(i);x2=r2(i);lpc=lpconst(xl x2)；if lpc=lz=mylp(xl x2)；if zz0z0=z;sol=xl x2;endendend五、模拟类型1.离散系统模拟1.1 离散时间模型的模拟离散时间模型中的时间表示为整数序列（代表某一时间单位的整数倍），只考虑系统在这些时刻上的状态变化。这种模型的一个典型模拟程序包括下列步骤：1）置 模拟时间T的初始值为to;2）置状态变量的初始值；3）给出当前模拟时间输入变量的值后，根据模型中的状态转移函数，确定在下一时刻T=t+h状