計算建模中的數(shù)學(xué)物理方法5

1、計算建模中的數(shù)學(xué)物理方法及其應(yīng)用 北京大學(xué) 胡俊峰2013-1-29先來看一道期末考題：問題的直觀描述Distribution of identical billiard ballsRegion ofinterestCenter ofmassMean ShiftvectorObjective : Find the densest regionMean Shift Distribution of identical billiard ballsRegion ofinterestCenter ofmassMean ShiftvectorObjective : Find the densest re

2、gionMean Shift Distribution of identical billiard ballsRegion ofinterestCenter ofmassMean ShiftvectorObjective : Find the densest regionMean ShiftDistribution of identical billiard ballsRegion ofinterestCenter ofmassMean ShiftvectorObjective : Find the densest regionMean ShiftDistribution of identic

3、al billiard ballsRegion ofinterestCenter ofmassMean ShiftvectorObjective : Find the densest regionMean ShiftDistribution of identical billiard ballsRegion ofinterestCenter ofmassMean ShiftvectorObjective : Find the densest regionMean ShiftDistribution of identical billiard ballsRegion ofinterestCent

4、er ofmassObjective : Find the densest region算法的定性表達：St = xi; where | mt - xi | 無自環(huán)的無向圖GG中的團： Gibbs Random Fields對于給定的點集S和鄰居系統(tǒng)NGibbs分布 T通常為1Z為常數(shù)，讓分布概率和為1， Vc為團c出現(xiàn)的可能性度量馬爾科夫 隨機場 與 Gibbs分布MRF刻畫的是局部依賴性，Gibbs則描述了一個全局的分布馬爾科夫-吉布斯等式（在同一鄰域系統(tǒng)下）一個MRF中隨機變量的聯(lián)合分布為Gibbs分布一個Gibbs分布滿足MRF的兩個條件對MRF的定義轉(zhuǎn)為對能量函數(shù)U的定義定義了U，聯(lián)

5、合分布便確定最大概率等價于最小能量舉一個簡單的應(yīng)用（by 徐源盛）圖像前背景分割手動標注某些像素為前景，某些像素為背景自動把前景圖像摳出圖像前背景分割給每個像素一個標號，0表示背景，1表示前景某些像素已經(jīng)有初始標號，不能改變目標：確定所有像素的標號假設(shè)：一個點的標號選取只依賴于與它相鄰的四個點（標號和色差）圖像前背景分割MRF每個像素為一個隨機變量，取值范圍0,1 假定象素屬于前景值為1 屬于背景值為0定義每個象素的鄰域系統(tǒng)四領(lǐng)域目標：找一個概率最大的聯(lián)合分布圖像前背景分割進一步定義Gibbs分布的能量根據(jù)鄰域系統(tǒng)假設(shè)V2與i,i的位置無關(guān)，且對稱定義像素自身能量 E1 E1刻畫了每個像素點的

6、概率分布： 一個像素越像前景E1(x=1)越小 K為聚類結(jié)果 定義像素間 E2E2刻畫對兩個隨機變量取值的可能性一般來說相鄰兩個隨機變量取值相同 相鄰兩個隨機變量取同樣值，能量為0Cij為RGB空間中兩個像素的曼哈頓距離平方網(wǎng)絡(luò)流算法求解最小能量構(gòu)圖一個S到T的割對應(yīng)一個能量最小割效果：算法分析設(shè)定系統(tǒng)為2態(tài)關(guān)系為鄰接2元關(guān)系能夠求得全局最優(yōu)解能量函數(shù)定義過分嚴格，不能很好的支持主題分類問題78金屬中氦泡的形成：KMC模擬薛建明北京大學(xué)物理學(xué)院2012-12-179主要內(nèi)容問題背景KMC算法介紹模擬He在BCC鐵中遷移和聚集行為總結(jié)80One cup of water = 20 gallons

