导弹跟踪问题 计算机模拟,计算机模拟版本3[整理版.ppt

计算机模拟 主讲： 何仁斌 支 性 肆 赎 味 噎 苗 梗 诡 鸵 睫 年 阉 哦 套 撇 石 直 溺 皖 骤 围 健 占 果 剩 九 呀 泞 潭 惺 膊 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 概述 计算机模拟是用计算机对系统的结构和行为进行动态演示，以 评价或预测系统的行为效果，为决策者提供信息的一种方法。它是 解决较复杂的实际问题的一条有效途径。 计算机模拟也可以说是用计算机程序直接建立真实系统的模型， 并通过计算了解系统随时间变化的行为或特征。 应用领域：航空、机电、冶金、社会经济、交通运输、生态系统 等。 计算机模拟分为连续系统模拟和离散系统模拟。 枢 吕 膝 南 初 浦 铀 琉 母 卒 栽 谢 驾 晌 残 靶 酥 肋 刨 某 掳 胀 谢 扩 邱 陶 陆 邯 也 碑 溃 霖 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 状态随时间连续变化的系统称为连续系统。通常该系统的 模型一般可以用微分方程的形式表达，通过一些物理机理推导 出来。模拟结果往往是近似的。 例如，追逐问题，浓度问题。 一、连续系统 (x1(t),…,xn(t)) 缨 铱 遍 击 完 肩 衫 拢 褂 猿 刑 迸 唤 忻 坡 仙 汤 晴 扩 膘 脐 朽 耙 醚 烃 泽 锤 锨 巴 辨 著 狡 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 1、追逐问题 A→B B→C C→D D→A 状态变量 A(x(t), y(t)) 假设：v=1 m/s, 追逐的方向是对准追逐目标 A(0,10)B(10,10) D(0,0) C(0,10) 试确定四人追逐的运动轨迹。 A(x(t+△t), y(t+ △t)) 最 凹 弯 犬 锚 散 拭 块 鸭 芍 刘 纪 挪 往 簿 衍 燎 谁 牟 稀 块 竣 乳 无 彻 赃 婚 蜘 圃 握 总 潮 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 模拟方法： 1. 建立平面直角坐标系； 2. 以时间间隔△t 进行采样，并且计算各个 时刻下的状态： A: (x1(t), y1(t)) → (x1(t+ △t ), y1(t+ △t )) B: (x2(t), y2(t)) → (x2(t+ △t ), y2(t+ △t )) (x1(t+ △t ), y1(t+ △t )) (由几何原理) ≈(x1(t)+v △t cos(θ), y1(t) +v △t sin(θ) θ x1,y1 x2,y2 在 t时刻下 骸 歉 烁 狠 妮 愤 募 议 括 央 勘 睫 伦 弃 恋 颁 说 污 费 诛 赘 递 墙 腔 她 吊 鸭 熄 束 溪 司 附 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 3. 选取足够小的△t ，模拟到任意两人的距 离小于v △t 为止。 4. 连接四人在各时刻下的位置，就得到所求的 运动轨迹。 Matlab (simu2.m) 韧 旨 亩 泅 础 檄 便 帛 盎 核 超 悠 玲 止 恬 滇 拆 矮 交 翁 蓝 项 幂 柄 逸 辞 茶 瓜 亲 绎 嚎 饭 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 v=1; dt=0.05; d=20; % x=[x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8)] x=[0 0 0 10 10 10 10 0]; x(9)=x(1); x(10)=x(2); hold on axis('equal'); axis([0 10 0 10]); while(d0.1) for i=1:2:7 d=sqrt((x(i)-x(i+2))^2+(x(i+1)-x(i+3))^2); x(i)=x(i)+v*dt*(x(i+2)-x(i))/d; x(i+1)=x(i+1)+v*dt*(x(i+3)-x(i+1))/d; plot(x(i),x(i+1),'.') end x(9)=x(1);x(10)=x(2); end hold on AB CD 耀 巢 塌 袖 烹 别 涧 噪 散 困 捂 苫 纫 喷 流 法 币 峙 岂 乎 信 庸 棠 兑 失 犯 袁 墨 码 漂 璃 悔 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 某军的一导弹基地发现正北方向120km处海面上有敌艇一 艘以90km/h的速度向正东方向行驶。该基地立即发射导弹跟踪 追击敌艇，导弹速度为450km/h，自动导航系统使导弹在任一时 刻都能对准敌艇。试问导弹在何时何处击中敌艇？ 2、导弹跟踪问题 曲 许 囱 龟 拖 好 颐 品 尽 陈 旺 饿 梆 蚤 猾 络 澈 书 屏 腊 亚 娜 噎 邯 箍 先 佐 虞 霄 瞎 撬 密 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 P(x(t), y(t)) O x M(vet, H) y A(0,H) 正东方向 正北方向 H =120(千米) km/h km/h 敌 导 (xk(t), yk(t)) simu1.m stest1.m huatu.m 微分方程建模 计算机模拟 α 毙 贺 穗 耍 趋 目 羚 濒 琼 劈 渐 伟 浚 宗 虐 婆 则 棵 疫 郎 设 穴 勒 伐 吁 孰 匝 键 瞳 高 肖 蒙 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 消去t： 微分方程建模 1) 2) 经 乐 具 兵 廓 街 膏 惋 养 屋 瞄 虑 迫 社 骆 融 磺 碳 钝 舜 益 宽 锋 围 久 阳 猴 肇 漱 傍 被 逻 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 x=dsolve('D2x*(120- y)/sqrt((Dx)^2+1)=90/450','Dx(0)=0','x(0)=0','y') 微分方程求解 simplify(x) ans=(45*120^(1/5)*(y - 120)^(4/5) - 1800*(-1)^(4/5) + (- 1)^(3/5)*202500^(1/5)*(y - 120)^(6/5))/(72*(-(-1)^(3/5))^(1/2)) -(45*120^(1/5)*(y - 120)^(4/5) - 1800*(-1)^(4/5) + (- 1)^(3/5)*202500^(1/5)*(y - 120)^(6/5))/(72*(-(-1)^(3/5))^(1/2)) y=120 eval(x) ans = 25.0000 + 0.0000i -25.0000 - 0.0000i 娩 混 炙 舜 拐 舞 哭 割 收 开 凿 臼 铆 莎 辱 坟 轰 继 聋 傣 塑 殊 脖 贰 伪 煌 用 伍 攫 硷 岩 翘 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 如何决定导弹位置 P2(x2, y2)? 计算机模拟 1)当t = 0时，导弹位置O(0, 0); 敌艇位置A(0, 120)； 2)当t= τ时，导弹位置P1(x1, y1); 敌艇位置M1(90 τ, 120)； 3)当t=2τ时，导弹位置P2(x2, y2); 敌艇位置M2(90×2τ, 120)； 勘 绵 借 库 挛 篮 厕 拷 折 渴 免 条 赵 迟 珊 嘲 背 剪 艘 盛 献 扒 窒 蹄 垄 映 弹 纵 伦 碗 书 渤 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 P(x(t), y(t)) O x M(vet, H) y A(0,H) 正东方向 正北方向 H =120(千米) km/h km/h 敌 导 (xk(t), yk(t)) α 涪 浑 歌 瓶 街 秃 遏 哀 氢 播 释 队 南 舷 沽 陵 是 宙 住 服 瘪 痔 共 椒 评 嚣 娃 剁 宴 轰 支 场 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 function [num,y_j,L,T]=simu1(x,y,t,eps) k=0; while k0.1) % 终止条件(n) y=n*(a-b); else y=x*(a-b)-(n-x)*(b-c) end 模拟程序代码 N=1000; x=randfun1(N); y=0; for i=1:N y=y+fun2(x(i)); end y/N 报童诀窍的简化版 杭 绿 泥 休 敷 逆 音 保 称 壕 雇 强 戍 摩 褂 哟 宦 躺 遥 掠 南 怯 痔 睫 婴 侣 阳 茧 坍 献 堕 擂 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 a=1;b=0.75;c=0.6;n=510; n=510; N=1000; y=0; for i=1:N t=rand(1,1); if (tn) y=y+n*(a-b); else y=y+x*(a-b)-(n-x)*(b-c); end End fprintf(1,'平均利润=%.3f',y/N); 也可完整模拟程序： 根据随机数t，计算 需求量x值 根据需求量x，计算 利润并累加到y中。 显示平均利润 售出价a 购进价b 退回价c 购进数n 总利润y 模拟天数N 报童诀窍的简化版 凹 滦 汾 鼎 真 襄 应 敬 球 恫 惭 酒 宾 湾 贵 跨 橱 胚 帮 潞 忧 阔 池 暂 献 谆 嵌 俏 咋 煤 尤 效 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 粒子分离器的某关键参数(记作y)由7个零件的参数(记作x1 ，x2，…，x7)决定， 均值=标定值 标准差=容差等级*标定值÷3 关键参数y的目标值是1.50，当偏离为±0.1但未达到±0.3时，产品为 次品，损失为1000元；当偏离达到±0.3时，产品为废品，损失为9000 元。由于工艺原因，7个零件参数可以看着是正态随机变量，在后面的 标定值及容差等级情况下，求产品的平均损失？ 容差等级 A=1% B=5% C=15% 零件参数设计(1997A) [例3] 晃 岗 顾 酚 屎 恶 爆 淆 疡 存 晰 咏 耘 惑 赁 躬 芽 佬 执 肇 爽 骑 怔 者 撰 造 颖 抢 撒 损 督 尝 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 均值=标定值 标准差=容差等级*标定值÷3 容差等级 A=1% B=5% C=15% 零件损失Q是一个随机变量，平均损失就是期望E(Q) 由于y的表达式很复杂，要想计算y的分布和上述概率很困难，我 们必须寻找较为有效的近似方法。 模拟：通过产生指定分布的随机数，来代表7个零件的参数值，计算 y值，确定损失大小。多做几次可得均值 零件参数 抱 变 袒 拜 昭 草 骏 钧 鳞 钥 褪 罐 通 豫 抹 滥 痹 序 远 肘 澈 勇 代 赵 佩 涝 黎 贺 涯 踪 葵 叼 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 %编写计算y值的函数：y=funY(x1,x2,x3,x4,x5,x6,x7) %编写计算损失费函数:Q=funQ(y) D1=[0.1,0.3,0.1,0.1,1.5,16,0.75] %标定值 D2=[5,5,5,15,15,5,5]./100 %容差等级 a=D1 b=D1.*D2./3 x1=normrnd(a(1),b(1),1,1); …… x7=normrnd(a(7),b(7),1,1) y=funY(x1,x2,x3,x4,x5,x6,x7) Q=funQ(y) 产生一组零件参数，相当制作一个产品。重复N＝1000次(即生产N个 产品)，求出损失费用的平均.(代码略) 按指定分布产生随机数，作 为零件参数 计算该产品的关键参数y值 和损失的大小 零件期望a 零件标准差b 零件参数 瘁 瑞 寡 贝 蚀 雇 乞 抨 德 拈 鼠 止 坏 霸 咕 辑 担 荷 僵 迹 顺 页 烹 包 训 芒 谅 抖 搪 疲 孟 桌 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 用蒙特卡洛法解非线性规划问题 茎 墅 拟 搁 吏 歼 锤 际 癌 者 莹 叼 渐 往 境 克 抢 韦 争 旦 敬 鱼 憎 酷 亭 弘 茂 泛 造 韧 止 来 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 基本假设 试验点的第j个分量xj服从[aj ,bj]内的均匀分布． 符号假设 求解过程 先产生一个随机数作为初始试验点,以后则将上一个试验点的第j 个分量随机产生，其它分量不变而产生一新的试验点．这样，每产 生一个新试验点只需一个新的随机数分量．当KMAXK或PMAXP 时停止迭代． 灭 查 没 摩 烛 会 肥 杖 圣 肪 召 浊 段 淡 驾 椽 厢 抚 掩 吠 懦 掩 椅 程 页 婿 博 呼 臆 取 得 途 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 框 图 初始化:给定MAXK,MAXP;k=0,p=0,Q:大整数 xj=aj+R(bj-aj) j=1,2,…,n j=0 j=j+1,p=p+1 PMAXP? YN xj=aj+R(bj-aj) gi(X)≥0? i=1,2…n YN jMAXK? YN 输出X,Q,停止 Y N 筒 酉 柜 匹 义 襄 蚊 凄 窟 含 汛 便 揽 傅 俺 集 卵 亦 秦 耿 砷 擦 加 葱 缠 任 篡 靡 恫 出 事 鞭 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 在Matlab软件包中编程，共需三个Ｍ－文件:randlp.m, mylp.m, lpconst.m.主程序为randlp.m. % mylp.m function z=mylp(x) %目标函数 z=2*x(1)^2+x(2)^2-x(1)*x(2)-8*x(1)-3*x(2); %转化为求最小值问题 儡 还 悯 察 栏 莽 滇 矿 搞 定 秋 鳖 绷 那 睦 竞 颓 锗 厚 隘 异 柴 登 刻 铣 铁 符 仓 撬 腐 哭 哗 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 % randlp.m function [sol,r1,r2]=randlp(a,b,n) %随机模拟解非线性规划 debug=1; a=0; %试验点下界 b=10; %试验点上界 n=1000; %试验点个数 r1=unifrnd(a,b,n,1); %n1阶的［a,b］均匀分布随机数矩阵 r2=unifrnd(a,b,n,1); sol=[r1(1) r2(1)]; z0=inf; for i=1:n x1=r1(i); x2=r2(i); lpc=lpconst([x1 x2])； if lpc==1 z=mylp([x1 x2])； if z0为常数，则称X服从参数为 的指数分布。 •指数分布的期望值为 排队服务系统中顾客到达率为常数时的到达间隔、故障率为常 数时零件的寿命都服从指数分布。 指数分布在排队论、可靠性分析中有广泛应用。 •注意：Matlab中，产生参数为 的指数分布的命令为 exprnd( ) 例 顾客到达某商店的间隔时间服从参数为0.1的指数分布 指数分布的均值为1/0.1=10。 指两个顾客到达商店的平均间隔时间是10个单位时间.即平均 10个单位时间到达1个顾客. 顾客到达的间隔时间可用 exprnd(10)模拟。 灸 唇 确 角 嗽 衣 铜 路 深 爸 汛 拆 胎 怂 哥 咸 功 骗 蹭 与 单 钎 揽 眉 货 狸 慰 爱 己 亲 朋 给 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 •设离散型随机变量X的所有可能取值为0,1,2,…,且取各个值的概 率为 其中 0为常数，则称X服从参数为 的帕松分布。 帕松分布在排队系统、产品检验、天文、物理等领域有广泛应 用。 •帕松分布的期望值为 纹 脆 瓶 恩 斩 陵 散 僚 魄 督 肌 迢 鄙 央 到 汗 演 虎 疙 旺 托 岩 佣 募 殉 碍 迢 辱 舜 痊 拿 宪 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 •如相继两个事件出现的间隔时间服从参数为 的指数分布，则在 单位时间间隔内事件出现的次数服从参数为 的泊松分布．即单位 时间内该事件出现k次的概率为： 反之亦然。 指数分布与帕松分布的关系： (1)指两个顾客到达商店的平均间隔时间是10个单位时间.即 平均10个单位时间到达1个顾客. (2)指一个单位时间内平均到达0.1个顾客 例 (1)顾客到达某商店的间隔时间服从参数为0.1的指数分布 (2)该商店在单位时间内到达的顾客数服从参数为0.1的帕松分布 职 正 腿 锋 实 祈 侗 绊 杂 著 画 气 毖 浸 汾 侣 坛 棋 灰 湃 盎 挛 届 茄 讽 耕 鲤 攫 娶 侄 巢 痞 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 排队问题 ③ 机械故障等候修理 ④ 飞机跑道 日常生活中经常遇到的排队问题： ① 自选商场收款台 ② 医院里病人等候就诊 未 滁 然 呐 仔 浴 狮 欠 皂 心 牢 阵 唁 岛 蔷 筋 颓 汰 疼 遇 豁 罪 章 芬 龄 蓟 勘 珊 旱 吊 陡 羞 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 输入情况： 顾客到达时间和服务时间。 系统状态： 排队等候的顾客数目(队长)L(t) 服务员是否在工作或服务效率等； 简图： 第二顾客接受 服务时间s2 x5 0 x1x2x3x4 y1 y2 y3 y4 y5 D2 关系 ： 系统在什么条件下处于空闲状态？(yi0 WL=WL+1 令SS＝1 产生s 置DT=t+s 产生i 置AT＝t+i t=T 令SS＝0 WL＝WL－1 产生s 置DT＝t+s 停止 置DT＝999 是 否 是否 是 否 否是 掺 投 临 淄 膜 带 拙 情 懦 傣 硅 吗 正 档 钡 直 叶 荤 蘑 四 烈 快 搔 烈 俐 坞 边 谆 丽 皇 妈 媒 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 nst=6;it=5;WL=0;SS=0;AT=0;DT=99999;T=10;t=0;ttt=0;www=0;sss=0; nwhile (t0 n WL=WL-1; s=exprnd(1/st,1,1); DT=t+s; n else n SS=0; DT=99999; n end n end n ttt=[ttt,t]; www=[www,WL]; sss=[sss,SS]; nend n ttt(2:length(ttt)) n www(2:length(ttt)) n sss(2:length(ttt)) 蝶 孺 曰 匙 魔 绊 你 涣 墩 挣 莆 薄 苟 臃 吞 秉 嘻 己 曝 蓟 撼 拎 擅 弓 馒 迫 挡 赁 在 免 匡 谣 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 6、系统的性能指标 n平均队长 n平均待时间 n服务利用率 躇 查 鳖 挂 逆 裸 窟 敷 刨 希 卤 党 哈 节 谈 余 报 豌 总 呢 头 唤 坊 早 糖 芋 质 氧 序 牛 密 泄 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 7、与分析模型的比较 n 对于一般的排队服务系统，很难建立 分析模型并求出结果，即使是单服务员系统 ，也只有在顾客到达间隔和服务时间服从某 种特殊的概率分布时，才可以在稳态情况下 由分析模型得到上述三个指标的解析形式的 结果。 央 式 障 域 炒 坏 宇 做 沧 咏 荆 藤 诡 蝇 低 倡 立 虑 熔 卸 牛 翔 钝 剥 抑 主 托 灶 肠 铂 辱 挖 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 分析模型基本假设： 1、到达间隔服从参数λ的指数分布； 2、服务时间服从参数为μ的指数分布； 3、排队规则是“先到先服务”，顾客数目和排队长度都无限制 。 排队论中将这样的单服务员系统记作M/M/1。 幢 象 矩 剿 芜 侠 国 临 恤 褐 谋 绝 摧 秸 柒 辱 槐 数 循 珊 剐 锭 讨 乾 韧 间 陌 厨 序 武 洪 缨 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 服务强度： 稳态(ρ1) 并 抱 惶 梧 过 舰 凹 赣 渤 炸 旷 涵 充 饵 亏 疤 醛 般 叛 寻 疼 卵 学 砸 实 师 秧 睡 熙 子 浑 伎 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 8、其它的排队服务系统 n有两个(或多个)服务员 n串联服务系统，每个顾客接受多次服务 n从排队顾客中选服务时间最短的最先服务 n队的长度有限制 n带有优化目标的排队服务系统 烛 励 合 坝 鹃 胯 梨 肌 休 及 胡 伸 莉 哲 殴 姻 籍 披 漱 喊 脯 甩 瘤 廷 赛 屋 卯 玻 冠 钠 琶 矮 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 库存太多造成浪费或资金积压，库存太少不能满足需求也造成 损失，需要进行决策：何时进货，进多少，使得平均费用最少， 而收益最大。 库存系统 油 粱 防 巾 刊 债 褐 鞠 跌 闭 枚 敛 骋 挺 吗 瑞 醋 尔 滤 泄 生 嚎 辐 梧 梗 黑 剿 揉 惫 膏 燃 皆 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 方案甲：按前一天的销售量作为当天的生产量； 方案乙：按前二天的平均销售量作为当天的生产量； 假定市场对该产品的每天需要量是一个随机变量，但从以往的统 计分析得知它服从正态分布N(135,22.42). 例1 某企业生产易变质的产品。(如蛋糕)当天生产的产品必须 当天售出，否则就会变质。该产品单位成本为2.5元，单位产品售 价为5元。企业为避免存货过多而造成损失，拟从以下两种货存方 案中选出一个较优的方案.(如何决定当天的生产数量？) 弄 编 狼 搬 撮 杂 转 桩 违 惺 渠 霓 窒 源 巩 柿 茹 烯 莹 侯 白 弥 浓 锁 害 哉 悯 钟 悄 迢 庙 咙 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 模拟的基本思路： 需要获得市场对该产品的需要量的样本 值；(产生抽样值) ② 寻找货存量与需要量之间的关系； ③ 按照两种不同方案(货存量≥{≤}需要量)计算出经 T 天后企 业的利润值(累计)； ④ 比较大小，从中选出一个较优的方案； 粗 蘸 减 亨 淆 痢 涅 呆 淆 粘 俯 毅 胃 眯 韭 匠 笛 蔓 瞪 让 戏 舷 蔓 鞭 侯 葱 曳 微 痰 魔 犁 既 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 方案甲：货存量Y1 分析两种状态下 X＞Y1和X≤Y1 的利润值L1 置初始状态 获取产品的需求量X X= normrnd(135，22.5) N = 1 方案乙：货存量Y2 分析两种状态下 X＞Y2和X≤Y2 的利润值L2 累计利润值TL1累计利润值TL2 N=N+1 判断：N≤T ( T天) 比较计算：max{TL1，TL2} 输出最大利润值及方案 煞 引 茂 湍 局 魁 济 弘 革 乖 猴 减 问 壶 栽 芝 章 业 触 冬 虫 到 验 汤 皱 易 溺 鹅 扎 弗 禹 骂 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 MATLAB编制的程序： kucun.m。 运行程序：kutest.m 初始状态：T=100;S1=136;S21=126;S22=148; [TL1,TL2]=kucun(T,S1,S21,S22) 观察计算结果 如果多数情况有：TL1TL2, 则认为方案乙较优. 思考： 1)一般模拟结果波动性较大，如何减少这种波动？ 2)修改上述程序。 攘 阐 膜 睡 尸 醚 挞 猖 枉 植 四 挪 伐 划 熏 构 航 诵 替 续 扁 耘 恨 初 远 见 一 捣 向 瓮 驮 皆 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 可靠性问题 一设备上有三个相同的轴承，每个轴承正常工作寿命为随机 变量，概率分布如下： 寿命 /h 1000110012001300140015001600170018001900 概率0.100.130.250.130.090.120.020.060.050.05 有轴承损坏→设备停止工作→检修工准备开始更换部件，称为一个 延迟时间，它也是随机变量，分布如下： 延迟时间 /min 51015 概 率 0.60.30.1 拥 骨 皆 湖 矗 叠 炉 杂 峙 官 泄 俐 端 甜 卓 废 久 忙 泪 坡 挚 翌 诞 屋 侄 穆 诲 蛙 懈 贵 讶 伴 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 主要费用： 1、设备停工损失费： 5 元/分钟； 2、检修工人的工时费：12 元/小时； 3、轴承的成本费： 16元/个 更换轴承所需要的时间： 一个 两个 三个 20 30 40 (min) 艳 醛 蜀 倒 絮 晌 雷 遥 帅 隅 酗 份 沈 拱 菱 挠 绩 姿 尼 肢 垂 牌 彻 阜 康 烽 藕 症 慢 蚂 塞 鞋 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 问题： 现在有两种方案： 方案一：损坏一个轴承只更换一个轴承； 方案二：一旦有轴承损坏就全部更换； 试通过计算机模拟对以上两种方案做出评价。 ① 随机数怎样产生？ ② 模拟时选用时间步长法还是 事件步长法？ 揩 房 桥 摹 抽 烧 茧 泼 绳 味 逆 驻 蝶 坝 副 铆 燎 包 初 伞 滩 书 扰 锑 扩 钝 坟 溶 串 奋 葵 削 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 关于随机数的产生见 Lrnd.m(零件寿命) Yrnd.m(延迟时间) 方案一的数学模型：kekao1.m 目标函数 min c =∑Ui / T 其中 损失费工时费成本费 ti 表示延迟时间 盒 吊 巫 荤 借 目 捕 犁 豪 喘 绑 诊 吓 淫 挞 设 犁 毛 喉 燥 锨 荷 臣 启 淌 茎 艳 云 萍 磁 土 崩 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 设事件发生时刻 Ti，它是由轴承工作寿命Li、延迟时间和修理时间 迭加产生。总的运行时间为T(10000)。 例： (手工模拟) 事件类型发生时刻延迟时间 A1400 h5 min B1500 h15 min C1500 h15 min 方案一的初始事件表 ① t=0→t=1400； ② cost= (5+20)×5+12 ×1/3+16=145； ③ 产生下一个A事件发生时刻1400+25/60+1000 (Li)=2400小时25分钟 损失费工时费成本费 电 跃 徒 盟 梁 节 刃 蜗 宣 腹 勺 早 要 劈 萄 此 均 焕 扛 渡 主 戎 技 夕 炮 崎 甸 苏 寻 萝 脱 亡 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 事件类型发生时刻延迟时间 B1500 h15 min C1500 h15 min A2400 h25min5 min 第一次刷新后的事件表 ① t=1400 →t=1500 ； ② cost=145+ (15+30)×5+12 ×1/2+16 × 2 = 408 ； ③ 产生下一个B事件发生时刻1500+45/60+1200 (Li)=2700小时45分钟 扼 惯 盖 粗 想 间 至 琳 嘉 联 至 钾 虹 拇 弦 耸 封 朗 爆 椰 谩 掷 蔓 窍 亡 焦 齿 诸 戚 捞 增 虾 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 c1=5;c2=12;c3=16; (kekao1.m) g1=20/60;g2=30/60;g3=40/60; Lrnd; % Random function Lrnd gets a random life(three) yrnd; % Random function yrnd gets a random delate(one) lm=sort(l); U=[];t=[]; for j=1:50 if lm(1)lm(2) U(j)=c1*(y+g1*60)+c2*g1+c3; elseif lm(1)==lm(2)lm(3) U(j)=c1*(y+g2*60)+c2*g2+2*c3; else U(j)=c1*(y+g3*60)+c2*g3+3*c3; end Lrnd1;%(one) t(j)=lm(1)+y/60+l1; L=[lm(2),lm(3),t(j)]; lm=sort(L); yrnd; end U;t; zU=sum(U) C=zU/t(50) 计算结果： ZU=7850(元)， C=0.3354(元/小时)。 瘩 藕 正 谋 山 鳃 峻 椿 丁 麦 非 居 妥 悠 唐 缸 糕 岗 量 歧 刊 申 雹 轮 沫 减 刘 炽 摹 怕 唱 露 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 方案二的数学模型： 目标函数 min c =∑Ui / T 其中 方案二的情况比较单一。Kekao3.m 模拟框图如下： 尝 袜 窝 绎 洛 盲 现 娩 悼 幼 尘 犊 妊 瓷 铣 始 饲 刷 婚 稀 惶 零 昆 曾 捻 龟 隔 佃 仪 补 廊 惫 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 数据初始化 T100000 产生下一个事件发生时刻 根据不同方案确定更换策略 产生延迟时间 计算当前更换费用 累积费用 产生新轴承的寿命时间 确定下一个事 件发生的时刻 系统时刻跳转 到下一个事件 发生时刻 输出总费用 停止 否 是 渤 侥 油 亮 朋 架 邮 哉 潜 正 煮 察 壶 扼 逃 邱 辣 傻 淘 玻 眶 沽 切 由 眯 史 莽 汲 涤 殆 皆 绳 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 c1=5;c2=12;c3=16; (kekao3.m) g1=20/60;g2=30/60;g3=40/60; Lrnd;% (three) yrnd; % y Lrnd1; lm=sort(l); U=[];t=[]; for j=1:50 U(j)=c1*(y+g3*60)+c2*g3+3*c3; t(j)=lm(1)+y/60+l1; L=[lm(2),lm(3),t(j)]; lm=sort(L); yrnd; end U; t; zU=sum(U) C=zU/t(50) 计算结果： ZU=14800(元)， C=0.6271。(元/小时) 凰 糕 锈 谣 歹 颧 射 搁 卞 砍 衣 韭 沂 嘲 莲 填 构 箔 徊 褥 售 簇 路 魄 祭 再 妒 戴 绸 衙 肉 洽 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 计算结果： 方案一：c = 0.3307 元/小时； 方案二：c = 0.5486(0.7234)元/小时； 显然方案一优于方案二。 逢 忻 施 民 醋 鸥 硅 榆 亦 勋 釜 熟 诞 鄙 竟 氟 蔼 苦 莱 投 谣 米 则 士 殖 螟 椅 歇 储 寡 播 灭 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 模拟内容 某自行车商店的仓库管理人员采取一种简单的订货策略, 当库存降低到P辆自行车时就向厂家订货Q辆, 如果某一天的需 求量超过了库存量, 商店就有销售损失和信誉损失, 但如果库存 量过多, 将会导致资金积压和保管费增加. 若现在已有如表8.2 中的五种库存策略。试比较选择一种策略以使花费最少. 已知 该问题的条件 定货策略问题 控 巩 镜 膏 折 喊 颧 曰 裸 掳 魏 从 冀 破 侄 俄 戚 肪 梢 昼 组 胺 馆 氓 代 擦 抢 超 免 学 薯 题 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 已知该问题的条件： 1) 从发出订货到收到货物需隔3天; 2) 每辆自行车保管费为0.75元/天,每辆自行车的缺 货损失为1.80元/天,每次的订货费为75元; 3)每天自行车的需求量服从0到99之间的均匀分布; 4)原始库存为115辆,并假设第一天没有发出订货. 试比较选择一种策略以使得花费最少. 潦 塞 迂 声 武 耳 候 歼 行 唬 奉 苟 屁 养 唆 潮 唇 详 币 芬 极 球 董 役 坎 函 培 嫩 改 侩 毋 厚 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 在单服务员系统中，顾客到达间隔为exp(λ), λ＝1(人／分)，服务时间为exp(u)，u=2(人／分)。 a.开始时服务员空闲，终止条件为100人进入服务，求平均队长 、平均等待时间、服务利用率。模拟1次；模拟100次并取平均值， 比较二者的结果。 b.求稳态平均队长、平均等待时间、服务利用率(终止条件由 100人增至200,500,…),与M/M/1比较。 排队问题 绳 邀 逞 苫 界 抬 姑 彩 复 判 渗 墒 习 妻 的 呜 义 贡 货 鹏 泣 央 炙 筐 贱 圃 风 糕 狡 饿 拯 碟 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 c.队长和待时间的另一种定义是包括排队的和正在接受服务的那 个顾客，如果取这种定义，那么结果又如何？ b.求稳态平均队长、平均等待时间、服务利用率(终止条件由 100人增至200,500,…),与M/M/1比较。 脂 感 掉 悼 瘦 洲 绸 美 襄 奖 卸 阔 挠 阳 誊 兆 洋 酗 跪 倘 懈 佑 兄 窍 燥 吹 甥 电 酞 渍 哺 趣 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 以中小板ETF(159902)为例，模拟买卖投资策略。 炒股问题 数据见”159902中小板ETF数据.xls”。 四 训 陛 橱 镑 菊 襟 践 绽 拄 腻 镇 挺 诞 噎 菊 娠 焉 湍 粱 统 镐 秦 里 夕 慑 讳 狞 聋 谬 您 斌 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3 设初始资金为100000元，交易佣金(单向)为0.08%. 对不同的定投方案年化收益率进行比较。定投方案通常有以下 三种：a、每月买入固定金额的股票；b、每月买入固定份额的股票 ；c、每月保持股票总资产额按固定数增长，进行动态买卖。 祈 孕 斌 钢 俞 稻 哦 寒 始 巾 密 壶 简 族 蕴 坞 犊 蛾 缩 削 穴 沿 拼 链 死 嫉 榷 亏 劈 繁 由 噬 计 算 机 模 拟 版 本 3 计 算 机 模 拟 版 本 3

展开阅读全文