第十五章 随机模拟技术第一节 引言一、模拟的定义模拟是一种数量技术,它利用计算机化的数学模型来表现在某些不确定的条件下所做出的实际决策,来评价一些根据事实及假设所建立的可供选择的行动方案。
二、模拟的意义
1.复杂的实际问题难于用解析理论处理,需要模拟方法提供数值解。
2.理论研究中的某些假设或结论需要经实际系统来检验,计算机模拟可代替费用昂贵的试验。
模拟法分类一、运筹对策法(主要用于军事对策和企业管理对策。如现代化战争的军事演习、新式武器的试验等。最早于 40年代末美国纽曼等人首先用运筹模拟法解决了核屏蔽实验问题。)
型问题。)当然也可用于解决随机很复杂及多变量积分,此法主要用于似的有
。易见,近,及矩形内总点数分的点数于曲线下阴影部均匀选随机点,计算落
,可在矩形内计算方法。例如求一种特殊的数值二、蒙特卡洛法(这是
f ( x )
)(
N
f ( x ) dx,
)(
f ( x ) dx
N
f ( x ) dx
b
a
b
a
b
a
abc
n
abc
n
Nn
a b
c)(xf
x
三、系统模拟法(是用数字对含有随机变量的系统进行模拟,可看作是蒙特卡洛法的应用。一般说来,蒙特卡洛法用于静态计算,而系统模拟法用于动态模型计算。
我们主要讨论此法。)
我们在排队论中讨论了 M/M/C,M/G/1等系统,并用解析方法得出了精确解。但对于到达与服务均为任意分布的排队系统的求解就不可能用那一套公式和方法。
例 1,设某商店顾客到达的时间间隔均匀分布在 1到 10分钟之间,而每一顾客所需要的服务时间均匀分布在 1到 6分钟之间。求顾客在商店所花费的平均时间和售货员空闲时间占全部工作时间的百分比。
分析,到达与服务皆为均匀分布,不能利用 M/M/C或 M/G/1
的公式。但由于问题的特性:
1- 10 1- 10
1- 10
1- 6 1- 6
到达间隔:在 分钟间均匀分布 在 中等可能取值在标有 的牌中任抽取服务时间:在 分钟间均匀分布 在 中等可能取值掷均匀的骰子可用人工方法模拟系统当时的真实情况从而求解。
(用标有 1-10的扑克牌及骰子分别得出用于模拟 20名顾客到达间隔与服务时间的一串数称为随机数,从而推知相关结果。具体怎样做?)
经考察开门后的 20名顾客的被服务情况可知,20名顾客在系统中的全部时间是 68分钟,售货员空闲时间是 55分钟,而售货员从 8点至 9点 57分在班上共 117分钟。于是可得,WS=68/20=3.4(分钟) P0=55/117=0.47
(空闲率过大,可加以调整)
由此例我们初步了解了系统模拟的方法。其中的重要步骤是得到一串关于系统中随机规律的随机数,用以模拟系统的真实情况(故模拟也称仿真),从而求解。而此例中均匀分布的随机数是采用人工方法得到的,即麻烦又不可靠,且局限性很大。所以我们还要寻求产生任意分布随机数的一般方法。
第二节 随机数的产生
RFRfx101
1,[ 0 1 ]
1 x [ 0,1]
( )
0
0 x 0
( ) 0 x 1
1 x 1
[ 0 1 ]
[ 0 1 ]
[ 0 1 ]
R
R
R
R f x
R F x x
R
定义:设 为,上均匀分布的随机变量,
即 的密度为,
其他的分布函数为则 的样本值,即以等概率取自,的一串数称为,上均匀分布的随机数。
一、,区间上均匀分布随机数的产生
2.产生方法
(1)物理方法:一是放射性物质随机蜕变;二是电子管回路的热噪声。(如可将热噪声源装于计算机外部,按其噪声电压的大小表示不同的随机数。此法产生的随机性最好,但产生过程复杂。)
(2)查随机数表 ----”Rand Table”( 1955年由美国兰德公司编制,有随机数 100万个。)随机数表中的数字具有均匀的随机性,没有周期性。使用时,可根据需要任取一段(横或竖)。如需 20个,便可从中取(顺次) 20个,需要几位取几位,随机数表无所谓位数,不能四舍五入。
(3)由递推公式(如同余数公式)在计算机内产生伪随机数。
由于第 i+1个随机数是由第 i个按一定公式推算出来的,故并非真正的随机数。但满足:
a)有较好的随机、均匀性。
b)周期长、重复性差。
c)算法过程不退化(即不能反复出现某一常数。)
d)算法可再现,速度快。
故这是目前最常用的方法。
二、任意概率分布的随机数的产生以上介绍了 R的随机数 r1,r2…… 的产生方法,那么任意分布 X的随机数如何产生?我们说,X的随机数 x1,
x2…… 可以利用 r1,r2…… 得到。那么 X与 R间必有一定关系。这种关系又是什么?
定理,设 R是服从 [0,1]区间上均匀分布的随机变量,X的分布函数为 FX(x),则 X=FX-1(R) 。
( X=FX-1(R) 即 FX(X)=R。即任意分布的随机变量 X被它自己的分布函数作用后所得的随机变量恰为 R。)
1
1 1 1
1
1
( ) ( ) ( )
( ) ( )
()
( ) { } { ( ) } { ( ) } ( ( ) )
0
( ) 0 1 ( )
1 1
(
X
X X X
YX
X
Y X X R X
RX
RX
RF
F R X F R X Y F R
F x F x
Y F R
F x P Y x P F R x P R F x F F x
x
F x x x F x
x
FF
当然与之等价的说法即 被 作用后所得的随机变量恰为 。要证,只要证明左边的分布函数 等于 。
证明:记,则
0
由于 =,而0 1,
( ) ) ( )
X
x F x
YX?
=
即 (分布函数的唯一性)
此定理说明:因为 x轴上的点经 FX(x)映射到 y轴的 [0,1]上便是的 R取值(如图)。反之,y轴的 [0,1]上 R的点经 FX-1映射到 x轴便是 X的取值。所以若知 R的随机数 r1,r2…… 便可得 X
的随机数 x1,x2…… 其中 xi= FX-1(ri)。
1r XF 1?XF0 0x 1x01 0 01 1)( xFX X的分布任 )( xFR R的分布注:若不加以说明则随机数即指 [0,1]上均匀分布的随机数。
例 2,利用 [0,1]区间均匀分布的随机数 r1,r2……
表示服从负指数分布的随机数。
( ) 1 ( 0 )
( ) 1 1
1
l n( 1 )
x
X
x i x i
Xi
ii
X F x e x
R F x r e r i e
xr
解:设 服从负指数分布,则,
由,,或即为所求。
例 3,已知 X的概率分布如右表,
试根据 [0,1]区间上均匀分布 R的随机数列 36,55,
70,38,36,98,
50,95,92,67。
产生 X的随机数列。
X P(x) FX(x) 对应的随机数
0
1
2
3
4
5
0.23
0.30
0.30
0.10
0.05
0.02
0.23
0.53
0.83
0.93
0.98
1.00
00— 23
24— 53
54— 83
84— 93
94— 98
99— 100
分析,所给概率分布如图,而由 R转到 X需用分布函数 FX(如图)。用数表表达即累积概率。
解,先求出 X的累积概率即 FX(x)如右表,然后由
X=FX-1(R)得 X的随机数。
( 注,当随机数落在交界点上,如 98,规定属于前一个范围,当然也可以规定属于后一个范围,只要一致即可。)
R X=FX-1(R)
36
55
70
38
36
98
50
95
92
67
1
2
2
1
1
4
1
4
3
2
第三节 系统的模拟一、排队系统的模拟例 4,有一单服务台的排队系统,根据经验资料知道到达的间隔时间和服务时间的概率分布如下表,其他条件符合标准情形。
到达间隔概率累积概率对应随机数 (a)
2 0.4
6 0.3
10 0.2
14 0.1
服务时间概率累积概率对应随机数 (b)
1 0.4
3 0.4
5 0.2
0.4
0.7
0.9
1.0
0.0-0.4
0.4-0.7
0.7-0.9
0.9-1.0
0.4
0.8
1.0
0.0-0.4
0.4-0.8
0.8-1.0
(1)今由随机数表任选两组随机数
RNa,902,321,211,021,198,383,107,799,439
RNb,612,484,048,605,583,773,054,853,313,200
试根据这两组随机数分别产生表示开始十位顾客的到达间隔时间的随机数 AT和表示服务时间的随机数 ST。
(2)模拟这个排队系统的运行情况,这一阶段共运行多少分钟?
(3)求系统空闲的概率 P0。
顾客
RNa AT 进入系统时刻
RNb ST 服务时间 系统空闲时间开始 结束
1 — 612
2 902 484
3 321 048
4 211 605
5 021 583
6 198 773
7 383 054
8 107 853
9 799 313
10 439 200
—
14
2
2
2
2
2
2
10
6
0
14
16
18
20
22
24
26
36
42
3
3
1
3
3
3
1
5
1
1
0
14
17
18
21
24
27
28
36
42
3
17
18
21
24
27
28
33
37
43
11
3
5
解,先求出到达间隔与服务时间的累积概率。
(见上表)。然后求出 AT,ST及求解 (2),(3)所需的有关数据。
求解结果:
(1)AT与 ST见表
(2)系统此阶段( 10位顾客)共运行 43分钟
(3)空闲率 P0=( 11+3+5) /43=19/43=0.44
二、存贮系统的模拟例 5,有某种货物的存贮系统,市场对这种货物的需求量和订货提前期都是随机的,它们的概率分布如下:
需求量 概率 累积概率
0 0.02
1 0.08
2 0.22
3 0.34
4 0.18
5 0.09
6 0.07
提前期 概率 累积概率
1 0.23
2 0.45
3 0.17
4 0.09
5 0.06
0.02
0.10
0.32
0.66
0.84
0.93
1.00
0.23
0.68
0.85
0.94
1.00
现在考虑订货、存贮、缺货损失三项费用:订货费用每次 25元,订货量每次 20单位,订货点为
15单位。(即存货低于 15单位时订货,但已订货未到前不再订)存贮费每件每周 10元,缺货损失费每件每周 500元。对于缺货,货到后不补,设开始时存货为 20单位。试利用所给随机数 R1(在下表内)模拟需求量,R2( 50,86,15…… )模拟订货提前期。模拟 14周的运行情况:
并求订货费用、存贮费用、缺货费用以及周平均费用。
周需求 到货量存贮量订货 缺货量随机数
R1
需求量 是否订货提前期
R2 提前期
0 20
1 68
2 52
3 90
4 59
5 08
6 72
7 44
—
4
3
5
3
1
4
3
—
16
13 √ 50 2
8
26
24
20
17
周需求 到货量存贮量订货 缺货量随机数 R1
需求量是否订货提前期
R2 提前期
8 95
9 81
10 94
11 28
12 89
13 60
14 03
6
4
6
2
5
3
1
11 √ 86 4
7
1
0 1
15
12 √ 15 1
31
可 求 得,订 货 费 用 25 3=7 5,存 贮 费 用 10 200 =20 00,缺 货 费 用 500 1=5 00
1周 平 均 费 用 ( 75+2000+500 ) =183.9 。
14
二、模拟的意义
1.复杂的实际问题难于用解析理论处理,需要模拟方法提供数值解。
2.理论研究中的某些假设或结论需要经实际系统来检验,计算机模拟可代替费用昂贵的试验。
模拟法分类一、运筹对策法(主要用于军事对策和企业管理对策。如现代化战争的军事演习、新式武器的试验等。最早于 40年代末美国纽曼等人首先用运筹模拟法解决了核屏蔽实验问题。)
型问题。)当然也可用于解决随机很复杂及多变量积分,此法主要用于似的有
。易见,近,及矩形内总点数分的点数于曲线下阴影部均匀选随机点,计算落
,可在矩形内计算方法。例如求一种特殊的数值二、蒙特卡洛法(这是
f ( x )
)(
N
f ( x ) dx,
)(
f ( x ) dx
N
f ( x ) dx
b
a
b
a
b
a
abc
n
abc
n
Nn
a b
c)(xf
x
三、系统模拟法(是用数字对含有随机变量的系统进行模拟,可看作是蒙特卡洛法的应用。一般说来,蒙特卡洛法用于静态计算,而系统模拟法用于动态模型计算。
我们主要讨论此法。)
我们在排队论中讨论了 M/M/C,M/G/1等系统,并用解析方法得出了精确解。但对于到达与服务均为任意分布的排队系统的求解就不可能用那一套公式和方法。
例 1,设某商店顾客到达的时间间隔均匀分布在 1到 10分钟之间,而每一顾客所需要的服务时间均匀分布在 1到 6分钟之间。求顾客在商店所花费的平均时间和售货员空闲时间占全部工作时间的百分比。
分析,到达与服务皆为均匀分布,不能利用 M/M/C或 M/G/1
的公式。但由于问题的特性:
1- 10 1- 10
1- 10
1- 6 1- 6
到达间隔:在 分钟间均匀分布 在 中等可能取值在标有 的牌中任抽取服务时间:在 分钟间均匀分布 在 中等可能取值掷均匀的骰子可用人工方法模拟系统当时的真实情况从而求解。
(用标有 1-10的扑克牌及骰子分别得出用于模拟 20名顾客到达间隔与服务时间的一串数称为随机数,从而推知相关结果。具体怎样做?)
经考察开门后的 20名顾客的被服务情况可知,20名顾客在系统中的全部时间是 68分钟,售货员空闲时间是 55分钟,而售货员从 8点至 9点 57分在班上共 117分钟。于是可得,WS=68/20=3.4(分钟) P0=55/117=0.47
(空闲率过大,可加以调整)
由此例我们初步了解了系统模拟的方法。其中的重要步骤是得到一串关于系统中随机规律的随机数,用以模拟系统的真实情况(故模拟也称仿真),从而求解。而此例中均匀分布的随机数是采用人工方法得到的,即麻烦又不可靠,且局限性很大。所以我们还要寻求产生任意分布随机数的一般方法。
第二节 随机数的产生
RFRfx101
1,[ 0 1 ]
1 x [ 0,1]
( )
0
0 x 0
( ) 0 x 1
1 x 1
[ 0 1 ]
[ 0 1 ]
[ 0 1 ]
R
R
R
R f x
R F x x
R
定义:设 为,上均匀分布的随机变量,
即 的密度为,
其他的分布函数为则 的样本值,即以等概率取自,的一串数称为,上均匀分布的随机数。
一、,区间上均匀分布随机数的产生
2.产生方法
(1)物理方法:一是放射性物质随机蜕变;二是电子管回路的热噪声。(如可将热噪声源装于计算机外部,按其噪声电压的大小表示不同的随机数。此法产生的随机性最好,但产生过程复杂。)
(2)查随机数表 ----”Rand Table”( 1955年由美国兰德公司编制,有随机数 100万个。)随机数表中的数字具有均匀的随机性,没有周期性。使用时,可根据需要任取一段(横或竖)。如需 20个,便可从中取(顺次) 20个,需要几位取几位,随机数表无所谓位数,不能四舍五入。
(3)由递推公式(如同余数公式)在计算机内产生伪随机数。
由于第 i+1个随机数是由第 i个按一定公式推算出来的,故并非真正的随机数。但满足:
a)有较好的随机、均匀性。
b)周期长、重复性差。
c)算法过程不退化(即不能反复出现某一常数。)
d)算法可再现,速度快。
故这是目前最常用的方法。
二、任意概率分布的随机数的产生以上介绍了 R的随机数 r1,r2…… 的产生方法,那么任意分布 X的随机数如何产生?我们说,X的随机数 x1,
x2…… 可以利用 r1,r2…… 得到。那么 X与 R间必有一定关系。这种关系又是什么?
定理,设 R是服从 [0,1]区间上均匀分布的随机变量,X的分布函数为 FX(x),则 X=FX-1(R) 。
( X=FX-1(R) 即 FX(X)=R。即任意分布的随机变量 X被它自己的分布函数作用后所得的随机变量恰为 R。)
1
1 1 1
1
1
( ) ( ) ( )
( ) ( )
()
( ) { } { ( ) } { ( ) } ( ( ) )
0
( ) 0 1 ( )
1 1
(
X
X X X
YX
X
Y X X R X
RX
RX
RF
F R X F R X Y F R
F x F x
Y F R
F x P Y x P F R x P R F x F F x
x
F x x x F x
x
FF
当然与之等价的说法即 被 作用后所得的随机变量恰为 。要证,只要证明左边的分布函数 等于 。
证明:记,则
0
由于 =,而0 1,
( ) ) ( )
X
x F x
YX?
=
即 (分布函数的唯一性)
此定理说明:因为 x轴上的点经 FX(x)映射到 y轴的 [0,1]上便是的 R取值(如图)。反之,y轴的 [0,1]上 R的点经 FX-1映射到 x轴便是 X的取值。所以若知 R的随机数 r1,r2…… 便可得 X
的随机数 x1,x2…… 其中 xi= FX-1(ri)。
1r XF 1?XF0 0x 1x01 0 01 1)( xFX X的分布任 )( xFR R的分布注:若不加以说明则随机数即指 [0,1]上均匀分布的随机数。
例 2,利用 [0,1]区间均匀分布的随机数 r1,r2……
表示服从负指数分布的随机数。
( ) 1 ( 0 )
( ) 1 1
1
l n( 1 )
x
X
x i x i
Xi
ii
X F x e x
R F x r e r i e
xr
解:设 服从负指数分布,则,
由,,或即为所求。
例 3,已知 X的概率分布如右表,
试根据 [0,1]区间上均匀分布 R的随机数列 36,55,
70,38,36,98,
50,95,92,67。
产生 X的随机数列。
X P(x) FX(x) 对应的随机数
0
1
2
3
4
5
0.23
0.30
0.30
0.10
0.05
0.02
0.23
0.53
0.83
0.93
0.98
1.00
00— 23
24— 53
54— 83
84— 93
94— 98
99— 100
分析,所给概率分布如图,而由 R转到 X需用分布函数 FX(如图)。用数表表达即累积概率。
解,先求出 X的累积概率即 FX(x)如右表,然后由
X=FX-1(R)得 X的随机数。
( 注,当随机数落在交界点上,如 98,规定属于前一个范围,当然也可以规定属于后一个范围,只要一致即可。)
R X=FX-1(R)
36
55
70
38
36
98
50
95
92
67
1
2
2
1
1
4
1
4
3
2
第三节 系统的模拟一、排队系统的模拟例 4,有一单服务台的排队系统,根据经验资料知道到达的间隔时间和服务时间的概率分布如下表,其他条件符合标准情形。
到达间隔概率累积概率对应随机数 (a)
2 0.4
6 0.3
10 0.2
14 0.1
服务时间概率累积概率对应随机数 (b)
1 0.4
3 0.4
5 0.2
0.4
0.7
0.9
1.0
0.0-0.4
0.4-0.7
0.7-0.9
0.9-1.0
0.4
0.8
1.0
0.0-0.4
0.4-0.8
0.8-1.0
(1)今由随机数表任选两组随机数
RNa,902,321,211,021,198,383,107,799,439
RNb,612,484,048,605,583,773,054,853,313,200
试根据这两组随机数分别产生表示开始十位顾客的到达间隔时间的随机数 AT和表示服务时间的随机数 ST。
(2)模拟这个排队系统的运行情况,这一阶段共运行多少分钟?
(3)求系统空闲的概率 P0。
顾客
RNa AT 进入系统时刻
RNb ST 服务时间 系统空闲时间开始 结束
1 — 612
2 902 484
3 321 048
4 211 605
5 021 583
6 198 773
7 383 054
8 107 853
9 799 313
10 439 200
—
14
2
2
2
2
2
2
10
6
0
14
16
18
20
22
24
26
36
42
3
3
1
3
3
3
1
5
1
1
0
14
17
18
21
24
27
28
36
42
3
17
18
21
24
27
28
33
37
43
11
3
5
解,先求出到达间隔与服务时间的累积概率。
(见上表)。然后求出 AT,ST及求解 (2),(3)所需的有关数据。
求解结果:
(1)AT与 ST见表
(2)系统此阶段( 10位顾客)共运行 43分钟
(3)空闲率 P0=( 11+3+5) /43=19/43=0.44
二、存贮系统的模拟例 5,有某种货物的存贮系统,市场对这种货物的需求量和订货提前期都是随机的,它们的概率分布如下:
需求量 概率 累积概率
0 0.02
1 0.08
2 0.22
3 0.34
4 0.18
5 0.09
6 0.07
提前期 概率 累积概率
1 0.23
2 0.45
3 0.17
4 0.09
5 0.06
0.02
0.10
0.32
0.66
0.84
0.93
1.00
0.23
0.68
0.85
0.94
1.00
现在考虑订货、存贮、缺货损失三项费用:订货费用每次 25元,订货量每次 20单位,订货点为
15单位。(即存货低于 15单位时订货,但已订货未到前不再订)存贮费每件每周 10元,缺货损失费每件每周 500元。对于缺货,货到后不补,设开始时存货为 20单位。试利用所给随机数 R1(在下表内)模拟需求量,R2( 50,86,15…… )模拟订货提前期。模拟 14周的运行情况:
并求订货费用、存贮费用、缺货费用以及周平均费用。
周需求 到货量存贮量订货 缺货量随机数
R1
需求量 是否订货提前期
R2 提前期
0 20
1 68
2 52
3 90
4 59
5 08
6 72
7 44
—
4
3
5
3
1
4
3
—
16
13 √ 50 2
8
26
24
20
17
周需求 到货量存贮量订货 缺货量随机数 R1
需求量是否订货提前期
R2 提前期
8 95
9 81
10 94
11 28
12 89
13 60
14 03
6
4
6
2
5
3
1
11 √ 86 4
7
1
0 1
15
12 √ 15 1
31
可 求 得,订 货 费 用 25 3=7 5,存 贮 费 用 10 200 =20 00,缺 货 费 用 500 1=5 00
1周 平 均 费 用 ( 75+2000+500 ) =183.9 。
14