第八講Monte-Carlo方法及應(yīng)用

MonteCarlo方法一、MonteCarlo方法簡介本章主要內(nèi)容二、蒙特卡羅方法的應(yīng)用MonteCarlo方法1.1統(tǒng)計(jì)力學(xué)層次的計(jì)算機(jī)分子模擬1、宏觀化學(xué)現(xiàn)象是~1024個分子（原子）的集體行為，固有統(tǒng)計(jì)屬性2、量子力學(xué)方法的局限性：對象為平衡態(tài)、單分子或幾個分子組成的體系；不適用于動力學(xué)過程和有溫度壓力變化的體系。一、MonteCarlo方法簡介分子模擬的兩種主要方法： ⑴ 分子動力學(xué)法

(MD，MolecularDynamics)

基于粒子運(yùn)動的經(jīng)典軌跡 ⑵ MonteCarlo法(MC)

基于概率和統(tǒng)計(jì)力學(xué)MonteCarlo方法MonteCarlo原為地中海沿岸Monaco(摩納哥)的一個城市的地名,是世界聞名的大賭場，MonteCarlo方法的隨機(jī)抽樣特征在它的命名上得到了反映。1.2MonteCarlo方法的發(fā)展歷史MonteCarlo方法

美國LosAlamos(洛斯·阿拉莫斯)實(shí)驗(yàn)室中子輸運(yùn)和輻射輸運(yùn)等物理過程用MonteCarlo來解決實(shí)際問題卻始于本世紀(jì)40年代。Metropolis等人在這一時(shí)期的工作主要就是對中子擴(kuò)散進(jìn)行隨機(jī)抽樣計(jì)算機(jī)模擬，得出所要求算的相關(guān)參數(shù)，并把這種隨機(jī)抽樣方法命名為MonteCarlo方法。

隨著電子計(jì)算機(jī)的迅速發(fā)展，人們開始有意識地、廣泛、系統(tǒng)地應(yīng)用隨機(jī)抽樣方法來解決大量的數(shù)學(xué)、物理和化學(xué)等方面的問題，并且將MonteCarlo方法作為一門獨(dú)立的計(jì)算方法進(jìn)行研究，并隨之向各個學(xué)科領(lǐng)域滲透。MonteCarlo方法MonteCarlo方法解決的問題：1、問題本身是確定性問題，要求我們?nèi)ふ乙粋€隨機(jī)過程，使該隨機(jī)過程的統(tǒng)計(jì)平均是所求問題的解；2、問題本身就是一個隨機(jī)過程，可根據(jù)問題本身的實(shí)際過程來進(jìn)行計(jì)算機(jī)模擬，并采用統(tǒng)計(jì)方法來求得問題的解。1.3MonteCarlo方法簡介基本思想：當(dāng)所求的問題是某種事件出現(xiàn)的概率，或者是某個隨機(jī)變量的期望值時(shí)，它們可以通過某種“隨機(jī)試驗(yàn)”的方法，得到這種事件出現(xiàn)的概率，或者得到這個隨機(jī)變量的統(tǒng)計(jì)平均值，并用它們作為問題的解。MonteCarlo方法1、問題本身是確定性問題MonteCarlo方法圓周率、定積分?2、問題本身就是一個隨機(jī)過程MonteCarlo方法MonteCarlo方法應(yīng)用：

1、數(shù)學(xué)：本身已形成計(jì)算數(shù)學(xué)的一個分支；

2、粒子物理：輸運(yùn)問題、屏蔽問題、核武器試驗(yàn)分析等；

3、統(tǒng)計(jì)物理、化學(xué)，材料、工程各領(lǐng)域；

4、其它：疾病傳播與免疫、系統(tǒng)工程與管理優(yōu)化等等。MonteCarlo方法應(yīng)用此方法求解工程技術(shù)問題可以分為兩類:確定性問題和隨機(jī)性問題。

二、蒙特卡羅方法的應(yīng)用MonteCarlo方法蒙特卡羅模擬流程圖MonteCarlo方法產(chǎn)生隨機(jī)數(shù)，進(jìn)行多次重復(fù)試驗(yàn)。

統(tǒng)計(jì)分析模擬試驗(yàn)結(jié)果，給出問題的概率解以及解的精度估計(jì)。

構(gòu)造一個簡單、適用的概率模型

上圖是葉片的橫截面，葉片外表面與高溫燃?xì)饨咏佑|，高溫燃?xì)獾臏囟韧ㄟ^熱障涂層傳遞到葉片的外表面，葉片內(nèi)表面用冷卻介質(zhì)進(jìn)行冷卻。由于熱障涂層厚度會在一個小范圍內(nèi)波動，求葉片表面溫度變化范圍。熱障涂層熱傳導(dǎo)的MC模擬MonteCarlo方法

