计算建模中的数学物理方法5

1、计算建模中的数学物理方法及其应用 北京大学 胡俊峰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 | 无自环的无向图GG中的团： Gibbs Random Fields对于给定的点集S和邻居系统NGibbs分布 T通常为1Z为常数，让分布概率和为1， Vc为团c出现的可能性度量马尔科夫 随机场 与 Gibbs分布MRF刻画的是局部依赖性，Gibbs则描述了一个全局的分布马尔科夫-吉布斯等式（在同一邻域系统下）一个MRF中随机变量的联合分布为Gibbs分布一个Gibbs分布满足MRF的两个条件对MRF的定义转为对能量函数U的定义定义了U，联

5、合分布便确定最大概率等价于最小能量举一个简单的应用（by 徐源盛）图像前背景分割手动标注某些像素为前景，某些像素为背景自动把前景图像抠出图像前背景分割给每个像素一个标号，0表示背景，1表示前景某些像素已经有初始标号，不能改变目标：确定所有像素的标号假设：一个点的标号选取只依赖于与它相邻的四个点（标号和色差）图像前背景分割MRF每个像素为一个随机变量，取值范围0,1 假定象素属于前景值为1 属于背景值为0定义每个象素的邻域系统四领域目标：找一个概率最大的联合分布图像前背景分割进一步定义Gibbs分布的能量根据邻域系统假设V2与i,i的位置无关，且对称定义像素自身能量 E1 E1刻画了每个像素点的

6、概率分布： 一个像素越像前景E1(x=1)越小 K为聚类结果 定义像素间 E2E2刻画对两个随机变量取值的可能性一般来说相邻两个随机变量取值相同 相邻两个随机变量取同样值，能量为0Cij为RGB空间中两个像素的曼哈顿距离平方网络流算法求解最小能量构图一个S到T的割对应一个能量最小割效果：算法分析设定系统为2态关系为邻接2元关系能够求得全局最优解能量函数定义过分严格，不能很好的支持主题分类问题78金属中氦泡的形成：KMC模拟薛建明北京大学物理学院2012-12-179主要内容问题背景KMC算法介绍模拟He在BCC铁中迁移和聚集行为总结80One cup of water = 20 gallons

7、 of oil聚变核反应From Dr.Finkle一、背景燃料充足！为什么是聚变能？81一、背景核聚变最重要问题：材料磁约束 （ITER）惯性约束(ICF/NIF)高强度中子辐照，导致目前所有材料都不满足未来商用聚变堆材料的应用要求（Nature materials，VOL 5，2006） ITERICF82材料关键问题之一：He bubble一、背景负面影响：疲劳断裂机械强度（n,He）嬗变反应He原子集积成泡物理机制？理想材料的构能？最关键问题83由于受到实验手段的限制，目前实验是难以实施的。计算机模拟技术显示出了独特的优越性，在极端环境下新材料研究中占据了非常重要的地位。 一、背景如何

8、解决上述问题？84多尺度计算机模拟技术一、背景multi-scale materials modeling不同模拟方法适合模拟不同规模的系统85KMC介绍二、算法介绍通过构造随机过程研究体系演化的方法即为动理蒙特卡洛方法(kinetic Monte Carlo, KMC)喝醉的鱼，不停随机跳动，每次跳动高度不同86问题：开始鱼均匀分布在各处，经过非常多的中间态，最后鱼到哪里？分析：各种状态；随机跳动；元素（鱼）有效跳动：越过大坝有效，其它忽略跳动高度不一，跳过的概率也就不同。其它事件：消失、出生、聚集（渔网）。KMC介绍二、算法介绍87不关注与体系穿越势垒无关的无效微小振动，而只着眼于体系的组

9、态变化。因此，不能描绘原子的运动轨迹，但是作为体系演化，其“组态轨迹”仍然是正确的。状态： 左（数量） 右（数量） 左（1） 右（9）左（2） 右（8）左（3） 右（7）。左（9） 右（1）状态1状态2状态3KMC特点二、算法介绍88组态变化的时间间隔很长，体系完成的连续两次演化是独立的，无记忆的，所以这个过程是一种典型的马尔可夫过程(Markov process)，即体系从组态到组态， 这一过程只与其跃迁速率有关。如果精确地知道随机过程的所有可能性，我们便可以构造一个随机过程，使得体系按照正确（结果）的轨迹演化。开始 - 状态1 状态2 状态3 - 结果KMC特点算法本身简单、直观，因此应用

