计算机模拟

1. 问题分析 需要模拟以下两件事：
[1] 观察所对目标的指示正确与否？
模拟试验有 2 种结果：每种结果概率都是 1/2。 与投掷一枚硬币方式相仿，出现硬币正面时为指示正 确，反之为不正确。
[2] 指示正确时，我方的射击效果如何？
模拟试验有 3 种结果：毁伤一门火炮的可能性为 1/3(即 2/6)，毁伤两门的可能性为 1/6，没能毁伤敌火炮 的可能性为 1/2(即 3/6)。 与投掷骰子的方式相仿，若出现是１、２、３点，则 认为没击中敌人；若出现是４、５点，则认为毁伤敌人 一门火炮；若出现是６点，则认为毁伤敌人两门火炮。
（二）蒙特卡洛方法
蒙特卡洛(Monte Carlo)方法是计算机模拟 的基础。 它得名于赌城——摩纳哥的蒙特卡洛。
基本思路： ⑴建立一个概率模型，使所求问题的解正 好是该模型的参数或特征量。 ⑵通过模拟一统计试验，统计某事件发生 的比率，当试验次数很大时比率就近似事件发 生的概率。 ⑶利用所建概率模型，求要估计的参数。
（四）实例：炮击效果的计算
在我方前沿防守地域，敌方以一个炮排(含两 门火炮)为单位对我方干扰和破坏。为躲避我方 打击，敌方的阵地进行伪装并变换射击地点。 经观察发现，我方指挥对敌方目标的指示 50 ％准确， 在指示正确时而我方火力有 1/3 的射击 效果能毁伤敌人一门火炮， 1/6 的射击效果能 有 全部消灭敌人。 现希望能用某种方式把我方将要对敌人实施 的 20 次打击结果显现出来，确定有效射击的比 率及毁伤敌方火炮的平均值。
2 2
。
当研究对象视为大量相互独立的随机变量之和，且其 中每一种变量对总和的影响都很小时，可认为该对象服 从正态分布。 例如，机械加工得到的零件尺寸的偏差、射击命中点 与目标的偏差、各种测量误差、人的身高、体重等，都 可近似看成服从正态分布。
6. exprnd , m, n 产生 m n 阶期望值为 的 指数分布的随机数矩阵；exprnd 产生期望值 为 的指数分布的随机数。 若连续型随机变量 X 的概率密度函数为
静态模拟和动态模拟
静态模拟(static simulation)：通过模拟方法 只要求求得问题的结果的模拟方法。 典型实例是用蒙特卡罗（Monte Carlo）方法 求数值积分的值。 动态模拟(dynamic simulation)：通过模拟方 法研究系统的状态随时间连续变化的情形。 动态模拟又可分为连续系统模拟和离散事件 系统模拟。
1. 符号假设
[1] 观察所对目标的指示正确与否？
模拟试验有 2 种结果：每种结果概率都是 1/2。 与投掷一枚硬币方式相仿，出现硬币正面时为指示正 确，反之为不正确。
[2] 指示正确时，我方的射击效果如何？
模拟试验有 3 种结果：毁伤一门火炮的可能性为 1/3(即 2/6)，毁伤两门的可能性为 1/6，没能毁伤敌火炮 的可能性为 1/2(即 3/6)。 与投掷骰子的方式相仿，若出现是１、２、３点，则 认为没击中敌人；若出现是４、５点，则认为毁伤敌人 一门火炮；若出现是６点，则认为毁伤敌人两门火炮。
E1 0 13 4 3 1 2 =0.5 20 20 20
E=7/20=0.35
5. 理论计算
确 0 观察所对目标指示不正 设： j 1 观察所对目标指示正确 A0：射中敌方火炮的事件；A1：射中敌方一门火炮的事件； A2：射中敌方两门火炮的事件． 则由全概率公式： E = P(A0) = P(j=0)P(A0∣j=0) + P(j=1)P(A0∣j=1) 1 1 1 = 0 0.25 2 2 2 P(A1) = P(j=0)P(A1∣j=0) + P(j=1)P(A1∣j=1) 1 1 1 1 = 0 2 2 3 6 P(A2) = P(j=0)P(A2∣j=0) + P(j=1)P(A2∣j=1) 1 1 1 1 = 0 2 2 6 12 1 1 E1 = 1 2 0.33 6 12
执行后得到 I 2.6681 。 N 值越大越精确。
（三）生成随机数的函数
在 Matlab 中可产生满足各种分布的随机数： 1．unidrnd(N,m,n) 产生服从 1,2,…,N 的等概 率的离散均匀分布的 mn 阶的随机数矩阵； unidrnd(N) 产生相应的随机数。
它的概率分布为 P k 1 N , k 0,1,2,..., N 。 2． binornd (k, p, m, n) 生成参数为 k,p 的 mn 阶的二项分布随机数矩阵；binornd (k, p)生成 k,p 的相应的随机数。 它的概率分布为 P 1 p, P 1 1 p 。
数学模拟-计算机模拟
在一定的假设条件下，运用数学运算模拟系 统的运行，称为数学模拟。现代数学模拟都在 计算机上进行，故也称计算机模拟。 优点：计算机模拟可反复进行，改变系统的 结构和系数都比较容易。 实际的一些带随机因素的复杂系统，用分析 方法建模常需要作许多简化假设，使模型与实 际问题可能相差甚远，导致解答根本应用。这 时，计算机模拟几乎成唯一选择。
7．poissrnd , m, n 产生 m n 阶的服从参数 为 的帕松分布的随机数矩阵； poissrnd 产生 服从参数为 的帕松分布的随机数。 若离散型随机变量 X 的所有可能取值为 0,1,2,…,且取各个值的概率为
P( X k )
k e
6. 结果比较
理论计算和模拟结果的比较
分类 项目 模 理 拟 论 无效射击 ０．６５ ０．７５ 有效射击 ０．３５ ０．２５ 平均值 ０．５ ０．３３
⑶产生在 内的 N 个服从均匀分布的随机 数，若 Nc 个落在 V 内，则 I 的体积 NV / N 。
算例 I
x2 y 2 1
1 x2 dxdy
曲顶柱体 V ： x2 y 2 1， x 2 z 2 1 ， z 0 。 包含 V 的立方体 ： x 1, y 1 ， 0 z 1 。
若顾客到达某商店的间隔时间服从参数为 0.1 的指数分 布，则指两个顾客到达商店的平均间隔时间是 10 个单位 时间，即平均 10 个单位时间到达 1 个顾客。 若该商店在单位时间内到达的顾客数服从参数为 0.1 的 帕松分布，则指一个单位时间内平均到达 0.1 个顾客。
k e
k!
, k 0,1, 2, 。
e x , x 0, f x 其中 0 为常数， x 0. 0,
则称 X 服从参数为 的指数分布。 参数为 的指数分布的期望值为 1 。
排队服务系统中，顾客到达率为常数时的到达间隔、 故障率为常数时零件的寿命都服从指数分布。指数分布 在排队论、可靠性分析中有广泛应用。
5. normrnd , , m, n 产生 m n 阶的服从均值
、方差 的正态分布的随机数矩阵；
normrnd , 产生服从均值 、 方差 的正 态分布的随机数。