7、 of oil聚變核反應(yīng)From Dr.Finkle一、背景燃料充足！為什么是聚變能？81一、背景核聚變最重要問題：材料磁約束 （ITER）慣性約束(ICF/NIF)高強度中子輻照，導(dǎo)致目前所有材料都不滿足未來商用聚變堆材料的應(yīng)用要求（Nature materials，VOL 5，2006） ITERICF82材料關(guān)鍵問題之一：He bubble一、背景負面影響：疲勞斷裂機械強度（n,He）嬗變反應(yīng)He原子集積成泡物理機制？理想材料的構(gòu)能？最關(guān)鍵問題83由于受到實驗手段的限制，目前實驗是難以實施的。計算機模擬技術(shù)顯示出了獨特的優(yōu)越性，在極端環(huán)境下新材料研究中占據(jù)了非常重要的地位。 一、背景如何

8、解決上述問題？84多尺度計算機模擬技術(shù)一、背景multi-scale materials modeling不同模擬方法適合模擬不同規(guī)模的系統(tǒng)85KMC介紹二、算法介紹通過構(gòu)造隨機過程研究體系演化的方法即為動理蒙特卡洛方法(kinetic Monte Carlo, KMC)喝醉的魚，不停隨機跳動，每次跳動高度不同86問題：開始魚均勻分布在各處，經(jīng)過非常多的中間態(tài)，最后魚到哪里？分析：各種狀態(tài)；隨機跳動；元素（魚）有效跳動：越過大壩有效，其它忽略跳動高度不一，跳過的概率也就不同。其它事件：消失、出生、聚集（漁網(wǎng)）。KMC介紹二、算法介紹87不關(guān)注與體系穿越勢壘無關(guān)的無效微小振動，而只著眼于體系的組

9、態(tài)變化。因此，不能描繪原子的運動軌跡，但是作為體系演化，其“組態(tài)軌跡”仍然是正確的。狀態(tài)： 左（數(shù)量） 右（數(shù)量） 左（1） 右（9）左（2） 右（8）左（3） 右（7）。左（9） 右（1）狀態(tài)1狀態(tài)2狀態(tài)3KMC特點二、算法介紹88組態(tài)變化的時間間隔很長，體系完成的連續(xù)兩次演化是獨立的，無記憶的，所以這個過程是一種典型的馬爾可夫過程(Markov process)，即體系從組態(tài)到組態(tài)， 這一過程只與其躍遷速率有關(guān)。如果精確地知道隨機過程的所有可能性，我們便可以構(gòu)造一個隨機過程，使得體系按照正確（結(jié)果）的軌跡演化。開始 - 狀態(tài)1 狀態(tài)2 狀態(tài)3 - 結(jié)果KMC特點算法本身簡單、直觀，因此應(yīng)用

10、非常的廣泛！二、算法介紹892.2 KMC算法：直接法算法主要思想：設(shè)體系處于狀態(tài)t，有n個躍遷事件可能發(fā)生，將每條躍遷途徑想象成長度與躍遷速率成正比的線段。將這些線段首尾相連，構(gòu)造隨機選擇的樣本空間。R=10R=4R=5R=11R=910141930390E1E2E3E4E5E1:事件1對應(yīng)的非歸一化概率10E2:事件1對應(yīng)的非歸一化概率4。E5:事件1對應(yīng)的非歸一化概率9902.2 KMC算法：直接法如果產(chǎn)生的隨機數(shù)落在某個線段中，這個線段所代表的躍遷途徑就被選中。 R=10R=4R=5R=11R=910141930390E1E2E3E4E5在0,39區(qū)間隨即均勻抽樣R，R如果在事件n對應(yīng)

11、的范圍內(nèi)，則事件n發(fā)生。912.2 KMC算法：直接法主要步驟如下：計算體系處于不同狀態(tài)時的總躍遷速率； 選擇在0,Rtotal區(qū)間的隨機數(shù)； 在事件概率對應(yīng)的線性排列中，查找上述隨機數(shù)對應(yīng)的時間E； 上述事件E發(fā)生，體系移動到下一個狀態(tài)；模擬時間前進； 重復(fù)上述過程。 92一般方法：首先計算躍遷概率總和T，確定抽樣范圍：在0，1范圍抽樣r，然后:R=r*T如果 0R10 則事件E1發(fā)生如果 10R14則事件E2發(fā)生如果 14R19則事件E3發(fā)生。R=10R=4R=5R=11R=910141930390E1E2E3E4E5選擇發(fā)生的事件93躍遷幾率計算TST理論當(dāng)體系處于穩(wěn)定狀態(tài)時，我們可以將

