用马尔可夫链的蒙特卡洛方法更新及预测连续油藏模拟模型_杨菁华

用马尔可夫链的蒙特卡洛方法更新及预测连续油藏模拟模型

编译:杨菁华　王奇(西安石油大学石油工程学院)审校:王青青(西安石油大学石油工程学院)

摘要　大多数油藏模拟研究都局限于静

态范畴,因时间和预算限制,未知参数明显

减少,这容易导致低估预测的不确定性或作

出不明智的决策。马尔可夫链的蒙特卡洛

(MCM C)方法已被用于静态研究以对预测

参数空间的不确定性进行严格的量化考察。

但这些方法在长期性和系列稳定性方面存在

局限。文中将M CMC应用于实时油藏建

模。较之传统方法,M CMC方法在某一特

殊时间点上应用更少的模型来实现合理的概

率预测,它也提供了一种随时校准不确定性

预测的机制。

关键词　连续油藏模拟　马尔可夫链　蒙特卡洛　概率预测　不确定性量化

DOI:10.3969/j.issn.1002-641X.2010.6.003

1　背景

1.1　不确定性量化技术

对产量预测中的不确定性进行量化完全依赖于历史拟合。进行历史拟合包含三个主要步骤:第一,根据描述油藏几何与流动特性的整套参数来定义油藏;第二,根据预估的概率分布对油藏未知参数进行赋值;第三,从预先的分配油藏模型可作为实际模拟例子,并在目标函数中量化模拟结果和观察数据的差异。许多方法用于寻求合理的模型,某些历史拟合的技术———梯度法技术,它们的目标就是最优化,即模拟结果最符合观测结果的模型就是所需要的模型。用梯度法进行的不确定性量化也许还不是最完美的方法,因为它只在目标函数的一个或多个极疏区的相邻区取样。M CMC方法已广泛应用于从一个复杂的分布函数中取样,特别是当不知道函数的精确形式的时候。

在油藏模型研究中,MCM C技术用于研究贝叶斯推断的后验分布。最后分布只不过产生于一套案例中油藏模型所用的样本。首先,一个随机模型取样于一个先验分布的样本,这是马尔可夫链的出发点。其次,第二个模型被随机选择但受限且关联于马尔可夫链中已存在的前一个模型。当这个马尔可夫链足够长后,可以用此链中的模型产生后验分布。生成未知分布的另一种相关方法叫作遗传算法(GA s),它有多个相关的应用软件。遗传算法是基于控制自然界遗传基因规则、类别涉及面广的一类最优化算法。

在遗传算法中的一个称为“繁殖”的过程中,特种油藏模型的子模型产生于先前运行的模型混合参数值。最后应用所有的子模型生成一种分布。较之GA s,M CMC方法在统计上更严格。M CM C也被当作是遗传算法的一个类型,因为下一个模型依赖于其前一个模型,这在GAs就是其母模型。

1.2　实时数据与综合最小方差线性递推估算滤波器

综合最小方差线性递推估算滤波器(卡尔曼滤波器)技术广泛应用于天气预报,近来应用于石油工程以进行概率预测。卡尔曼滤波器是一种蒙特卡洛手段,是一种通过连续模型校正而进行不确定性量化的技术。较之传统的历史拟合中进行同步数据同化,EnKF(综合卡尔曼滤波器)要顺序进行。EnKF从受静态数据(如岩心、测井及构造资料)约束的整体油藏模型开始运行。当数据可用时, EnKF用来自于综合模型的静态信息和模型预测数据更新每个模型的实现。因此,EnKF生成一套受生产历史约束的拟模型实现,而且从理论上讲,它应该尊重先前的静态或地质信息。然而,由于受依赖于样本数量或模型实现的所有统计测量中固有的杂波影响,对于小的整体规模,EnKF更新能导致地质不协调实现。因此,如果最终的模型实现尊重历史生产数据,那么模型可能就不符合先前的地质信息。在EnKF过程中使用大的整体规模可以减少运行的困难,但对于油田规模计算费用可能比较昂贵。EnKF仍是为达到历史拟合而进行积极研究的主题,并且该技术正致力于解决运行中的难题。

2　连续模拟过程

用连续MCM C方法进行历史拟合以及生成概率预测需要几个部分相结合:首先,不确定性油藏的参数以及先验分布必须是确定的。其次,既然将来的油藏动态的不确定性通常要通过一套油藏模型的模拟动态来评估,那么必须有在参数空间取样以

及生成油藏模型的方法。M CMC方法在这里被用来探究不确定性参数空间,并自动生成运行模拟的编码。针对每一示例模型,随着历史拟合的进行,预测也一并进行,同时运行结果也被存储起来。在连续的历史拟合过程中,油田现场的新数据也一直被加入目标函数中。被更新的目标函数致力于后来的模拟运行。最后,各个运行结果被合并入概率预测。

2.1　参数化

在作任何模拟之前都必须首先决定考虑哪些不确定性参数。一旦确定了相关的参数,就要指定先验分布(这里通常是连续的)以在数量上描述参数的不确定性。正确的分布类型通常是建立在油藏特征数据上的。后期的分布就在贝叶斯框架内定义,而先前的分布根据观测的动态数据作出修正,表示如下:

P(m|d obs)∝P(d obs|m)P(m)(1)式中　d obs———从油田现场上观测的动态数据,如

含水率;

m———不确定性参数;

P(m)———不确定性参数的先验概率分布;

P(d ob s|m)———与观测数据有关的似然函

数;

P(m|d obs)———后期的分布函数。

如果假设先验模型和数据误差服从正态分布,那么这里的后验分布P(m|d ob s)就成为如下形式(Ho w ard,2005):

P(m|d obs)∝ex p-1

2

{(m-μ)T C-1m(m-μ)+

[g(m)-d ob s]T C-1D[g(m)-d obs]}(2)式中　g(m)———符合不确定性参数的先验概率分

布m的模拟油藏响应(曲线);

C m———参数协方差;

C D———资料协方差。

2.2　目标函数

目标函数是从数量上评估一个独立模型再现油田现场观测数据的满意程度的函数。这一项被定义为如下后验分布函数的一部分:

O(m)=(m-μ)T C-1m(m-μ)+

[g(m)-d obs]T C-1D[g(m)-d obs](3) 在(3)式中,(m-μ)T C-1m(m-μ)是先验项; [g(m)-d obs]T C-1D[g(m)-d obs]是似然项;μ表示先验平均值。目标函数的结果就是先验信息与观测信息的结合,是一个在贝叶斯框架内后验分布结构的推论。

2.3　Metropolis-H asting M CMC算法

该MCM C方法的主要目标是构造一个平稳分布与后验分布相匹配的马尔科夫链。这种后验分布典型地定义在多维参数空间上并常有多路系统式。这里所使用的Metropolis-H asting M CMC方法(H asting s,1970)常被用来从复杂后验分布中采样,其随机步进M-H采样过程如下:

(1)从先验分布中随机采样一套参数,表示为m t1;

(2)从状态t i时刻到t i+1,m t i+1=m t i+σε,其中,ε是一个标准的正态随机变量;σ是一个比例因子;

(3)R=P(m

t

i+1|d t i+1obs)

P(m t i|d t i obs)

;

(4)从0～1间的均匀分布中随机取一个数y。如果y≤R,则接受链中的m t i+1;否则又将m t i放入链中;

(5)回到第(2)步。

