数学建模社区-数学中国
标题:
排队论模型(八):Matlab 生成随机数、排队模型的计算机模拟
[打印本页]
作者:
浅夏110
时间:
2020-6-13 09:35
标题:
排队论模型(八):Matlab 生成随机数、排队模型的计算机模拟
1 产生给定分布的随机数的方法
0 S( l0 c) m7 M) k# X
Matlab 可以产生常用分布的随机数。下面我们介绍按照给定的概率分布产生随机数的一般方法,这些方法都以U(0,1) 分布的随机变量为基础。
5 ^4 M) G1 M) `2 u; L
T0 u! l2 Z# L% P0 R5 g! s
(i)反变换法
; e5 L# t# N) f5 s/ \) N
定理 设 X 是一个具有连续分布函数 F(x) 的随机变量,则 F(X ) 在 [0,1] 上服 从均匀分布。
' e1 ~* j7 i1 H" b0 Z
& {$ h# E% i# i6 d$ B
: z" v; t- n& {& ^! f- J0 x# w
8 D- c& {% ^, i( V. b2 @ m
8 T* ]% m5 u+ k. ^0 [
! C, K& [8 ~! G
(ii)卷积法
2 r( K+ v6 V" t0 }- s8 C1 P% ^
* U1 z0 j, m- |0 s; m
7 X! x4 \9 x1 e
0 a4 |/ V% D( ?3 O/ P% X- k
(iii)取舍法
3 y m6 o3 g* m7 J, D, N! T5 m S H
若随机变量 X 在有限区间(a,b) 内变化,但概率密度 f (x)具有任意形式(甚至没 有解析表达式),无法用前面的方法产生时,可用取舍法。一种比较简单的取舍法的步 骤是:
$ {8 I: [ k2 \8 M2 Z' X
. r, I& w( ~+ |
% t, X" x @$ C5 O- B% |9 t
6 B- W, j9 S% J" t
2 排队模型的计算机模拟
$ ~6 \" n) B( Y2 E H; A
2.1 确定随机变量概率分布的常用方法
% q/ H% s* X1 [$ L' `' t
在模拟一个带有随机因素的实际系统时,究竟用什么样的概率分布描述问题中的随 机变量,是我们总是要碰到的一个问题,下面简单介绍确定分布的常用方法:
4 a- Y. f% S$ C. ~
/ h8 m$ P4 E P7 D3 d
【1 】根据一般知识和经验,可以假定其概率分布的形式,如顾客到达间隔服从指数 分布 Exp(λ) ;产品需求量服从正态分布 ;订票后但未能按时前往机场登机 的人数服从二项分布 B(n, p) 。然后由实际数据估计分布的参数 λ,μ,σ 等,参数估计 可用极大似然估计、矩估计等方法。
0 h4 R5 t2 S! N) A) X' _, t
% H" d1 j$ o F; }; q3 \& {2 w! S' D. @
【2】 直接由大量的实际数据作直方图,得到经验分布,再通过假设检验,拟合分布 函数,可用 检验等方法。 3 o 既缺少先验知识,又缺少数据时,对区间(a,b) 内变化的随机变量,可选用 Beta 分布(包括均匀分布)。先根据经验确定随机变量的均值 μ 和频率最高时的数值(即密度函数的最大值点)m ,则 Beta 分布中的参数 可由以下关系求出:
- M4 m ~: U( L) G! x
+ r& T# v- c8 M0 ^; i' S" `, I0 s. f% Q
N# V) B: v# e
" @- K5 M: O5 N/ R) O
2 .2 计算机模拟
0 x! Q! {4 D6 V: ^
当排队系统的到达间隔时间和服务时间的概率分布很复杂时,或不能用公式给出 时,那么就不能用解析法求解。这就需用随机模拟法求解,现举例说明。
8 a9 n: V5 N6 d. _: p# D
. ?+ w4 X1 v3 ~( b2 @5 w& N
例 14 设某仓库前有一卸货场,货车一般是夜间到达,白天卸货,每天只能卸货 2 车,若一天内到达数超过 2 车,那么就推迟到次日卸货。根据表 3 所示的数据,货车到 达数的概率分布(相对频率)平均为 1.5 车/天,求每天推迟卸货的平均车数。
4 C$ R9 E d3 Q& r' @' i
, \( O; J& W, C7 E
- {7 Z3 }8 s t' O6 S. L
( ?' x! a; r u/ Z2 ^8 N1 R
解 这是单服务台的排队系统,可验证到达车数不服从泊松分布,服务时间也不服 从指数分布(这是定长服务时间)。 随机模拟法首先要求事件能按历史的概率分布规律出现。模拟时产生的随机数与事 件的对应关系如表 4。
+ g/ X( O9 g7 C, v. m
0 A( I: Y2 E1 S, j; [
0 q. N4 |" O. R( T# e
3 M( c! Z! Q/ s2 L; W; G
我们用 a1 表示产生的随机数,a2 表示到达的车数,a3 表示需要卸货车数,a4 表 示实际卸货车数,a5 表示推迟卸货车数。编写程序如下:
* |% i7 t% W' ]! |" ^
1 D( w! z& s2 J7 u4 ?
clear
$ w/ U& _6 q. R, R" H2 b
rand('state',sum(100*clock));
" t. I) D Q/ c. [7 T; G. J- v
n=50000;
; Z% O' ~: E* Q ]8 L6 ^8 V+ X
m=2
, G4 ~; i$ A0 O3 ^
a1=rand(n,1);
1 ?0 S _! W% |& d4 {4 p# H: h
a2=a1; %a2初始化
8 P/ ^* o& Q& d
a2(find(a1<0.23))=0;
# }0 p% ]- H4 j8 u
a2(find(0.23<=a1&a1<0.53))=1;
; h% T. r6 z, p$ ~+ \# {
a2(find(0.53<=a1&a1<0.83))=2;
0 n9 \ Q4 q+ Y/ l0 P
a2(find(0.83<=a1&a1<0.93),1)=3;
* |2 w0 J0 z; G3 C. j) N6 E
a2(find(0.93<=a1&a1<0.98),1)=4;
7 v* |/ e( v1 `# f3 F; s
a2(find(a1>=0.98))=5;
$ A4 B. y6 e/ s! _3 ]1 n
a3=zeros(n,1);a4=zeros(n,1);a5=zeros(n,1); %a2初始化
1 `& h$ T$ d, F! W2 h! M: E8 j
a3(1)=a2(1);
/ a/ e! c9 S6 B" E' ~
if a3(1)<=m
" b8 l$ c4 ^$ B6 ?$ O+ [
a4(1)=a3(1);a5(1)=0;
* b' w |( W5 ^2 p* @6 h: P2 ]
else
2 H2 F$ @; N& E
a4(1)=m;a5(1)=a2(1)-m;
( {3 N4 k% p% D9 z7 t- T7 t' P
end
$ v8 a& u" q1 m
for i=2:n
# u$ ]. @: O: e
a3(i)=a2(i)+a5(i-1);
2 N1 M C6 v8 n# C# f
if a3(i)<=m
" A( q& g0 S& Z9 c- |
a4(i)=a3(i);a5(i)=0;
6 t9 m$ a2 L _- m) \
else
8 K3 D% E3 l7 O5 A
a4(i)=m;a5(i)=a3(i)-m;
: ?! B% b+ h6 N
end
5 c# Y0 u; _: V
end
2 q! P8 y, f1 a: k* ^" Q
a=[a1,a2,a3,a4,a5];
1 g3 X: o! W k0 K X
sum(a)/n
+ Q2 G2 L: h s
例 15 银行计划安置自动取款机,已知 A 型机的价格是 B 型机的 2 倍,而 A 型机 的性能—平均服务率也是 B 型机的 2 倍,问应该购置 1 台 A 型机还是 2 台 B 型机。 为了通过模拟回答这类问题,作如下具体假设,顾客平均每分钟到达 1 位, A 型 机的平均服务时间为 0.9 分钟, B 型机为 1.8 分钟,顾客到达间隔和服务时间都服从 指数分布,2 台 B 型机采取 M / M / 2 模型(排一队),用前 100 名顾客(第 1 位顾客到 达时取款机前为空)的平均等待时间为指标,对 A 型机和 B 型机分别作 1000 次模拟, 进行比较。
6 {2 j5 s7 L3 @- S$ R& Z0 S
# I6 [ |8 C; |0 Q
" S) k! Q0 K0 C) R: j& w. f
0 z/ g( ^' U& B4 g& A6 e/ l
在模拟 A 型机时,我们用cspan表示到达间隔时间,sspan表示服务时间,ctime 表示到达时间,gtime表示离开时间,wtime表示等待时间。我们总共模拟了m 次, 每次n 个顾客。程序如下:
# _: R3 d& j6 S9 h9 u1 s
# }8 v) A- m* Z7 P
tic
7 A% `' \- C% A1 {. a9 A
rand('state',sum(100*clock));
5 I9 Z c4 D( b; C6 f) L. i
n=100;m=1000;mu1=1;mu2=0.9;
4 h p( w/ v# b3 s
for j=1:m
8 A* }2 e. W( K4 p- U; a
cspan=exprnd(mu1,1,n);sspan=exprnd(mu2,1,n);
3 i. T$ L8 R0 p c( L( m S
ctime(1)=cspan(1);
- z5 S7 ]7 N) n% _8 {
gtime(1)=ctime(1)+sspan(1);
" }6 X# m& S) B0 _9 l) h! H$ i
wtime(1)=0;
0 D. K: i Z' ]6 ~! @( P
for i=2:n
/ Y9 H, Q; E4 a! w! U& [
ctime(i)=ctime(i-1)+cspan(i);
- ^7 \* |- v# ?
gtime(i)=max(ctime(i),gtime(i-1))+sspan(i);
. Z* ^0 l! ]: \0 x1 b
wtime(i)=max(0,gtime(i-1)-ctime(i));
6 l$ j8 t T* M" z% a5 l2 X
end
& L* Z' j, ]2 I* I& l) v
result1(j)=sum(wtime)/n;
6 l o2 ~! q& Q% P9 O) u3 Y
end
+ e" C* m' a& q U3 c( {
result_1=sum(result1)/m
( c5 u# o$ F: i8 p& c4 V' a' d
toc
) d+ C9 R5 J9 L: }+ V, G6 T( s5 w( n
类似地,模拟 B 型机的程序如下:
$ C" z- b' P0 p+ y# K
( V/ s3 h% U& l$ a
tic
7 Z% A. J5 L3 `4 N
rand('state',sum(100*clock));
' X T" C9 Y( U$ e* f
n=100;m=1000;mu1=1;mu2=1.8;
7 b% P+ f) W8 f# ]0 v
for j=1:m
9 w A0 N# }( A3 Y
cspan=exprnd(mu1,1,n);sspan=exprnd(mu2,1,n);
/ z/ Y+ u+ M1 H1 h* P
ctime(1)=cspan(1);ctime(2)=ctime(1)+cspan(2);
. z# v( z- M; `6 S
gtime(1:2)=ctime(1:2)+sspan(1:2);
3 \3 v \0 y0 j3 G, J1 O" n
wtime(1:2)=0;flag=gtime(1:2);
" J3 G, x. x5 W& ?1 C% h7 ^, a- ^1 R& Q
for i=3:n
/ \6 f5 n' ~" `" o
ctime(i)=ctime(i-1)+cspan(i);
# r2 V* L- Z( O; T4 K1 Z
gtime(i)=max(ctime(i),min(flag))+sspan(i);
) I3 S: {3 S& t% Q& s# \
wtime(i)=max(0,min(flag)-ctime(i));
3 X9 A) g, y) Q! L9 t% M& Q$ K
flag=[max(flag),gtime(i)];
8 b m! A5 o. Y1 ]" p
end
( j0 x9 ^2 s! {! L# B
result2(j)=sum(wtime)/n;
+ w5 n2 e. ?! [( a5 C
end
# k+ t- E3 q5 I8 S
result_2=sum(result2)/m
5 c* G. Q: a4 [2 Q4 Z$ l% }
toc
7 P( Q# f7 r1 ?
读者可以用下面的程序与上面的程序比较了解编程的效率问题。
T/ i( G5 e6 T7 Y% u% v7 `1 N$ |
' k9 g0 S* M3 x, l
tic
( R8 Q8 j$ S5 b$ ?# j4 _1 p2 Y4 l
clear
; R9 }# n6 y, q0 d. `4 n& y
rand('state',sum(100*clock));
; t9 o3 {% ^* u5 }
n=100;m=1000;mu1=1;mu2=0.9;
5 i9 @% W* a% N5 v( i5 P
for j=1:m
' B! X. z4 o. |! Q: P4 f2 ^1 H
ctime(1)=exprnd(mu1);
- l# a. `6 q8 Q* m6 P/ {3 p5 r
gtime(1)=ctime(1)+exprnd(mu2);
( n& K; m" C& m' a
wtime(1)=0;
2 n; [+ l$ k9 X* v2 G8 ?
for i=2:n
& Y0 A; a: H5 R! p- X$ e" s
ctime(i)=ctime(i-1)+exprnd(mu1);
8 d% V$ O4 { D4 H5 a; D
gtime(i)=max(ctime(i),gtime(i-1))+exprnd(mu2);
% y7 `4 d5 O ^) Q# A, l1 x: P
wtime(i)=max(0,gtime(i-1)-ctime(i));
6 N/ l3 |9 J( K6 i5 W& }
end
3 H: J" ?6 ] O7 i( ^9 x
result(j)=sum(wtime)/n;
( K& Q& w% F" I1 |- i
end
C+ @6 J- l- t4 d0 w) ~8 y9 {8 B
result=sum(result)/m
& U: Y# d2 Y- }" |
toc
' _6 s) e! i1 n9 F, j0 O/ v, D0 F
1. 一个车间内有10台相同的机器,每台机器运行时每小时能创造4元的利润,且平 均每小时损坏一次。而一个修理工修复一台机器平均需4小时。以上时间均服从指数分 布。设一名修理工一小时工资为6元,试求:
4 _8 E/ y4 _! j! U, F- r
. ^4 ~3 n# W+ a( u
(i)该车间应设多少名修理工,使总费用为最小;
- F8 x2 H! L# V5 V% J
1 `! L: [% _$ a$ H8 O
(ii)若要求不能运转的机器的期望数小于4台,则应设多少名修理工;
$ @$ s* [2 k& `+ B* i
0 g( O: o" v! A1 s0 ]
(iii)若要求损坏机器等待修理的时间少于4小时,又应设多少名修理工。
6 q$ g5 v: D; Y
4 x4 D. [# w* i8 c6 l5 M
2. 到达某铁路售票处顾客分两类:一类买南方线路票,到达率为λ1 /小时,另一 类买北方线路票,到达率为λ2 /小时,以上均服从泊松分布。该售票处设两个窗口,各窗口服务一名顾客时间均服从参数 μ = 10 的指数分布。试比较下列情况时顾客分别等 待时间Wq :
& V: Y* q' W/ Y
& y2 r& B) G. s
(i)两个窗口分别售南方票和北方票;
5 x9 f" U/ O E+ f Z4 F* R5 T4 V& A
6 Q5 ]: E$ m+ Q( {' I1 \$ C
(ii)每个窗口两种票均出售。(分别比较 λ1 = λ2 = 2,4,6,8 时的情形)
7 [! W. w. Y: D- O
% U0 n7 {- ^$ ~9 i9 h- k% ~
3. 一名修理工负责5台机器的维修,每台机器平均每2h损坏一次,又修理工修复一 台机器平均需时18.75min,以上时间均服从负指数分布。试求:
F# d& b! p {3 p" y0 R* ]5 `
. D$ w5 y5 @) R# A! }0 k2 _
(1)所有机器均正常运转的概率;
6 h9 ?5 [: N' J: c+ J
( _0 @6 i. K' [+ F" b/ T' j0 u5 v
(2)等待维修的机器的期望数;
* x1 F! u7 W* |' \$ |0 Y* S( K
0 h+ g4 O `: ~
(3)假如希望做到有一半时间所有机器都正常运转,则该修理工最多看管多少台 机器。
5 t4 [( m/ h% j) T- f
6 M4 U. @3 F7 ^' }& `. p
(4)假如维修工工资为8元/h,机器不能正常运转时的损失为40元/h,则该修理工 看管多少台机器较为经济合理。
4 D" h# X$ Z: D1 M8 l4 d
————————————————
, L. s, o# p5 [9 x0 C( ^# O
版权声明:本文为CSDN博主「wamg潇潇」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
+ i, ~; ~$ k2 B4 m- s
原文链接:https://blog.csdn.net/qq_29831163/java/article/details/89738145
9 E- \+ y/ e# x& \8 ^! z
, q, B! W$ G, G( \# R
- g' N% Y8 s/ A1 j4 E9 h' c; ?" o
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5