12、其描述為處于3N維勢能函數(shù)面的一個局域極小值（阱底）處。有限溫度下，雖然體系內(nèi)的原子不停的進行熱運動，但是絕大部分時間內(nèi)原子都是在勢能阱底附近振動。偶然情況下體系會越過不同勢阱間的勢壘從而完成一次“演化”，這類小概率事件才是決定體系演化的重點。94跳躍幾率D計算躍遷速率：過渡態(tài)理論(TST) 其中Em 為遷移活化能Em , f為熱振動頻率, 一般取1014 s- 1.T為系統(tǒng)溫度，Kb是Boltzmann常數(shù)。躍遷幾率計算Au(100)表面的勢壘分布95KMC適用范圍非常廣泛，包括經(jīng)濟、金融。KMC在物理上的其它應(yīng)用 雪花 外延生長 ？96三、 金屬中He行為物理模型He原子在樣品中隨機產(chǎn)生，

13、在熱運動驅(qū)動下沿不同的方向遷移，不同方向的遷移能不同，不同He原子遇到之后結(jié)合長大，并不再移動。最后He泡的大小？分布？He的遷移路徑？影響因素：He濃度溫度遷移能。物理問題色彩代表勢能大小97簡化物理過程He原子遷移遇到He bubble結(jié)合遷移能相同He atom假定He bubble 中He原子數(shù)大于1后，靜止不再移動;He可以看做只有一個原子的bubble;分析：各種狀態(tài)；隨機跳動(熱運動，動能正態(tài))；有效跳動：MST理論KMC方法三、 金屬中He行為物理模型983.1 物理模型BCC metalMigration ways物理模型N0個He原子初開空間分布He的所有遷移路徑不同遷移路

14、徑概率T時刻發(fā)生事件的抽樣結(jié)果統(tǒng)計He復(fù)合判斷：距離小于一個特定值六個遷移路徑：六個事件99操作添加He原子：InsertAtom()He/He bubble 合并: AtomMerging()事件對應(yīng)的躍遷概率:TransPoss()線性化所有事件概率：SortPoss()選擇事件：PickupEvent()產(chǎn)生隨機數(shù)0,1: RandomNumber()原子位置移動: MovingAtom()系統(tǒng)演化時間：ElapsedTime().3.1 程序：操作100程序結(jié)構(gòu)示意圖構(gòu)造樣品：包括物理參數(shù)的讀入以及He原子分布、樣品特殊性描述（缺陷、邊界等）。進行反應(yīng)：事件對應(yīng)的計算和操作，包括合并、

15、分解、遷移，以及由此導(dǎo)致的其它物理性質(zhì)的變化。結(jié)束時間：最長模擬時間3.2 程序程序結(jié)構(gòu)101數(shù)據(jù)更新事件查找隨機數(shù)問題3.3關(guān)鍵算法在模擬過程中，He原子數(shù)目可能高達百萬，每個He原子對應(yīng)時間可能超過10個；因此，有效數(shù)據(jù)在108。如何有效的組織數(shù)據(jù)，如何設(shè)計高效算法成為模擬的關(guān)鍵！102模擬系統(tǒng)描述：在BCC結(jié)構(gòu)單晶鐵中，中子輻照產(chǎn)生間隙原子（SIA）、空位和He三種缺陷，它們之間的遷移和復(fù)合是互相影響還是可以近似獨立？研究原因：目前的實驗多數(shù)是三個因素分別研究，并將獲得的結(jié)果用來評估實際核設(shè)施中材料的性能，但該方法的誤差一直是個迷。3.4 模擬結(jié)果基于當(dāng)前的實驗研究，一個核能裝置中的壽命為100年，真實條件下是？1000年？1年？103HeVV+HeOSIAHeHenVm+He+HeVHenVm-SIAsSIAsVm-SIAsVm+V+I+V+I可能發(fā)生事件示意圖+V3.4 模擬結(jié)果104結(jié)果分布3.4 模擬結(jié)果105He濃度影響濃度越高