2.4　连续数据的同化作用

在连续模拟过程中的不同时间点上,油田现场新数据就显得很有用了。由于通常会设想油田现场的信息越丰富则导致预测和不确定性评估的效果就越好,因此尽快在模拟过程中引入新数据是有利的。当更多的数据加入,则(3)式中将包括更多的观测及模拟数据点。这样,观测数据的限定不符合目标函数中相关项的状况就将随着时间的推移而改变。尽管持续变化的目标函数在统计上是不严格的,但可以假设在这种违规情况下能产生合理的概率预测。

2.5　概率预测

模拟过程的最后一步是将为单独的马尔可夫链样本所做的生产预测合并入概率预测中。在连续的MCMC过程中,通过使用足够多的最新采样模型可在任何时候生成概率预测。要么用少量的模型以使改变目标函数的影响最小,要么用大量的模型以使获取可代表后验分布的平稳系列分布的机会最大,这二者之间需要寻求一个可以接受的平衡。

3　PU N Q-S3油藏实例研究

PUNQ-S3人造油藏作为在生产预测中量化不确定性的测试实例已应用了很多次。PUN Q-S3油藏模型是一个基于实际油田的5层3相人造油藏(Bo s,1999)。该模型包括13300(19×28×25)个网格块,其中1761个网格块起作用。油藏被东部和南部的断层所限制并与西和北边相当厚的含水

层有关系。带井位的结构如图1所示,图中有一个小气顶在中央并以红色显示,6口生产井表示为黑色的圆点。石油行业和学术理论上的几个合作伙伴用不同的方法来检验历史拟合与不确定性量化技术。这个PUNQ 项目的目标就是确立一种方法———可以通过适当方式将油藏模型、油藏参数、油井观察结果的不确定性与生产预测的不确定性结合起来

。

图1　P U NQ 人造油藏结构图

用来确定不确定性参数及指定先验分布的过程

与传统上模拟研究中评估输入不确定性时所做的工作是一致的。在PUNQ -S3模型中,假设孔隙度呈正态分布,而渗透率呈对数正态分布。另外,研究中的不确定性参数不直接使用孔隙度和渗透率值而是孔隙度和渗透率的乘积,在运行模拟过程中,这些乘子应用于孔隙度和渗透率底图。这样做的效果与直接使用孔隙度和渗透率值的效果相同,但是这种方法简化了实施过程。对每一层使用统一的属性值绘制底图。基于Barker 2001年提供的平均孔隙度常量值,垂向和纵向渗透率通过(4)、(5)式(Gu and Oliver ,2004)计算得到。表1列出了底图中使用的平均孔隙度和渗透率值。

lo g 10(k h )=9.02 +0.77(4)k v =0.31k h +3.12

(5)

式中,k v 表示垂向渗透率;k h 表示横向渗透率。

用每层6个同源区对PUNQ -S3模型进行参数化。另外,通过分裂油藏的方式将各区视为相互独立的区域,这一点是基于实时地质描述的,指示油藏广泛分布着东南走向的优质夹层(Imperial Co l -lege ,2007)。这样,在先验分布中,5层6区3个

属性(孔隙度、垂向渗透率和横向渗透率)就产生总共90个未确定的乘子。孔隙度乘子被指定服从中值为1、标准偏差为0.3的正态分布(Barker et al ,2001)。基于(4)、(5)式,垂向和横向渗透率乘子被指定服从中值为1、标准偏差为1.35的对数正态分布。为防止出现极端和不实际的渗透率值,此乘子分布被冠以上不过4而下不过0的限制。同样,孔隙度范围在0～2.28。

表1　P U N Q -S3油藏每层的平均孔隙度、渗透率值

平均值第1层第2层第3层第4层第5层孔隙度

0.170.080.170.160.19水平渗透率/m D 43233432196654垂直渗透率/m D 137

13

137

64

205

　1mD =1.02×10-3μm 2

3.1　目标函数

由于在研究中确定的不确定性参数是孔隙度和渗透率,目标函数变为:

O (X )=(X -μ)T C -1x (X -μ)+[g (X )-d obs ]T C -1

D [g (X )-d obs ]

(6)

式中,X 代表不确定的乘子。

在实验过程中,在第5、6、7、8、9年初顺序加入历史资料。因此,当加入新数据时,目标函数

中的观测资料值的数量成倍数增加。

3.2　参数空间研究

在连续历史模拟过程中,模拟运行过程是与历时8.04年的观测资料进行拟合。进行的预测是16.5年的生产经历。PUNQ -S3模型连续模拟开始于第4.5年,持续到第9年底,并作出16.5年的预测。在第一个半年(从第4.5年到第5年),用MCMC 采样方法的模拟中有4500个模型作为样品,并且在其余的年份里每年样本模型数为9000,最后模拟样品总数是49500。所有这些都假设9000个模型的运行都发生在一年中的实际时间,即每小时运行一个模型,这不是今天石油行业中用于油田模拟的非典型的运行频率。图2显示模拟运行的累积次数与油藏生产时期的对比。包括在目标函数计算中的额外数据也导致目标函数更大且引起

目标函数曲线中太多的数值起落跳跃(图3)。3.3　预测

一个独立的预测是与每一次历史拟合一起运行。在第5、6、7、8、9年和第10年,利用在链中可用的样品模型在相应的时间进行概率预测。通过过去一年(或第5年中的半年)的模拟运行可以生成这些概率预测,这些概率预测的累积分布一起显示在图4中。

图2　

连续测试中模型运行数与油藏生产各时期的关系曲线

图3　目标函数与连续实例中M CM C

链模型的对应关系

图4　累积分布函数的比较

3.4　不确定性估计的校准

从图5中可以看出,前3个预测没能反映真实情形,在4.5～5年间偏离明显。因此,在前几年里不确定性明显淡化了,这可能归因于先验分布中对不确定性的低估,或是对观测资料误差的低估。即使没有大量的可用动态资料,人们也总是要去量化不确定性。为了增加不确定性,可以在先验分布或观测资料或两者中使用更大的标准差。为了克服

这些缺点,增大了先验乘子的标准差。渗透率乘子标准差从1.35增大到20,而孔隙度乘子标准差从0.3增大到0.5。伴随着标准差的增大,先前的后

验分布被重新建立,直到它们在本质上支持所有继发的后验分布

。

图5　连续实例预测与P U NG 油藏分布的预测值的对比

4　结论

MCMC 方法是一种应用于历史拟合与不确定性量化的强有力工具。通过对油藏的历史拟合及其整个生命周期的连续预测,随着时间的推移和观测

资料的增加,预测不确定性范围逐渐变小。较之MCMC 方法传统上应用于过去一时间上的模拟研究,连续M CMC 方法可用很少的模型在特别的时间点上进行合理的概率预测,这是因为连续方法得益于最新观测资料同化前的模型运行。在油藏的长期生产中持续建模可考虑以下方面:在参数化过程中获得更多的不确定性参数,更多的样本以及由此造成的参数空间大得多的片断取样。这应该最终导致更可靠的概率预测。

连续的模拟方法也提供一种从始至终的不确定性估计校准机制。如果观测到后验分布除了变窄以外始终发生的跳跃起伏,那么可以放大先验分布中的标准差或者增加观测资料以增加不确定性。应该采取调整措施直到继发的后验分布与所有先前的分布相一致。

资料来源于美国《SPE 119197》(收稿日期　2009-03-23)