10、非常的广泛！二、算法介绍892.2 KMC算法：直接法算法主要思想：设体系处于状态t，有n个跃迁事件可能发生，将每条跃迁途径想象成长度与跃迁速率成正比的线段。将这些线段首尾相连，构造随机选择的样本空间。R=10R=4R=5R=11R=910141930390E1E2E3E4E5E1:事件1对应的非归一化概率10E2:事件1对应的非归一化概率4。E5:事件1对应的非归一化概率9902.2 KMC算法：直接法如果产生的随机数落在某个线段中，这个线段所代表的跃迁途径就被选中。 R=10R=4R=5R=11R=910141930390E1E2E3E4E5在0,39区间随即均匀抽样R，R如果在事件n对应

11、的范围内，则事件n发生。912.2 KMC算法：直接法主要步骤如下：计算体系处于不同状态时的总跃迁速率； 选择在0,Rtotal区间的随机数； 在事件概率对应的线性排列中，查找上述随机数对应的时间E； 上述事件E发生，体系移动到下一个状态；模拟时间前进； 重复上述过程。 92一般方法：首先计算跃迁概率总和T，确定抽样范围：在0，1范围抽样r，然后:R=r*T如果 0R10 则事件E1发生如果 10R14则事件E2发生如果 14R19则事件E3发生。R=10R=4R=5R=11R=910141930390E1E2E3E4E5选择发生的事件93跃迁几率计算TST理论当体系处于稳定状态时，我们可以将

12、其描述为处于3N维势能函数面的一个局域极小值（阱底）处。有限温度下，虽然体系内的原子不停的进行热运动，但是绝大部分时间内原子都是在势能阱底附近振动。偶然情况下体系会越过不同势阱间的势垒从而完成一次“演化”，这类小概率事件才是决定体系演化的重点。94跳跃几率D计算跃迁速率：过渡态理论(TST) 其中Em 为迁移活化能Em , f为热振动频率, 一般取1014 s- 1.T为系统温度，Kb是Boltzmann常数。跃迁几率计算Au(100)表面的势垒分布95KMC适用范围非常广泛，包括经济、金融。KMC在物理上的其它应用 雪花 外延生长 ？96三、 金属中He行为物理模型He原子在样品中随机产生，

13、在热运动驱动下沿不同的方向迁移，不同方向的迁移能不同，不同He原子遇到之后结合长大，并不再移动。最后He泡的大小？分布？He的迁移路径？影响因素：He浓度温度迁移能。物理问题色彩代表势能大小97简化物理过程He原子迁移遇到He bubble结合迁移能相同He atom假定He bubble 中He原子数大于1后，静止不再移动;He可以看做只有一个原子的bubble;分析：各种状态；随机跳动(热运动，动能正态)；有效跳动：MST理论KMC方法三、 金属中He行为物理模型983.1 物理模型BCC metalMigration ways物理模型N0个He原子初开空间分布He的所有迁移路径不同迁移路

14、径概率T时刻发生事件的抽样结果统计He复合判断：距离小于一个特定值六个迁移路径：六个事件99操作添加He原子：InsertAtom()He/He bubble 合并: AtomMerging()事件对应的跃迁概率:TransPoss()线性化所有事件概率：SortPoss()选择事件：PickupEvent()产生随机数0,1: RandomNumber()原子位置移动: MovingAtom()系统演化时间：ElapsedTime().3.1 程序：操作100程序结构示意图构造样品：包括物理参数的读入以及He原子分布、样品特殊性描述（缺陷、边界等）。进行反应：事件对应的计算和操作，包括合并、

15、分解、迁移，以及由此导致的其它物理性质的变化。结束时间：最长模拟时间3.2 程序程序结构101数据更新事件查找随机数问题3.3关键算法在模拟过程中，He原子数目可能高达百万，每个He原子对应时间可能超过10个；因此，有效数据在108。如何有效的组织数据，如何设计高效算法成为模拟的关键！102模拟系统描述：在BCC结构单晶铁中，中子辐照产生间隙原子（SIA）、空位和He三种缺陷，它们之间的迁移和复合是互相影响还是可以近似独立？研究原因：目前的实验多数是三个因素分别研究，并将获得的结果用来评估实际核设施中材料的性能，但该方法的误差一直是个迷。3.4 模拟结果基于当前的实验研究，一个核能装置中的寿命为100年，真实条件下是？1000年？1年？103HeVV+HeOSIAHeHenVm+He+HeVHenVm-SIAsSIAsVm-SIAsVm+V+I+V+I可能发生事件示意图+V3.4 模拟结果104结果分布3.4 模拟结果105He浓度影响浓度越高