上圖是葉片局部放大圖，

Tgas-燃?xì)鉁囟?，TTBC-熱障涂層溫度，TMH-金屬葉片外表面溫度，TMC-金屬葉片內(nèi)表面溫度，Tcool-冷卻介質(zhì)溫度MonteCarlo方法葉片冷卻邊金屬葉片熱障涂層涂層接觸燃?xì)釺gasTTBCTMHTcoolTMC模型建立MonteCarlo方法模型建立MonteCarlo方法求解線性方程組的程序MonteCarlo方法function[Ttbc,Tmh,Tmc,q]=mtkl1(hgas,Tgas,ktbc,Ltbc,km,Lm,hcool,Tcool)%計(jì)算矩陣K=[-hgas, 0, 0, -1;ktbc/Ltbc, -ktbc/Ltbc, 0, -1;0, km/Lm, -km/Lm, -1;0, 0, hcool, -1;];%計(jì)算等式右側(cè)量b=[-hgas*Tgas;0;0;hcool*Tcool];u=K\b;%輸出計(jì)算結(jié)果Ttbc=u(1);Tmh=u(2);Tmc=u(3);q=u(4);模型求解MonteCarlo方法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

涂層厚度單參數(shù)均勻分布的蒙特卡羅法模擬MonteCarlo方法其概率密度函數(shù)如右圖涂層厚度單參數(shù)均勻分布的蒙特卡羅法模擬程序MonteCarlo方法clearall;%參數(shù)的標(biāo)稱值hgas=3000; %TBC-氣體熱傳遞系數(shù)(W/m^2)Tgas=1300; %混合氣體溫度(c)ktbc=1;%TBC熱傳導(dǎo)。(W/mK)km=21.5%金屬溫度傳導(dǎo)。（W/mK）Lm=0.003%金屬厚度（m）hcool=1000; %冷卻液—金屬熱傳遞系數(shù)(W/m^2)Tcool= 200; %冷卻液溫度(c)%蒙特卡洛試驗(yàn)次數(shù)Ntrial=1000;forn=1:Ntrial,Ltbc(n)=0.00025+0.0005*rand;%用均勻分布產(chǎn)生Ltbc的值。

[Ttbc,Tmh(n),Tmc,q]=mtkl1(hgas,Tgas,ktbc,Ltbc(n),km,Lm,hcool,Tcool);if(n>1),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');涂層厚度與葉片外表面溫度在不同區(qū)間分布數(shù)MonteCarlo方法模擬實(shí)驗(yàn)次數(shù)與葉片外表面溫度及方差的關(guān)系MonteCarlo方法

隨著模擬次數(shù)增大，葉片外表面溫度偏離程度減少。多參數(shù)均勻分布的蒙特卡羅法模擬MonteCarlo方法

對于多個輸入的情況，和單變量類似，只是需要對每個輸入都產(chǎn)生一系列隨機(jī)數(shù)。問題：隨著隨機(jī)輸入的個數(shù)增加，TMH的分布將會發(fā)生怎樣的變化？

多參數(shù)均勻分布的蒙特卡羅法模擬程序MonteCarlo方法clearall;%輸入?yún)?shù)的上下界hgas=[1500,4500]; %TBC-氣體熱傳遞系數(shù)(W/m^2)Tgas=[1200,1400]; %混合氣體溫度(C)ktbc =[0.9,1.1]; %TBC熱傳導(dǎo)。(W/mK)Ltbc =[0.00025,0.00075];%TBC厚度(m)km =[20.0,23.0]; %金屬溫度傳導(dǎo)。（W/mK）Lm =[0.002, 0.004]; %金屬厚度（m）hcool=[500,1500]; %冷卻液—金屬熱傳遞系數(shù)(W/m^2)Tcool=[150,250]; %冷卻液溫度(C)%將參數(shù)邊界放進(jìn)一個向量中Pbound=[hgas;Tgas;ktbc;Ltbc;km;Lm;hcool;Tcool];Ntrial=1000;forn=1:Ntrial,%用均勻分布產(chǎn)生各個參數(shù)的值。P(:,n)=Pbound(:,1)+(Pbound(:,2)-Pbound(:,1)).*rand;%求解熱傳導(dǎo)問題[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(n>1),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');葉片外表面溫度在不同區(qū)間分布數(shù)模擬實(shí)驗(yàn)次數(shù)與葉片外表面溫度及方差的關(guān)系MonteCarlo方法多參數(shù)非均勻分布的蒙特卡羅法模擬MonteCarlo方法

怎樣處理非均勻分布的問題

在MATLAB（和很多模擬軟件包）中