蒙特卡罗方法1课件

,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,计算物理,125.217.162.13/lesson/ComputationalPhysics,蒙特卡罗方法,计算物理125.217.162.13/less,1,蒙特卡罗方法,蒲丰投针,收敛性、误差和优缺点,任意分布的随机数,粒子输运问题,随机过程模拟,梅氏抽样,蒙特卡罗方法蒲丰投针,2,蒲丰投针(1/5),蒙特卡罗方法,又称随机抽样技巧或统计试验方法,以概率统计理论为基础的,能够比较逼真地描述事物的特点及物理实验过程,解决一般数值方法难以解决的问题,随着电子计算机的发展而发展,首先在核武器的试验与研制中得到了应用,蒲丰投针,法国数学家蒲丰的1777年出版的著作：“在平面上画有一组间距为,d,的平行线，将一根长度为,l,(,l,d,),的针任意掷在这个平面上，求此针与平行线中任一条相交的概率。”,蒲丰投针(1/5)蒙特卡罗方法蒲丰投针,3,蒲丰投针(2/5),步骤,在桌面上画出,间距为,2,d,的平行线,准备长度为,2,l,(,l,d,),的针,向,桌面,随机,投,针,如果针与平行线相交，则计数器,n,加,1,计算：计数器,n,与总,投,针数,N,的比例(视作相交概率,P,),概率分析,P,=,?,各条平行线地位等同，仅考虑某条平行线附近的情况,平行线方向的,x,坐标,对概率没影响,针的中点的,y,坐标在线之间等概率落入(均匀分布在,0,d,)，仅当,y,l,才可能,针-,线,相交,针-,线的夹角,q,均匀分布在,0,p,，,q,与,y,独立,x,y,蒲丰投针(2/5)步骤概率分析 P=?xy,4,蒲丰投针(3/5),x,y,概率,P,=,2,l,/,(,p,d,),，可计算圆周率,实验者,时间,针长,总投数,相交数,p,值,Wolf,1850,0.8,5,000,2,532,3.1596,Smith,1855,0.6,3,204,1,218,3.1554,De Morgan,C,1860,1.0,600,382,3.137,Fox,1884,0.75,1,030,489,3.1595,Lezzerini,1901,0.83,3,408,1,808,3.1415929,Reina,1925,0.5419,2,520,859,3.1795,蒲丰投针(3/5)xy概率 P=2l/(p d)，,5,蒲丰投针(4/5),关于蒙特卡罗方法的分析和总结,基本思想,确定所求问题的解是某事件的概率(或某随机变量的数学期望、或与概率/数学期望有关的量，如,p,),通过试验方法，得出事件发生的频率(或该随机变量若干个具体观察值的算术平均值，如,P,)，求解,数学期望与概率：当随机变量的取值仅,为 1 或 0 时,，它的数学期望就是事件的概率；反之亦然,数学期望与算术平均值,用随机试验的方法计算积分，将积分看作服从分布,密度函数,f,(,r,),的随机变量,g,(,r,),的,数学期望,通过试验，得到,N,个观察值,r,1,r,2,r,N,(,从,f,(,r,),中抽取,N,个子样,r,1,r,2,r,N,),，求,g,(,r,),的算术平均值,蒲丰投针(4/5)关于蒙特卡罗方法的分析和总结通过试验，得到,6,蒲丰投针(5/5),试验方法和次数,试验方法不一定可行,精确的近似解需要巨量的试验次数，但人工方法作大量的试验相当困难，甚至是不可能的,用计算机模拟随机试验过程，完成巨量的次数,抽象,x,的分布密度函数：,q,的分布密度函数：,产生任意的,(,x,q,),=,由,f,1,(,x,),抽样,x,+,由,f,2,(,q,),抽样,q,对应的,随机变量：,数学期望与算术平均值,新的问题：误差？收敛？,蒲丰投针(5/5)试验方法和次数q 的分布密度函数：产生任意,7,收敛性、误差和,优缺点,(1/4),收敛性,求解：以,随机变量,X,的简单,子样,X,1,X,2,X,N,的算术平均值，作为求解的近似值,近似值的收敛性,大数定理：当试验次数足够多时，事件发生的频率无穷接近于该事件发生的概率,如果,X,1,X,N,独立分布，且期望值有限(,E,(,X,),)，那么,当,随机变量,X,的简单子样数,N,充分大时，其算术平均值以概率,1,收敛于期望值,E,(,X,),收敛性、误差和优缺点(1/4)收敛性近似值的收敛性 当随机,8,收敛性、误差和,优缺点,(2/4),误差,概率论的中心极限定理：如果随机变量序列,X,1,X,2,X,N,独立分布，且具有有限非零的方差,s,2,，即,其中的,f,(,x,),是,X,的分,布密度函数，则,当,N,足够大时,(蒙特卡罗方法),不等式 的概率约,1-,a,误差定义为 ，收敛速度为,a,0.500,0.050,0.003,l,a,0.674,1.960,2.968,收敛性、误差和优缺点(2/4)误差 其中的 f(x)是,9,收敛性、误差和,优缺点,(3/4),蒙特卡罗方法的误差为概率误差,均方差,s,是未知的，估计值为,减少误差的技巧,(在确定的置信度,a,前提,下),误差,e,与试验次数的开根号,N,1/2,成反比：精度,一个数量级，次数,N,两个数量级巨大的代价,误差,e,与均方差,s,成正比：精度,一个数量级，均方差,s,一个数量级可接受的代价,效率,降低方差,增加观察子样的时间,固定时间内样本数减少,代价增加,蒙特卡罗方法中的效率：由均方差,s,和观察一个子样的费用(计算机时),c,衡量,=,s,2,c,收敛性、误差和优缺点(3/4)蒙特卡罗方法的误差为概率误差减,10,收敛性、误差和,优缺点,(4/4),优点,较逼真地描述具有随机性质的事物的特点及物理实验过程,受几何条件限制小,收敛速度与问题的维数无关,具有同时计算多个方案与多个未知量的能力,误差容易确定,程序结构简单，易于实现,缺点,收敛速度慢,误差具有概率性,在某些问题中，计算结果与系统大小有关,主要应用范围,粒子输运，统计物理，典型数学，真空技术，激光技术，医学，生物，探矿等,收敛性、误差和优缺点(4/4)优点,11,任意分布的随机数,(1/11),随机数(,蒙特卡罗方法的关键),基本的：,均匀,地在,(0,1),内分布,任意的：,非均匀,地在,(,a,b,),内分布,方法：直接的(离散、连续)、舍选的(简单、乘、乘加),直接抽样方法,离散型：产生随机量,x,的抽样值,x,i,，概率为,P,i,(,i,=,1,2,),方法：先计算,x,的分布函数,对于产生的,R,，如果满足,F,i,-1,R,F,i,，则令抽样,x,=,x,i,证明：,任意分布的随机数(1/11)随机数(蒙特卡罗方法的关键)方法,12,任意分布的随机数,(2/11),例：二项分布,直接抽样方法,：,产生,R,，如果满足,，则令,x,=,n,例：泊松分布,直接抽样方法,：,产生,R,，如果满足,，则令,x,=,n,任意分布的随机数(2/11)例：二项分布 直接抽样,13,任意分布的随机数,(3/11),连续型：产生随机量,x,的抽样值，概率密度函数为,f,(,x,),方法：先计算,x,的分布函数,对于随机数,R,，解方程,R,=,F,(,x,)，得到,x,=,F,-1,(,R,),，则令,x,=,x,证明：,任意分布的随机数(3/11)连续型：产生随机量 x 的抽样值,14,任意分布的随机数,(4/11),例：均匀介质中，粒子运动的自由程,S,是随机量，其概率密度函数为,f,(,S,),=,S,e,-S,S,；其中的,S,=,n,s,是宏观总截面，,s,是原子截面，,n,为介质中的原子数密度,分布函数：,介质中粒子的碰撞过程,随机量,R,中的抽样过程,例：散射方位角余弦分布,任意分布的随机数(4/11)例：均匀介质中，粒子运动的自由程,15,舍,任意分布的随机数,(5/11),舍选抽样方法,直接法的特点：简单，但分布函数的反函数不一定有解析解,简单分布,定理：如果,Z,是,a,b,上均匀分布的随机数，那么利用条件,f,(,Z,)/,M,R,(,M,是,f,(,Z,),的上界)选出的,Z,将是,a,b,上概率密度函数为,f,(,Z,),的随机数,证明：如右图所示,d,Z,(,Z,R,),R,Z,f,(,Z,),/,M,a,b,1,0,Z,为横坐标，,R,为纵坐标，实曲线为函数,f,(,Z,)/,M,R,在,(,0,1,),内、,Z,在,a,b,内均匀分布,随机点,(,Z,R,),在虚框内,均匀分布,随机点落入窄条的概率,=,两面积之比,选,已知概率密度函数,抽样方法,舍任意分布的随机数(5/11)舍选抽样方法 证明：如右图所,16,任意分布的随机数,(6/11),简单分布抽样方法的流程,产生随机数,x,和,R,计算,Z,=,a,+(,b,-,a,),x,f,(,Z,),/,M,R,?,Z,=,Z,Y,N,效率(有一部分随机数被舍弃),舍,d,Z,(,Z,R,),R,Z,f,(,Z,),/,M,a,b,1,0,选,例：受限的散射方位角余弦分布,效率,任意分布的随机数(6/11)简单分布抽样方法的流程产生随机数,17,任意分布的随机数,(7/11),乘分布：,f,(,x,),有锐锋时效率很低，需要改进的简单分布,方法：将函数写成,f,(,Z,),=,h,(,Z,),f,1,(,Z,),，由容易抽样的,f,1,(,Z,),抽样出,Z,，代入,h,(,Z,),，,如果满足条件,h,(,Z,),/,M,R,(,M,是,h,(,Z,),的上界)，那么得到概率密度函数为,f,(,Z,),的抽样值,舍,证明：如右图所示,d,F,1,(,F,1,R,),R,F,1,(,Z,),h,(,Z,),/,M,1,1,0,f,1,(,Z,),的分布函数,F,1,(,Z,),为横坐标，,R,为纵坐标，实曲线为函数,h,(,Z,),/,M,R,和,F,1,在,(,0,1,),内均匀分布,随机点,(,F,1,R,),在虚框内,均匀分布,随机点落入窄条的概率,=,两面积之比,选,效率,任意分布的随机数(7/11)乘分布：f(x)有锐锋时效率很,18,任意分布的随机数,(8/11),乘分布抽样方法的流程,产生随机数,F,1,和,R,利用,F,1,，,由,f,1,(,Z,),抽样,Z,h,(,Z,),/,M,R,?,Z,=,Z,Y,N,例：半正态分布，概率密度函数,任意分布的随机数(8/11)乘分布抽样方法的流程产生随机数利,19,任意分布的随机数,(9/11),乘加分布,方法：如果随机量,x,具有以下乘加形式,其中,f,1,(,x,),和,f,2,(,x,),满足密度函数的要求，即,并且,h,1,(,x,),和,h,2,(,x,),满足,那么从,f,(,x,),的以下变形，可以得到,乘加分布,其中,和,为加权因子,证明：参照直接法和舍选法(略),任意分布的随机数(9/11)乘加分布 其中 f1(x)和,20,任意分布的随机数,(10/11),乘加分布抽样方法的流程,产生随机数,g,F,R,利用,F,由,f,1,(,x,),抽样,x,1,g,g,1,?,Y,N,h,1,(,x,1,)/,M,1,R,?,利用,F,由,f,2,(,x,),抽样,x,2,h,2,(,x,2,)/,M,2,R,?,Y,Y,N,N,x,=,x,1,x,=,x,2,效率,任意分布的随机数(10/11)乘加分布抽样方法的流程产生随机,21,任意分布的随机数,(11/11),例：散射光子能量抽样。能量为,E,0,的入射光子，经原子散射后的能量,E,按某种概率分布。如果令,x,=,E,0,/,E,，那么概率密度函数,f,(,x,),为(,其中,K,=,K,(,E,0,),为归一因子,),任意分布的随机数(11/11)例：散射光子能量抽样。能量为,22,粒子输运问题(1/7),蒙特卡罗模拟,粒子输运是随机过程，运动规律是大量粒子运动的统计,蒙特卡罗模拟：模拟一定数量粒子在介质中的运动，再现粒子运动的统计规律,例(平板介质模型)：由单一物质组成的均匀介质，厚度为,H,，能量为,E,00,的平行光子束垂直射入板内，求光子对板的投射率,H,抽样：自由程、作用类型、散射能量、散射方向,自由程抽样(PPT15)：,粒子运动的自由程,S,是随机量，其概率密度函数为,f,(,S,),=,S,e,-S,S,；其中的,S,=,n,s,是宏观总截面，,s,是原子截面，,n,为介质中的原子数密度,分布函数：,F,(,S,),=,1,-,e,-S,S,粒子的碰撞过程,随机量,R,中的,抽样过程：,S,=,-,ln,R,/,S,粒子输运问题(1/7)蒙特卡罗模拟例(平板介质模型)：由单一,23,粒子输运问题(2/7),作用类型抽样,光子-介质：散射(康普顿散射)和吸收(光电效应),散射截面,吸收截面,总截面：,S,t,=,S,s,+,S,a,光子随机地被散射或吸收的过程,随机量,R,中的,抽样过程：如果,S,s,/,S,t,R,，则散射，否则为吸收,粒子输运问题(2/7)作用类型抽样 吸收截面 总截面：S,24,粒子输运问题(3/7),散射能量抽样(PPT22)：能量为,E,0,的入射光子，经原子散射后的能量,E,按某种概率分布。如果令,x,=,E,0,/,E,，那么概率密度函数,f,(,x,),为(,其中,K,=,K,(,E,0,),为归一因子,),粒子输运问题(3/7)散射能量抽样(PPT22)：能量为 E,25,粒子输运问题(4/7),散射方向抽样,坐标系：以入射方向 (由,q,和,f,确定)为参考系，散射方向,由,散射角,q,和,散射方位角,f,确定,x,y,z,流程图,输入,g,0,E,0,E,计算,m,=,1,+1/,E,0,-1/,E,计算,g,产生随机数,R,计算,f,=,2,p,R,粒子输运问题(4/7)散射方向抽样坐标系：以入射方向 (,26,粒子输运问题(5/7),直接模拟方法,模拟粒子在介质中的真实物理过程,粒子在介质中的状态：空间位置，能量和运动方向,碰撞点的状态参数,(,Z,m,E,m,g,m,),表示从源出发的粒子在介质中经过,m,次碰撞后的状态,Z,H,O,Z,0,E,0,q,0,Z,1,E,1,q,1,Z,2,E,2,q,2,粒子的运动过程,=,碰撞点的状态序列,模拟粒子的运动过程,=,确定,状态序列的问题,粒子输运问题(5/7)直接模拟方法碰撞点的状态参数(Zm,27,粒子输运问题(6/7),记录和计算的内容,穿透率和误差估计(设,N,为入射粒子数，,N,1,为透射数),透射粒子的能量和方向分布,离散化能量：,E,min,=,E,0,E,1,E,i,E,I,离散化角度：,0,=,q,0,q,1,q,j,M,?,N,Y,i,N,?,N,Y,Z,H,?,N,Y,粒子输运问题(7/7)输入N,H,光子 i 的初态m=0,29,随机过程模拟,(1/3),蒙特卡罗模拟,随机过程,的两种情况,已知随机过程的概率分布函数,随机过程的统计特征,建立与,随机过程的模型，形成随机量，并使其数字特征(概率、平均值、方差等)是问题的解,由已知的概率分布函数进行大量的抽样,统计处理抽样结果，给出问题的解及其误差,已知随机现象的观测数据,概率分布函数,分析现象的特征，假设随机量服从某种分布，建立概率模型,由观测数据推断分布中的参数,按推断的分布进行大量的抽样,比较抽样值和观测值，根据误差，判断假设的准确性,随机过程模拟(1/3)蒙特卡罗模拟随机过程的两种情况,30,随机过程模拟,(2/3),例：,a,粒子衰变的蒙特卡罗模拟,随机现象的观测数据,概率分布函数,观测数据：每隔,D,t,Ng,?,Y,N,m,m,+1,计算平均值,方差,误差,主程序,结束,梅氏抽样(1/1)特别适合于多维随机量的系统梅氏1步输入产生,33,作业,用蒲丰投针在计算机上计算,p,值，取,d,=,4,l,=,3,分子热运动的速率概率密度函数是,取,f,1,(,x,),=,2,e,-2,x,/3,/3,，找出,h,(,x,),和,M,，用直接法对,f,1,(,x,),抽样并写出程序,作业用蒲丰投针在计算机上计算 p 值，取 d=4,l,34,