正态分布 N ,
2
的密度函数 f x
1 e 2
x 2
P k C p 1 p
k n k nk
, k 0,1,..., n 。
3．unifrnd(a,b,m, n) 产生 mn 阶在 a, b 上服 从均匀分布 U a, b 的随机数矩阵；unifrnd(a,b) 产生在 a, b 上服从均匀分布 U a, b 的随机数。
当只知道一个随机变量取值在 a, b 内，但不知道（也 没理由假设）它在何处取值的概率大，在何处取值的概 率小，就只好用 U a, b 来模拟。
4． rand (m, n) 生成 mn 阶的 0,1 上服从均匀 分布的随机数矩阵； rand 生成在 0,1 上服从均匀 分布的随机数。
计算机数学模拟
（一）模拟(Simulation)概念
1. 模拟：利用物理的、数学的模型来类比、 模仿现实系统及其演变过程，以寻求过程规律 的一种方法。 2. 思路： 建立一个包含所研究系统的主要特 点试验模型。通过对这个试验模型的运行，获 得所要研究系统的必要信息。
物理模拟
形式：对实际系统及其过程用功能相似的实 物系统去模仿。例如，军事演习、船艇实验、 沙盘作业等。 优点：直观。 缺点：花费较大、周期较长，在模型上改变 系统结构和系数都较困难。且许多系统无法进 行物理模拟，如社会经济系统、生态系统等。
的体积为 2 2 ：
clear; N = 1000000; MATLAB 程序 x = unifrnd(-1,1,N,1); p_jifen.m y = unifrnd(-1,1,N,1); z = rand(N,1); v1 = (x.^2+y.^2)<1; v2 = (x.^2+z.^2)<1; Nv = sum(v1&v2); I = Nv/N*4
连续系统动态模拟：研究系统的状态随时间 连续变化的情形，其模型一般是微分方程模 型。建模时首先确定系统的连续状态变量，然 后将它在时间上进行离散化处理, 并由此模拟 系统的运行状态。
离散事件系统模拟：离散事件系统模拟讨论 的是系统状态只在一些离散时间点上，由于随 机事件驱动而发生变化，其模型一般用流程图 或网络来表示。
k!
, k 0,1, 2, ， 其中 0 为常数， 则
则称 X 服从参数为 的帕松分布。 参数为 的帕松分布的期望值为 。
帕松分布在排队系统、产品检验、天文、物理等领域 有广泛应用。
指数分布与帕松分布的关系
如相继两个事件出现的间隔时间服从参数为 的指数分布，则在单位时间间隔内事件出现的次 数服从参数为 的泊松分布。即单位时间内该事 件出现 k 次的概率为 P( X k ) 反之亦然。
用蒙特卡洛方法计算二重积分 计算二重积分 I f x, y dxdy ，f x, y 0 ，