第八讲Monte-Carlo方法及应用案例

MonteCarlo方法一、MonteCarlo方法简介本章主要内容二、蒙特卡罗方法的应用MonteCarlo方法1.1统计力学层次的计算机分子模拟1、宏观化学现象是~1024个分子（原子）的集体行为，固有统计属性2、量子力学方法的局限性：对象为平衡态、单分子或几个分子组成的体系；不适用于动力学过程和有温度压力变化的体系。一、MonteCarlo方法简介•分子模拟的两种主要方法：⑴分子动力学法(MD，MolecularDynamics)基于粒子运动的经典轨迹⑵MonteCarlo法(MC)基于概率和统计力学MonteCarlo方法MonteCarlo原为地中海沿岸Monaco(摩纳哥)的一个城市的地名,是世界闻名的大赌场，MonteCarlo方法的随机抽样特征在它的命名上得到了反映。1.2MonteCarlo方法的发展历史MonteCarlo方法美国LosAlamos(洛斯·阿拉莫斯)实验室中子输运和辐射输运等物理过程用MonteCarlo来解决实际问题却始于本世纪40年代。Metropolis等人在这一时期的工作主要就是对中子扩散进行随机抽样计算机模拟，得出所要求算的相关参数，并把这种随机抽样方法命名为MonteCarlo方法。随着电子计算机的迅速发展，人们开始有意识地、广泛、系统地应用随机抽样方法来解决大量的数学、物理和化学等方面的问题，并且将MonteCarlo方法作为一门独立的计算方法进行研究，并随之向各个学科领域渗透。MonteCarlo方法MonteCarlo方法解决的问题：1、问题本身是确定性问题，要求我们去寻找一个随机过程，使该随机过程的统计平均是所求问题的解；2、问题本身就是一个随机过程，可根据问题本身的实际过程来进行计算机模拟，并采用统计方法来求得问题的解。1.3MonteCarlo方法简介基本思想：当所求的问题是某种事件出现的概率，或者是某个随机变量的期望值时，它们可以通过某种“随机试验”的方法，得到这种事件出现的概率，或者得到这个随机变量的统计平均值，并用它们作为问题的解。MonteCarlo方法1、问题本身是确定性问题MonteCarlo方法圆周率、定积分?2、问题本身就是一个随机过程MonteCarlo方法MonteCarlo方法应用：1、数学：本身已形成计算数学的一个分支；2、粒子物理：输运问题、屏蔽问题、核武器试验分析等；3、统计物理、化学，材料、工程各领域；4、其它：疾病传播与免疫、系统工程与管理优化等等。MonteCarlo方法蒙特卡罗模拟流程图MonteCarlo方法产生随机数，进行多次重复试验。统计分析模拟试验结果，给出问题的概率解以及解的精度估计。构造一个简单、适用的概率模型上图是叶片的横截面，叶片外表面与高温燃气接接触，高温燃气的温度通过热障涂层传递到叶片的外表面，叶片内表面用冷却介质进行冷却。由于热障涂层厚度会在一个小范围内波动，求叶片表面温度变化范围。热障涂层热传导的MC模拟MonteCarlo方法上图是叶片局部放大图，Tgas-燃气温度，TTBC-热障涂层温度，TMH-金属叶片外表面温度，TMC-金属叶片内表面温度，Tcool-冷却介质温度MonteCarlo方法叶片冷却边金属叶片热障涂层涂层接触燃气TgasTTBCTMHTcoolTMC模型建立MonteCarlo方法金属的厚度；热障涂层的厚度；金属的热导率；热障涂层的热导率；；冷却介质表面传热系数燃气表面传热系数；热流密度；MTBCMTBCcoolgascoolMCcoolMCMHMMMHTBCTBCTBCTBCgasgasLLqTTqTTLqTTLqTTq)()()()(模型建立MonteCarlo方法coolcoolgasgasMCMHTBCcoolMMMMTBCTBCTBCTBCgasTT-qTTT1-01-L-L1-0L-L1-00-0000求解线性方程组的程序MonteCarlo方法function[Ttbc,Tmh,Tmc,q]=mtkl1(hgas,Tgas,ktbc,Ltbc,km,Lm,hcool,Tcool)%计算矩阵K=[-hgas,0,0,-1;ktbc/Ltbc,-ktbc/Ltbc,0,-1;0,km/Lm,-km/Lm,-1;0,0,hcool,-1;];%计算等式右侧量b=[-hgas*Tgas;0;0;hcool*Tcool];u=K\b;%输出计算结果Ttbc=u(1);Tmh=u(2);Tmc=u(3);q=u(4);模型求解MonteCarlo方法CTCTmLmLCmWCmWCmWCmWTTLLcoolgasMTBCMTBCcoolgascoolgasMTBCMTBCcoolgas2001300003.00005.0)/(5.21),/(1)/(1000),/(300022，，，，，，，，，输入量：MonteCarlo方法hgas=3000;hcool=1000;ktbc=1;km=21.5;Ltbc=0.0005;Lm=0.003;Tcool=200;Tgas=1300;[Ttbc,Tmh,Tmc,q]=mtkl1(hgas,Tgas,ktbc,Ltbc,km,Lm,hcool,Tcool)Ttbc=1.1141e+03Tmh=835.3635Tmc=757.5639q=5.5756e+05涂层厚度单参数均匀分布的蒙特卡罗法模拟MonteCarlo方法范围内均匀分布在设涂层厚度m.LTBC00075.0000250其概率密度函数如右图涂层厚度单参数均匀分布的蒙特卡罗法模拟程序MonteCarlo方法clearall;%参数的标称值hgas=3000;%TBC-气体热传递系数(W/m^2)Tgas=1300;%混合气体温度(c)ktbc=1;%TBC热传导。(W/mK)km=21.5%金属温度传导。（W/mK）Lm=0.003%金属厚度（m）hcool=1000;%冷却液—金属热传递系数(W/m^2)Tcool=200;%冷却液温度(c)%蒙特卡洛试验次数Ntrial=1000;forn=1:Ntrial,Ltbc(n)=0.00025+0.0005*rand;%用均匀分布产生Ltbc的值。[Ttbc,Tmh(n),Tmc,q]=mtkl1(hgas,Tgas,ktbc,Ltbc(n),km,Lm,hcool,Tcool);if(n1),mTmh(n-1)=mean(Tmh);sTmh(n-1)=std(Tmh);end%plot(Ltbc(n),Tmh(n),’*’);holdon;%drawnow;endfprintf('MeanTmh=%f\n',mTmh(Ntrial-1));fprintf('StdTmh=%f\n',sTmh(Ntrial-1));subplot(211)hist(Ltbc,20);xlabel('Ltbc(m)');subplot(212);hist(Tmh,20);xlabel('Tmh(C)');figure;subplot(211);plot(mTmh);ylabel('MeanTmh');subplot(212);plot(sTmh);ylabel('StdTmh');涂层厚度与叶片外表面温度在不同区间分布数MonteCarlo方法2345678x10-4020406080Ltbc (m)760780800820840860880900920940020406080Tmh (C)模拟实验次数与叶片外表面温度及方差的关系MonteCarlo方法01002003004005006007008009001000780800820840860MeanTmh01002003004005006007008009001000020406080StdTmh随着模拟次数增大，叶片外表面温度偏离程度减少。多参数均匀分布的蒙特卡罗法模拟MonteCarlo方法对于多个输入的情况，和单变量类似，只是需要对每个输入都产生一系列随机数。问题：随着随机输入的个数增加，TMH的分布将会发生怎样的变化？多参数均匀分布的蒙特卡罗法模拟程序MonteCarlo方法clearall;%输入参数的上下界hgas=[1500,4500];%TBC-气体热传递系数(W/m^2)Tgas=[1200,1400];%混合气体温度(C)ktbc=[0.9,1.1];%TBC热传导。(W/mK)Ltbc=[0.00025,0.00075];%TBC厚度(m)km=[20.0,23.0];%金属温度传导。（W/mK）Lm=[0.002,0.004];%金属厚度（m）hcool=[500,1500];%冷却液—金属热传递系数(W/m^2)Tcool=[150,250];%冷却液温度(C)%将参数边界放进一个向量中Pbound=[hgas;Tgas;ktbc;Ltbc;km;Lm;hcool;Tcool];Ntrial=1000;forn=1:Ntrial,%用均匀分布产生各个参数的值。P(:,n)=Pbound(:,1)+(Pbound(:,2)-Pbound(:,1)).*rand;%求解热传导问题[Ttbc,Tmh(n),Tmc,q]=mtkl1(P(1,n),P(2,n),P(3,n),P(4,n),P(5,n),P(6,n),P(7,n),P(8,n));if(n1),mTmh(n-1)=mean(Tmh);sTmh(n-1)=std(Tmh);end%plot(Ltbc(n),Tmh(n),'*');holdon;%drawnow;endfprintf('MeanTmh=%f\n',mTmh(Ntrial-1));fprintf('StdTmh=%f\n',sTmh(Ntrial-1));subplot(311);hist(Tmh,20);xlabel('Tmh(C)');subplot(312);plot(mTmh);ylabel('MeanTmh');subplot(313);plot(sTmh);ylabel('StdTmh');叶片外表面温度在不同区间分布数模拟实验次数与叶片外表面温度及方差的关系MonteCarlo方法800810820830840850860870880050100Tmh (C)01002003004005006007008009001000820840860MeanTmh01002003004005006007008009001000152025StdTmh多参数非均匀分布的蒙特卡罗法模拟MonteCarlo方法怎样处理非均匀分布的问题在MATLAB（和很多模拟软件包）中，有现成的函数产生服从正态分布的输入随机数。MATLAB中“randn”就产生一个均值为零，方差为1的正态随机变量。人有了知识，就会具备各种分析能力，明辨是非的能力。所以我们要勤恳读书，广泛阅读，古人说“书中自有黄金屋。”通过阅读科技书籍，我们能丰富知识，培养逻辑思维能力；通过阅读文学作品，我们能提高文学鉴赏水平，培养文学情趣；通过阅读报刊，我们能增长见识，扩大自己的知识面。有许多书籍还能培养我们的道德情操，给我们巨大的精神力量，鼓舞我们前进。