蒙特卡洛(Monte Carlo)模拟法

当科学家们使用计算机来试图预测复杂的趋势和事件时, 他们通常应用一类需要长串的随机数的复杂计算。设计这种用来预测复杂趋势和事件的数字模型越来越依赖于一种称为蒙特卡罗模似的统计手段, 而这种模拟进一步又要取决于可靠的无穷尽的随机数目来源。 蒙特卡罗模拟因摩纳哥著名的赌场而得名。它能够帮助人们从数学上表述物理、化学、工程、经济学以及环境动力学中一些非常复杂的相互作用。数学家们称这种表述为“模式”, 而当一种模式足够精确时, 他能产生与实际操作中对同一条件相同的反应。但蒙特卡罗模拟有一个危险的缺陷: 如果必须输入一个模式中的随机数并不像设想的那样是随机数, 而却构成一些微妙的非随机模式, 那么整个的模拟(及其预测结果)都可能是错的。 最近, 由美国佐治亚大学的费伦博格博士作出的一分报告证明了最普遍用以产生随机数串 的计算机程序中有5个在用于一个简单的模拟磁性晶体中原子行为的数学模型时出现错误。科学家们发现, 出现这些错误的根源在于这5个程序产生的数串其实并不随机, 它们实际上隐藏了一些相互关系和样式, 这一点只是在这种微小的非随机性歪曲了晶体模型的已知特 性时才表露出来。贝尔实验室的里德博士告诫人们记住伟大的诺伊曼的忠告:“任何人如果相信计算机能够产生出真正的随机的数序组都是疯子。” 蒙特卡罗方法（MC） 蒙特卡罗（Monte Carlo）方法： 蒙特卡罗（Monte Carlo）方法，又称随机抽样或统计试验方法，属于计算数学的一个分支，它是在本世纪四十年代中期为了适应当时原子能事业的发展而发展起来的。传统的经验方法由于不能逼近真实的物理过程，很难得到满意的结果，而蒙特卡罗方法由于能够真实地模拟实际物理过程，故解决问题与实际非常符合，可以得到很圆满的结果。这也是我们采用该方法的原因。 蒙特卡罗方法的基本原理及思想如下： 当所要求解的问题是某种事件出现的概率，或者是某个随机变量的期望值时，它们可以通过某种“试验”的方法，得到这种事件出现的频率，或者这个随机变数的平均值，并用它们作为问题的解。这就是蒙特卡罗方法的基本思想。蒙特卡罗方法通过抓住事物运动的几何数量和几何特征，利用数学方法来加以模拟，即进行一种数字模拟实验。它是以一个概率模型为基础，按照这个模型所描绘的过程，通过模拟实验的结果，作为问题的近似解。可以把蒙特卡罗解题归结为三个主要步骤：构造或描述概率过程；实现从已知概率分布抽样；建立各种估计量。 蒙特卡罗解题三个主要步骤： 构造或描述概率过程： 对于本身就具有随机性质的问题，如粒子输运问题，主要是正确描述和模拟这个概率过程，对于本来不是随机性质的确定性问题，比如计算定积分，就必须事先构造一个人为的概率过程，它的某些参量正好是所要求问题的解。即要将不具有随机性质的问题转化为随机性质的问题。 实现从已知概率分布抽样： 构造了概率模型以后，由于各种概率模型都可以看作是由各种各样的概率分布构成的，因此产生已知概率分布的随机变量（或随机向量），就成为实现蒙特卡罗方法模拟实验的基本手段，这也是蒙特卡罗方法被称为随机抽样的原因。最简单、最基本、最重要的一个概率分布是(0,1)上的均匀分布（或称矩形分布）。随机数就是具有这种均匀分布的随机变量。随机数序列就是具有这种分布的总体的一个简单子样，也就是一个具有这种分布的相互独立的随机变数序列。产生随机数的问题，就是从这个分布的抽样问题。在计算机上，可以用物理方法产生随机数，但价格昂贵，不能重复，使用不便。另一种方法是用数学递推公式产生。这样

基于马尔可夫排队模型的行程时间预测方法

第34卷　第4期吉林大学学报(工学版)　Vol.34　No.4 2004年10月Journal of Jilin University(Engineering and Technology Edition)　Oct.2004 文章编号:1671-5497(2004)04-0671-04 基于马尔可夫排队模型的行程时间预测方法 杨志宏1,杨兆升2,于德新2,陈　林2 (1.宝路集团,吉林长春　130022;2.吉林大学交通学院,吉林长春　130022) 摘　要:针对城市交通流诱导系统(U TF GS)亟待解决的综合路段行程时间预测这一关键问题,利用马尔可夫排队模型给出了车辆路段(含信号交叉口)实时行程时间预测的基本公式,并结合实际工程项目对公式中的一些参数进行了简化,提高了模型的实用性。人工调查数据验证表明该模型具有较高的精度。同时给出了相对误差图。 关键词:交通运输工程;城市交通流诱导系统(U TF GS);马尔可夫排队模型;排队等待时间;实时动态行程时间 中图分类号:U491.2 文献标识码:A T ravel time prediction method based on Malcov queuing model YAN G Zhihong1,YAN G Zhaosheng2,YU Dexin2,CHEN Lin2 (1.China B aolu Com pany,Changchun130022,China;2.College of T ransportation,Jilin U niversity,Changchun 130022,China) Abstract:Aiming at the key problem of synthetic Link travel time prediction in Urban Traffic Flow Guidance System(U TF GS).A Vehicle link travel time prediction algorithm based on Malcov Queuing model was presented.With a quantity of traffic measurement data,some model parameters were simplized and confirmed,thus getting a high precision and also making the model more become applicable. K ey w ords:traffic engineering;U TF GS;Malcov queuing model;queuing wait time;real2time dynamic travel time 0　引　言 交通流诱导以交通流预测和实时动态交通分配(D TA)为基础,应用现代通信技术、电子技术、计算机技术等为路网上的出行者提供必要的交通信息,为其指出当前的最佳行驶路线,从而避免盲目出行造成的交通阻塞,到达路网畅通、高效运行的目的[1,2]。交通流诱导的方式一般分为路边显示板式和车内显示屏式两种。前者主要适用于高速公路以及城市路网集体车辆诱导,后者主要适用于城市路网中的个体车辆诱导[2]。 为了准确、快速地给出路网的最佳行驶路线,需要估计路网中各路段的行程时间。路网中的路段均指含一个相邻的下游交叉口(有信号灯控制)的路段。当车辆进入路段后,其行程时间随交通流量的变 收稿日期:2004205219. 基金项目:“十五”国家智能交通重大科技攻关项目(2002BA404A22B). 作者简介:杨志宏(1971-),男,工程师.E2mail:yangzhihong0527@https://www.wendangku.net/doc/ee14974782.html, 通讯联系人:杨兆升(1938-),男,教授,博士生导师.E2mail:yangzs@https://www.wendangku.net/doc/ee14974782.html,

马尔可夫链模型简介

马尔可夫链模型简介 设考察对象为一系统,若该系统在某一时刻可能出现的事件集合为,}{N N E E E E E E ??????,2,1,2,1,两两互斥,则陈i E 为状态。N i ???=,2,1。称该系统从一种状态i E 变化到另一状态j E 的过程称为状态转移，并把整个系统不断实现状态转移的过程称为马尔可夫过程。 定义1 具有下列两个性质的马尔可夫过程称为马尔可夫链： （1）无后效性，即系统的第n 次实验结果出现的状态，只与第1-n 次有关，而与它以前所处的状态无关； （2）具有稳定性，该过程逐渐趋于稳定状态，而与初始状态无关。 定义2 向量),,,(21n u u u u ???= 成为概率向量，如果u 满足： ?? ???=???=≥∑=n j j j u n j u 11,,2,10 定义3 如果方阵P 的每行都为概率向量，则称此方阵为概率矩阵。 如果矩阵A 和B 皆为概率矩阵，则AB ，k A ，k B 也都是概率矩阵（k 为正整数）。 定义4 系统由状态i E 经过一次转移到状态j E 的概率记为ij P ，称矩阵 ????????????????????????=32 12222111211N N N N N P P P P P P P P P P 为一次（或一步）转移矩阵。 转移矩阵必为概率矩阵，且具有以下两个性质： 1、P P P k k )1()(-=； 2、k k P P =)(

其中)(k P 为k 次转移矩阵。 定义5 对概率矩阵P ，若幂次方)(m P 的所有元素皆为正数，则矩阵P 称为正规概率矩阵。（此处2≥m ） 定理1 正规概率矩阵P 的幂次方序列P ，2P ，3P ，…趋近于某一方阵T ，T 的每一行均为同一概率向量t ，且满足t tP = 。 马尔可夫链模型如下： 设系统在0=k 时所处的初始状态 ),,() 0()0(2)0(1)0(N S S S S ???=为已知，经过k 次转移后的状态向量 ),,()()(2)(1)(k N k k k S S S S ???=),2,1(???=k ，则 ??????? ?????? ?????????????=NN N N N N k P P P P P P P P P S S 212222111211)0() ( 此式即为马尔可夫链预测模型。 由上式可以看出，系统在经过k 次转后所处的状态)(k S 取决与它的初始状态)0(S 和转移矩阵P 。 马尔可夫引例 例1：市场占有率预测 设有甲、乙、丙三家企业，生产同一种产品，共同供应1000家用户，各用户在各企业间自由选购，但不超出这三家企业，也无新的用户，假定在10月末经过市场调查得知，甲，乙，丙三家企业拥有的客户分别是：250户，300户，450户，而11月份用户可能的流动情况如下表所示：

蒙特卡洛模拟方法作业及答案(附程序)

蒙特卡洛习题 1.利用蒙特卡洛计算数值积分 () ()() 1280ln 1tan x x x xe dx +++? clear all ;clc;close all ; n=1000; count=0; x=0:0.01:1; y=log((1+x).^2+(tan(x).^8)+x.*exp(x)); plot(x,y,'linewidth',2) hold on for i=1:n x1=rand; y1=rand*y(end); plot(x1,y1,'g*') pause(0.01) if y1

2.分别用理论计算和计算机模拟计算，求连续掷两颗骰子，点数之和大于6且第一次掷出的点数大于第二次掷出点数的概率。 clear all;clc;close all; count=0; n=100000; for i=1:n x=floor(rand*6+1); y=ceil(rand*6); if x+y>6&&x>y count=count+1; end end P=count/n 3.

clear all;clc;close all; count=0; n=2000; ezplot('x^2/9+y^2/36=1'); hold on ezplot('x^2/36+y^2=1'); hold on ezplot('(x-2)^2+(y+1)^2=9') for i=1:n x=rand*12-6; y=rand*12-6; plot(x,y,'gh','linewidth',2) pause(0.01) if x^2/9+y^2/36<1&&x^2/36+y^2<1&&(x-2)^2+(y+1)^2<9

蒙特卡洛方法

蒙特卡洛方法 1、蒙特卡洛方法的由来 蒙特卡罗分析法（Monte Carlo method），又称为统计模拟法，是一种采用随机抽样（Random Sampling）统计来估算结果的计算方法。由于计算结果的精确度很大程度上取决于抽取样本的数量，一般需要大量的样本数据，因此在没有计算机的时代并没有受到重视。 第二次世界大战时期，美国曼哈顿原子弹计划的主要科学家之一，匈牙利美藉数学家约翰·冯·诺伊曼（现代电子计算机创始人之一）在研究物质裂变时中子扩散的实验中采用了随机抽样统计的手法，因为当时随机数的想法来自掷色子及轮盘等赌博用具，因此他采用摩洛哥著名赌城蒙特卡罗来命名这种计算方法，为这种算法增加了一层神秘色彩。 蒙特卡罗方法提出的初衷是用于物理数值模拟问题, 后来随着计算机的快速发展, 这一方法很快在函数值极小化、计算几何、组合计数等方面得到应用, 于是它作为一种独立的方法被提出来, 并发展成为一门新兴的计算科学, 属于计算数学的一个分支。如今MC方法已是求解科学、工程和科学技术领域大量应用问题的常用数值方法。 2、蒙特卡洛方法的核心—随机数 蒙特卡洛方法的基本理论就是通过对大量的随机数样本进行统计分析，从而得到我们所需要的变量。因此蒙特卡洛方法的核心就是随机数，只有样本中的随机数具有随机性，所得到的变量值才具有可信性和科学性。

在连续型随机变量的分布中, 最基本的分布是[0, 1]区间上的均匀分布, 也称单位均匀分布。由该分布抽取的简单子样ξ1，ξ2ξ3……称为随机数序列, 其中每一个体称为随机数, 有时称为标准随机数或真随机数, 独立性和均匀性是其必备的两个特点。真随机数是数学上的抽象, 真随机数序列是不可预计的, 因而也不可能重复产生两个相同的真随机数序列。真随机数只能用某些随机物理过程来产生, 如放射性衰变、电子设备的热噪音、宇宙射线的触发时间等。 实际使用的随机数通常都是采用某些数学公式产生的，称为伪随机数。真随机数只是一种数学的理想化概念，实际中我们所接触到的和使用的都是伪随机数。要把伪随机数当成真随机数来使用, 必须要通过随机数的一系列的统计检验。 无论伪随机数用什么方法产生,它的局限性都在于这些随机数总是一个有限长的循环集合, 而且序列偏差的上确界达到最大值。所以若能产生低偏差的确定性序列是很有用的,产生的序列应该具有这样的性质, 即任意长的子序列都能均匀地填充函数空间。 人们已经产生了若干种满足这个要求的序列,如Halton序列、Faure序列、Sobol序列和Niederreiter序列等。称这些序列为拟随机数序列。伪随机序列是为了模拟随机性, 而拟随机序列更致力于均匀性。 3、蒙特卡洛方法的原理 当问题可以抽象为某个确定的数学问题时，应当首先建立一个恰当的概率模型，即确定某个随机事件A或随机变量X，使得待求的解等

蒙特卡罗方法简介

第三章蒙特卡罗方法简介 3.1 Monte Carlo方法简介 Monte Carlo方法是诺斯阿拉莫斯实验室在总结其二战期间工作（曼哈顿计划）的基础上提出来的。Monte Carlo的发明，主要归功于Enrico Fermi、Von Neumann和Stanislaw Ulam等。自二战以来，Monte Carlo方法由于其在解决粒子输运问题上特有的优势而得到了迅速发展，并在核物理、辐射物理、数学、电子学等方面得到了广泛的应用。Monte Carlo的基本思想就是基于随机数选择的统计抽样，这和赌博中掷色子很类似，故取名Monte Carlo。 Monte Carlo方法非常适于解决复杂的三维问题，对于不能用确定性方法解决的问题尤其有用，可以用来模拟核子与物质的相互作用。在粒子输运中，Monte Carlo技术就是跟踪来自源的每个粒子，从粒子产生开始，直到其消亡（吸收或逃逸等）。在跟踪过程中，利用有关传输数据经随机抽样来决定粒子每一步的结果[6]。 3.2 Monte Carlo发展历程 MCNP程序全名为Monte Carlo Neutron and Photon Transport Code (蒙特卡罗中子-光子输运程序)。Monte Carlo模拟程序是在1940年美国实施“发展核武器计划”时，由洛斯阿拉莫斯实验室（LANL）提出的，为其所投入的研究、发展、程序编写及参数制作超过了500人年。1950年Monte Carlo方法的机器语言出现, 1963年通用性的Monte Carlo方法语言推出，在此基础上，20世纪70年代中期由中子程序和光子程序合并,形成了最初的MCNP程序。自那时起，每2—3年MCNP更新一次, 版本不断发展，功能不断增加，适应面也越来越广。已知的MCNP程序研制版本的更新时间表如下：MCNP-3：1983年写成，为标准的FORTRAN-77版本，截面采用ENDF /B2III。 MCNP-3A：1986年写成，加进了多种标准源，截面采用ENDF /B2I V[20]。

数学建模之马尔可夫预测

马尔可夫预测 马尔可夫过程是一种常见的比较简单的随机过程。该过程是研究一个系统的 状况及其转移的理论。它通过对不同状态的初始概率以及状态之间的转移概率的研究,来确定状态的变化趋势,从而达到对未来进行预测的目的。 三大特点： （1）无后效性 一事物的将来是什么状态，其概率有多大，只取决于该事物现在所处的状态如何，而与以前的状态无关。也就是说，事物第n 期的状态，只与第n 期内的变化和第n-1期状态有关,而与第n-1期以前的状态无关。 （2）遍历性 不管事物现在所处的状态如何，在较长的时间内马尔可夫过程逐渐趋于稳定状态，而与初始状态无关。 （3）过程的随机性。 该系统内部从一个状态转移到另一个状态是，转变的可能性由系统内部的原先历史情况的概率值表示。 1.模型的应用， ①水文预测， ②气象预测， ③地震预测， ④基金投资绩效评估的实证分析， ⑤混合动力车工作情况预测， ⑥产品的市场占有情况预测。 2.步骤 ①确定系统状态 有的系统状态很确定。如：机床工作的状态可划分为正常和故障，动物繁殖后代可以划分为雄性和雌性两种状态等。但很多预测中，状态需要人为确定。如：根据某种产品的市场销售量划分成滞销、正常、畅销等状态。这些状态的划分是依据不同产品、生产能力的大小以及企业的经营策略来确定的，一般没有什么统一的标准。在天气预报中，可以把降水量划分为旱、正常和涝等状态。 ②计算初始概率()0i S 用i M 表示实验中状态i E 出现的总次数，则初始概率为 ()()0 1 1,2,i i i n i i M S F i n M =≈= =∑L ③计算一步转移概率矩阵

令由状态i E 转移到状态j E 的概率为()|ij j i P P E E =，则得到一步转移概率矩阵为： 1112121 2221 2n n n n nn p p p p p p P p p p ??????=??????L L M M M M L ④计算K 步转移概率矩阵 若系统的状态经过了多次转移，则就要计算K 步转移概率与K 步转移概率矩阵。 K 步转移概率矩阵为： 11121212221 2()k n n k n n nn p p p p p p P k p p p p ??????==??????L L M M M M L ⑤预测及分析 根据转移概率矩阵对系统未来所处状态进行预测，即： () ()111210212221 2K n K n n n nn p p p p p p S S p p p ??????=??????L L M M M M L 例题： 设某企业生产洗涤剂为A 型,市场除A 型外,还有B 型、C 型两种。为了生产经营管理上的需要,某企业要了解本厂生产的A 型洗涤剂在未来三年的市场占有倩况。为此,进行了两项工作,一是进行市场调查,二是利用模型进行预测。 市场调查首先全面了解各型洗涤剂在市场占有情况。年终调查结果:市场洗涤剂目前总容量为100万件,其中A 型占40万,B 型和C 型各占30万。 再者,要调杏顾客购买各型洗涤剂的变动情况。调查发现去年购买A 型产品的顾客,今年仍购A 型产品24万件,转购B 型和C 型产品备占8万件,去年购买B 型产品顾客,今年仍购B 型产品9万件,转购A 型15万件,转购C 型6万件,去年购买C 型产品的顾客,今年仍购C 型产品9万件,转购A 型15万件,转购B 型6万件。计算各型产品保留和转购变动率。 模型的建立： ①计算初始概率 用i M 表示i E 型产品出现的总次数，则初始概率为 ()()0 1 1,2,i i i n i i M S F i n M =≈= =∑L （1） ②计算各类产品保留和转购变动率

蒙特卡洛模拟原理及步骤

二、蒙特卡洛模拟原理及步骤 （一）蒙特卡洛模拟原理：经济生活中存在大量的不确定与风险问题，很多确定性问题实际上是不确定与风险型问题的特例与简化，财务管理、管理会计中同样也存在大量的不确定与风险型问题，由于该问题比较复杂，一般教材对此问题涉及较少，但利用蒙特卡洛模拟可以揭示不确定与风险型问题的统计规律，还原一个真实的经济与管理客观面貌。 与常用确定性的数值计算方法不同，蒙特卡洛模拟是用来解决工程和经济中的非确定性问题，通过成千上万次的模拟，涵盖相应的可能概率分布空间，从而获得一定概率下的不同数据和频度分布，通过对大量样本值的统计分析，得到满足一定精度的结果，因此蒙特卡洛模拟是进行不确定与风险型问题的有力武器。 1、由于蒙特卡洛模拟是以实验为基础的，因此可以成为财务人员进行风险分析的“实验库”，获得大量有关财务风险等方面的信息，弥补确定型分析手段的不足，避免对不确定与风险决策问题的误导； 2、财务管理、管理会计中存在大量的不确定与风险型问题，目前大多数教材很少涉及这类问题，通过蒙特卡洛模拟，可以对其进行有效分析，解决常用决策方法所无法解决的难题，更加全面深入地分析不确定与风险型问题。 （二）蒙特卡洛模拟步骤以概率型量本利分析为例，蒙特卡洛模拟的分析步骤如下： 1、分析评价参数的特征，如企业经营中的销售数量、销售价格、产品生产的变动成本以及固定成本等，并根据历史资料或专家意见，确定随机变量的某些统计参数； 2、按照一定的参数分布规律，在计算机上产生随机数，如利用EXCEL提供的RAND函数，模拟量本利分析的概率分布，并利用VLOOKUP寻找对应概率分布下的销售数量、销售价格、产品生产的变动成本以及固定成本等参数； 3、建立管理会计的数学模型，对于概率型量本利分析有如下关系式，产品利润=产品销售数量×（产品单位销售价格-单位变动成本）-固定成本，这里需要说明的是以上分析参数不是确定型的，是依据某些概率分布存在的； 4、通过足够数量的计算机仿真，如文章利用RAND、VLOOKUP等函数进行30000次的模拟，得到30000组不同概率分布的各参数的排列与组合，由于模拟的数量比较大，所取得的实验数据具有一定的规律性； 5、根据计算机仿真的参数样本值，利用函数MAX、MIN、A VERAGE等，求出概率型量本利分析评价需要的指标值，通过对大量的评价指标值的样本分析，得到量本利分析中的利润点可能的概率分布，从而掌握企业经营与财务中的风险，为财务决策提供重要的参考。三、概率型量本利分析与比较 （一）期望值分析方法假设某企业为生产与销售单一产品的企业，经过全面分析与研究，预计未来年度的单位销售价格、销售数量、单位变动成本和固定成本的估计值及相应的概率如表1，其中销售数量单位为件，其余反映价值的指标单位为元，试计算该企业的生产利润。表1概率型量本利分析参数 项目概率数值 单位销售价格0.3 40 0.4 43 0.3 45 单位变动成本0.4 16 0.2 18 0.4 20 固定成本0.6 28000 0.4 30000

蒙特卡罗也称统计模拟方法

蒙特卡罗也称统计模拟方法，是二十世纪四十年代中期由于科学技术的发展和电子计算机的发明，而被提出的一种以概率统计理论为指导的一类非常重要的数值计算方法。是指使用随机数（或更常见的伪随机数）来解决很多计算问题的方法。蒙特卡罗方法的名字来源于摩纳哥的一个城市蒙地卡罗，该城市以赌博业闻名，而蒙特·罗方法正是以概率为基础的方法。与它对应的是确定性算法。 蒙特卡罗方法在金融工程学，宏观经济学，计算物理学(如粒子输运计算、量子热力学计算、空气动力学计算)等领域应用广泛。 基本思想 当所求解问题是某种随机事件出现的概率，或者是某个随机变量的期望值时，通过某种“实验”的方法，以这种事件出现的频率估计这一随机事件的概率，或者得到这个随机变量的某些数字特征，并将其作为问题的解。有一个例子可以使你比较直观地了解蒙特卡罗方法：假设我们要计算一个不规则图形的面积，那么图形的不规则程度和分析性计算（比如，积分）的复杂程度是成正比的。蒙特卡罗方法是怎么计算的呢？假想你有一袋豆子，把豆子均匀地朝这个图形上撒，然后数这个图形之中有多少颗豆子，这个豆子的数目就是图形的面积。当你的豆子越小，撒的越多的时候，结果就越精确。在这里我们要假定豆子都在一个平面上，相互之间没有重叠。 工作过程 在解决实际问题的时候应用蒙特卡罗方法主要有两部分工作： 用蒙特卡罗方法模拟某一过程时，需要产生各种概率分布的随机变量。 用统计方法把模型的数字特征估计出来，从而得到实际问题的数值解。 计算步骤 使用蒙特卡罗方法进行分子模拟计算是按照以下步骤进行的： ① 使用随机数发生器产生一个随机的分子构型。 ②对此分子构型的其中粒子坐标做无规则的改变，产生一个新的分子构型。 ③计算新的分子构型的能量。 ④比较新的分子构型于改变前的分子构型的能量变化，判断是否接受该构型。 若新的分子构型能量低于原分子构型的能量，则接受新的构型，使用这个构型重复再做下一次迭代。 若新的分子构型能量高于原分子构型的能量，则计算玻尔兹曼常数，同时产生一个随机数。

马尔可夫链模型

马尔可夫链模型 马尔可夫链模型（Markov Chain Model） 目录 [隐藏] ? 1 马尔可夫链模型概述 ? 2 马尔可夫链模型的性质 ? 3 离散状态空间中的马尔可夫链 模型 ? 4 马尔可夫链模型的应用 o 4.1 科学中的应用 o 4.2 人力资源中的应用 ? 5 马尔可夫模型案例分析[1] o 5.1 马尔可夫模型的建 立 o 5.2 马尔可夫模型的应 用 ? 6 参考文献 [编辑] 马尔可夫链模型概述 马尔可夫链因安德烈·马尔可夫（Andrey Markov，1856－1922）得名，是数学中具有马尔可夫性质的离散时间随机过程。该过程中，在给定当前知识或信息的情况下，过去（即当期以前的历史状态）对于预测将来（即当期以后的未来状态）是无关的。 时间和状态都是离散的马尔可夫过程称为马尔可夫链, 简记为。 马尔可夫链是随机变量的一个数列。这些变量的范围，即他们所有可能 取值的集合，被称为“状态空间”，而Xn的值则是在时间n的状态。如果Xn + 1对于过去状态的条件概率分布仅是Xn的一个函数，则 这里x为过程中的某个状态。上面这个恒等式可以被看作是马尔可夫性质。

马尔可夫在1906年首先做出了这类过程。而将此一般化到可数无限状态空间是由柯尔莫果洛夫在1936年给出的。 马尔可夫链与布朗运动以及遍历假说这两个二十世纪初期物理学重要课题是相联系的，但马尔可夫寻求的似乎不仅于数学动机，名义上是对于纵属事件大数法则的扩张。 马尔可夫链是满足下面两个假设的一种随机过程： 1、t+l时刻系统状态的概率分布只与t时刻的状态有关，与t时刻以前的状态无关； 2、从t时刻到t+l时刻的状态转移与t的值无关。一个马尔可夫链模型可表示为=(S，P，Q)，其中各元的含义如下： 1）S是系统所有可能的状态所组成的非空的状态集，有时也称之为系统的状态空间，它可以是有限的、可列的集合或任意非空集。本文中假定S是可数集(即有限或可列)。用小写字母i,j(或S i,S j)等来表示状态。 2）是系统的状态转移概率矩阵，其中P ij表示系统在时刻t处于状态i，在下一时刻t+l处于状态i的概率，N是系统所有可能的状态的个数。对于任意i∈s，有 。 3）是系统的初始概率分布，q i是系统在初始时刻处于状态i的概率， 满足。 [编辑] 马尔可夫链模型的性质 马尔可夫链是由一个条件分布来表示的 P(X n + 1 | X n) 这被称为是随机过程中的“转移概率”。这有时也被称作是“一步转移概率”。二、三，以及更多步的转移概率可以导自一步转移概率和马尔可夫性质：

中天会计事务所马尔可夫模型例题(最完整的例题分析)

中天会计事务所马尔可夫模型例题一、问题分析 中天会计事务所由于公司业务日益繁忙，常造成公司事务工作应接不暇，解决该公司出现的这种问题的有效办法是要实施人力资源的供给预测技术。根据对该公司材料的深入分析，可采用马尔可夫模型这一供给预测方法对该事务所的人力资源状况进行预测。 马尔可夫分析法是一种统计方法，其方法的基本思想是：找出过去人力资源变动的规律，用以来推测未来人力变动的趋势。马尔可夫分析法适用于外在环境变化不大的情况下，如果外在环境变化较大的时候这种方法则难以用过去的经验情况预测未来。马尔可夫分析法的分析过程通常是分几个时期来收集数据，然后在得出平均值，利用这些数据代表每一种职位的人员变动频率，就可以推测出人员的变动情况。 二、项目策划 （一）第一步是编制人员变动概率矩阵表。 根据公司提供的内部资料：公司的各职位人员如下表1所示。 表1：各职位人员表 职位代号人数 合伙人P 40 经理M 80 高级会计师S 120 会计员 A 160 制作一个人员变动概率矩阵表，表中的每一个元素表示从一个时期到另一个时期（如从某一年到下一年）在两个工作之间调动的雇员数量的历年平均百分比（以小数表示）。（注：一般以3—5年为周期来估计年平均百分比。周期越长，根据过去人员变动所推测的未来人员变动就越准确。） 表2：历年平均百分比人员变动概率矩阵表 职位合伙人 P 经理M 高级会计师S 会计员A 职位年度离职升为 合伙 人 离职升为经 理 降为 会计 员 离职升为高级 会计师 离职 2005 0.20 0.08 0.13 0.07 0.05 0.11 0.12 0.11 2006 0.23 0.07 0.27 0.05 0.08 0.12 0.15 0.29 2007 0.17 0.13 0.20 0.08 0.03 0.10 0.17 0.20 2008 0.21 0.12 0.21 0.03 0.07 0.09 0.13 0.19 2009 0.19 0.10 0.19 0.02 0.02 0.08 0.18 0.21 平均0.20 0.10 0.20 0.05 0.05 0.10 0.15 0.20

人力供给预测之马尔科夫模型

人力供给预测之马尔科夫模型 马尔科夫模型是根据历史数据，预测等时间间隔点上的各类人员分布状况。此方法的基本思想是根据过去人员变动的规律，推测未来人员变动的趋势。因此，运用马尔科夫模型时假设——未来的人员变动规律是过去变动规律的延续。既是说，转移率要么是一个固定比率，要么可以通过历史数据以某种方式推算出。 步骤： （1）根据历史数据推算各类人员的转移率，得出转移率的转移矩阵； （2）统计作为初始时刻点的各类人员分布状况； （3）建立马尔科夫模型，预测未来各类人员供给状况。 运用马尔科夫模型可以预测一个时间段后的人员分布，虽然这个时间段可以自由定义，但较为普遍的是以一年为一个时间段，因为这样最为实用。在确定转移率时，最粗略的方法就是以今年的转移率作为明年的转移率，这种方法认为最近时间段的变化规律将继续保持到下一时间段。虽然这样很简便，但实际上一年的数据过于单薄，很多因素没有考虑到，一个数据的误差可能非常大。因为以一年的数据得出的概率很难保证稳定，最好运用近几年的数据推算。在推算时，可以采用简单移动平均法、加权移动平均法、指数平滑法、趋势线外推法等，可以在试误的过程中发现哪种方法推算的转移率最准确。尝试用不同的方法计算转移率，然后用这个转移率和去年的数据来推算今年的实际情况，最后选择与实际情况最相符的计算方法。转移率是一类人员转移到另一类人员的比率，计算出所有的转移率后，可以得到人员转移率的转移矩阵。 转移出i类人员的数量 i类人员的转移率 = （3-1） i类人员原有总量 人员转移率的转移矩阵： P11 P12 (1) P21 P22 (2) P = P31 P32 (3) （3-2） ┇┇┇ P K1 P K2 ……P KK 一般是以现在的人员分布状况作为初始状况，所以只需统计当前的人员分布情况即可。这是企业的基本信息，人力资源部门可以很容易地找到这些数据。 建立模型前，要对员工的流动进行说明。流动包括外部到内部、内部之间、内部到外部的流动，内部之间的流动可以是提升、降职、平级调动等。由于推测的是整体情况，个别特殊调动不在考虑之内。马尔科夫模型的基本表达式为：

(定价策略)期权定价中的蒙特卡洛模拟方法

期权定价中的蒙特卡洛模拟方法 期权作为最基础的金融衍生产品之一，为其定价一直是金融工程的重要研究领域，主要使用的定价方法有偏微分方程法、鞅方法和数值方法。而数值方法又包括了二叉树方法、有限差分法和蒙特卡洛模拟方法。 蒙特卡洛方法的理论基础是概率论与数理统计，其实质是通过模拟标的资产价格路径预测期权的平均回报并得到期权价格估计值。蒙特卡洛方法的最大优势是误差收敛率不依赖于问题的维数，从而非常适宜为高维期权定价。 §1. 预备知识 ◆两个重要的定理：柯尔莫哥洛夫(Kolmogorov)强大数定律和莱维一林德贝格(Levy-Lindeberg)中心极限定理。 大数定律是概率论中用以说明大量随机现象平均结果稳定性的一系列极限定律。在蒙特卡洛方法中用到的是随机变量序列同分布的Kolmogorov 强大数定律： 设12,,ξξL 为独立同分布的随机变量序列，若 [],1,2,k E k ξμ=<∞=L 则有1 1(lim )1n k n k p n ξμ→∞===∑ 显然，若12,,,n ξξξL 是由同一总体中得到的抽样，那么由 此大数定律可知样本均值1 1n k k n ξ=∑当 n 很大时以概率1收敛于

总体均值μ。 中心极限定理是研究随机变量之和的极限分布在何种情形下是正态的，并由此应用正态分布的良好性质解决实际问题。 设12,,ξξL 为独立同分布的随机变量序列，若 2 [],[],1,2,k k E D k ξμξσ=<∞=<∞=L (0,1)n k d n N ξ μ -??→∑ 其等价形式为2 1 1lim ()exp(),2n x k k n t n P x dt x ξμσ =→∞ -∞ -≤= --∞<<∞∑?。 ◆Black-Scholes 期权定价模型 模型的假设条件： 1、标的证券的价格遵循几何布朗运动 dS dt dW S μσ=+ 其中，标的资产的价格S 是时间t 的函数，μ为标的资产 的瞬时期望收益率，σ为标的资产的波动率，dW 是维纳过程。 2、证券允许卖空、证券交易连续和证券高度可分。 3、不考虑交易费用或税收等交易成本。 4、在衍生证券的存续期内不支付红利。 5、市场上不存在无风险的套利机会。 6、无风险利率r 为一个固定的常数。 下面，通过构造标的资产与期权的资产组合并根据无套利定价原理建立期权定价模型。首先，为了得到期权的微分形式，先介绍随机微积分中的最重要的伊藤公式。

蒙特卡罗(Monte Carlo)方法简介

蒙特卡罗(Monte Carlo)方法简介

蒙特卡罗(Monte Carlo)方法简介 蒙特卡罗(Monte Carlo)方法，也称为计算机随机模拟方法，是一种基于"随机数"的计算方法。 一起源 这一方法源于美国在第二次世界大战进研制原子弹的"曼哈顿计划"。Monte Carlo方法创始人主要是这四位：Stanislaw Marcin Ulam, Enrico Fermi, John von Neumann（学计算机的肯定都认识这个牛人吧）和Nicholas Metropolis。 Stanislaw Marcin Ulam是波兰裔美籍数学家，早年是研究拓扑的，后因参与曼哈顿工程，兴趣遂转向应用数学，他首先提出用Monte Carlo方法解决计算数学中的一些问题，然后又将其应用到解决链式反应的理论中去，可以说是MC方法的奠基人；Enrico Fermi是个物理大牛，理论和实验同时都是大牛，这在物理界很少见，在“物理大牛的八卦”那篇文章里提到这个人很多次，对于这么牛的人只能是英年早逝了（别说我嘴损啊，上帝都嫉妒！）；John von Neumann可以说是计算机界的牛顿吧，太牛了，结果和Fermi一样，被上帝嫉妒了；Nicholas Metropolis，希腊裔美籍数学家，物理学家，计算机科学家，这个人对Monte Carlo方法做的贡献相当大，正式由于他提出的一种什么算法（名字忘了），才使得Monte Carlo方法能够得到如此广泛的应用，这人现在还活着，与前几位牛人不同，Metropolis很专一，他一生主要的贡献就是Monte Carlo方法。 蒙特卡罗方法的名字来源于摩纳哥的一个城市蒙地卡罗，该城市以赌博业闻名，而蒙特?罗方法正是以概率为基础的方法。与它对应的是确定性算法。 二解决问题的基本思路 Monte Carlo方法的基本思想很早以前就被人们所发现和利用。早在17世纪，人们就知道用事件发生的"频率"来决定事件的"概率"。19世纪人们用投针试验的方法来决定圆周率π。本世纪40年代电子计算机的出现，特

蒙特卡罗方法及应用实验讲义2016

蒙特卡罗方法及应用 实验讲义 东华理工大学核工系 2016.8

实验一 蒙特卡罗方法基本思想 一、实验目的 1、了解蒙特卡罗方法方法的基本思想； 2、掌握蒙特卡罗方法计算面积、体积的方法； 3、掌握由已知分布的随机抽样方法。 二、实验原理 Monte Carlo 方法，又称统计模拟方法或计算机随机模拟方法，是一种基于“随机数”进行数值模拟的方法，一种采用统计抽样理论近似求解物理或数学问题的方法。 如待求量可以表述成某些特征量的期望值、某些事件出现的概率或两者的函数形式，那么可采用蒙特卡罗方法求解。在求解某些特征量的期望值或某些事件出现的概率时，必须构建合符实际的数学模型。例如采用蒙特卡罗方法计算某函数所围面积时，构建的数学模型是构造一已知面积的可均匀抽样区域，在该区域投点，由伯努利定理大数定理可知，进入待求区域投点的频率依概率1收敛于该事件出现的概率（面积之比）。 由已知分布的随机抽样方法指的是由已知分布的总体中抽取简单子样。具体方法很多，详见教材第三章。 三、实验内容 1、安装所需计算工具（MATLAB 、fortran 、C++等）； 2、学习使用rand(m,n)、unifrnd(a,b,m,n)函数 3、求解下列问题： 3.0、蒲丰氏投针求圆周率。 3.1、给定曲线y =2 – x 2 和曲线y 3 = x 2，曲线的交点为：P 1( – 1，1 )、P 2( 1，1 )。曲线围成平面有限区域，用蒙特卡罗方法计算区域面积； 3.2 、计算1z z ?≥??≤??所围体积 其中{(,,)|11,11,02}x y z x y z Ω=-≤≤-≤≤≤≤。 4、对以下已知分布进行随机抽样：

蒙特卡洛模拟法

蒙特卡洛模拟法 蒙特卡洛(Monte Carlo)模拟是一种通过设定随机过程，反复生成时间序列，计算参数估计量和统计量，进而研究其分布特征的方法。具体的，当系统中各个单元的可靠性特征量已知，但系统的可靠性过于复杂，难以建立可靠性预计的精确数学模型或模型太复杂而不便应用时，可用随机模拟法近似计算出系统可靠性的预计值；随着模拟次数的增多，其预计精度也逐渐增高。由于涉及到时间序列的反复生成，蒙特卡洛模拟法是以高容量和高速度的计算机为前提条件的，因此只是在近些年才得到广泛推广。 这个术语是二战时期美国物理学家Metropolis执行曼哈顿计划的过程中提出来的。 蒙特卡洛模拟方法的原理是当问题或对象本身具有概率特征时，可以用计算机模拟的方法产生抽样结果，根据抽样计算统计量或者参数的值；随着模拟次数的增多，可以通过对各次统计量或参数的估计值求平均的方法得到稳定结论。 蒙特卡洛模拟法的应用领域 蒙特卡洛模拟法的应用领域主要有： 1.直接应用蒙特卡洛模拟:应用大规模的随机数列来模拟复杂系统，得到某些参数或重要指标。 2.蒙特卡洛积分:利用随机数列计算积分，维数越高，积分效率越高。 3.MCMC:这是直接应用蒙特卡洛模拟方法的推广，该方法中随机数的产生是采用的马尔科夫链形式。 （也叫随机模拟法）当系统中各个单元的可靠性特征量已知，但系统的可靠性过于复杂，难以建立可靠性预计的精确数学模型或模型太复杂而不便应用则可用随机模拟法近似计算出系统可靠性的预计值。随着模拟次数的增多，其预计精度也逐渐增高。由于需要大量反复的计算，一般均用计算机来完成。 应用此方法求解工程技术问题可以分为两类:确定性问题和随机性问题。解题步骤如下： 1.根据提出的问题构造一个简单、适用的概率模型或随机模型，使问题的解对应于该模型中随机变量的某些特征（如概率、均值和方差等），所构造的模型在主要特征参量方面要与实际问题或系统相一致 2 .根据模型中各个随机变量的分布，在计算机上产生随机数，实现一次模拟过程所需的足够数量的随机数。通常先产生均匀分布的随机数，然后生成服从某一分布的随机数，方可进行随机模拟试验。

蒙特卡洛方法及其在风险评估中的应用(1)

蒙特卡洛方法及其应用 1风险评估及蒙特卡洛方法概述 １.１蒙特卡洛方法。 蒙特卡洛方法，又称随机模拟方法或统计模拟方法，是在20世纪40年代随着电子计算机的发明而提出的。它是以统计抽样理论为基础，利用随机数，经过对随机变量已有数据的统计进行抽样实验或随机模拟，以求得统计量的某个数字特征并将其作为待解决问题的数值 解。 蒙特卡洛模拟方法的基本原理是：假定随机变量X1、X2、X3……X n、Y，其中X1、X2、X3……X n 的概率分布已知，且X1、X2、X3……X n、Y有函数关系：Y=F（X1、X2、X3……X n），希望求得随机变量Y的近似分布情况及数字特征。通过抽取符合其概率分布的随机数列X1、X2、X3……X n带入其函数关系式计算获得Y的值。当模拟的次数足够多的时候，我们就可以得到与实际情况相近的函数Y的概率分布和数字特征。 蒙特卡洛法的特点是预测结果给出了预测值的最大值,最小值和最可能值，给出了预测 值的区间范围及分布规律。 １.２风险评估概述。 风险表现为损损益的不确定性，说明风险产生的结果可能带来损失、获利或是无损失也无获利，属于广义风险。正是因为未来的不确定性使得每一个项目都存在风险。对于一个公司而言，各种投资项目通常会具有不同程度的风险，这些风险对于一个公司的影响不可小视，小到一个项目投资资本的按时回收，大到公司的总风险、公司正常运营。因此，对于风险的 测量以及控制是非常重要的一个环节。 风险评估就是量化测评某一事件或事物带来的影响的可能程度。根据“经济人”假设，收益最大化是投资者的主要追求目标，面对不可避免的风险时，降低风险，防止或减少损失， 以实现预期最佳是投资的目标。 当评价风险大小时，常有两种评价方式：定性分析与定量分析法。定性分析一般是根据风险度或风险大小等指标对风险因素进行优先级排序，为进一步分析或处理风险提供参考。这种方法适用于对比不同项目的风险程度，但这种方法最大的缺陷是在于，在多个项目中风险最小者也有可能亏损。而定量分析法则是将一些风险指标量化得到一系列的量化指标。通过这些简单易懂的指标，才能使公司的经营者、投资者对于项目分风险有正确的评估与判断，