小论文-Ti-Si-N表面生长KMC并行仿真.doc

Ti-Si-N薄膜仿真田伟伟（内蒙古科技大学 机械工程学院 05机械2班）摘要：本文从薄膜生长理论，仿真的方法，KMC程序实现等几个方面论述了Ti-Si-N薄膜仿真技术，然后对KMC仿真程序的结果作出总结。关键词：薄膜生长理论；计算机模拟；KMCAbstract：This article provides of Ti-Si-N film mulation technology ，Form theory of thin film growth，simulation method and KMC procedures to achieve，And KMC simulation program on the summary.Key Words: Theory of thin film growth; Simulation;KMC;随制造业的快速发展，一次性设计必须在越来越短的时间内进行重新设计。数字工厂彻底改变传统设计方式，从概念的提出，技术支持到工艺的安排及工艺参数的选取等一系列的过程在现实实现之前得到虚拟环境中模拟、优化。数字工厂可以制造出质量更好的产品，还降低了成本，缩短了时间，提高了效率。薄膜工艺在现代社会有着非常重要的地位，随着现代微电子及光电子工业向着集成化和微型化的发展趋势，预计21世纪分子/电子器件尺寸将是在nm量级。这些材料都是在非平衡态下通过薄膜生长获得的。Ti-Si-N纳米复合薄膜以其硬度高、抗高温氧化性能好、摩擦系数小、弹性模量高、和基体的结合力强、热稳定性优良等性能成为超硬材料研究领域里的研究热点。实验研究还未能给出Ti-Si-N薄膜制备的可靠工艺；理论研究未能给出超硬机理。利用计算机对真实的物质系统进行模拟，并通过模拟结果指导新材料研究与开发是材料设计的有效方法之一。1 薄膜生长理论薄膜生长在原子尺度主要表现为粒子的迁移与扩散，粒子如何迁移有两个重要的概念：过渡态和势能面1.1过渡态理论过渡态理论是1935年由艾林(Eyring)和波兰尼(Polany)等人在统计热力学和量子力学的基础上提出来的。他们认为由反应物分子变成生成物分子，中间一定要经过一个过渡态，而形成这个过渡态必须吸取一定的活化能，这个过渡态就称为活化络合物，所以又称为活化络合物理论。用该理论，只要知道分子的振动频率、质量、核间距等基本物性，就能计算反应的速率系数，所以又称为绝对反应速率理论(absolute rate theory)。该理论认为反应物分子间相互作用的势能是分子间相对位置的函数，活化络合物与反应物之间建立化学平衡，总反应速率由活化络合物转化成产物的概率决定。过渡态理论是在统计热力学和量子力学的基础上提出来的。 系统从稳定态到稳定态中间必须经过一个过渡态，而形成这个过渡态必须吸取一定的活化能。系统相互作用是原子间相对位置的函数，反应速率由过渡态转化为稳定态的概率决定。1.2 势能面 过渡态理论实现的主要难点在于势能面的计算及其反应路径的发现。蒙特卡罗方法计算系统构型能量时更为关注形成过渡态的激活能，具体实现上根据使用模型不同实现不同。在无格点KMC仿真中为在势能面上寻找鞍点的位置。2 仿真方法原子尺度模拟计算，所用的方法主要是分子动力学和蒙特卡罗方法。2.1 分子动力学所谓分子动力学，是指对原子核和电子所构成系统中，使用经典牛顿物理学，求解微粒的运动情况，采用Verlet或其它各种数值算法得到系统的解，取统计结果得到系统宏观的结构和性质。分子动力学一般方法是：1，先设定模型。例如在一个分子系统中，假定两个分之间的相互作用势为硬球势，确定其势函数。2，给定初始条件。分子动力学模型的过程进入对系统微分方程组成做数值求解时，需要知道粒子的初始位置和速度数值。不同的算法要求不同的初始条件。3，趋于平衡。按照给出的运动方程、边界条件和初始条件，就可以进行分子动力学模拟计算。但是，这样计算出的系统不会具有所要求的系统能量，并且这个状态本身也还不是一个平衡态。为了使系统达到平衡，模拟中需要一个趋衡过程。在这个过程中，我们增加或从系统中移出能量，直到系统具有所要求的能量。然后，再对运动方程中的时间向前积分若干步，使系统持续给出确定能量值。我们称这时系统已经达到平衡态。这段达到平衡所需的时间称为弛豫时间。4，宏观物理量的计算。实际计算宏观物理量是在分子动力学模拟的最后阶段进行的。它是沿着相空间轨迹求平均来计算得到的。分子动力学计算的基本思想是赋予分子体系初始运动状态之后利用粒子在力的作用下运动，然后在相空间中抽取样本进行统计计算，时间步长就是抽样的间隔，因而时间步长的选取对动力学模拟非常重要。太长的时间步长会造成分子间的激烈碰撞，体系数据溢出；太短的时间步长会降低模拟过程搜索相空间的能力，因此一般选取的时间步长为体系各个自由度中最短运动周期的十分之一。但是通常情况下，体系各自由度中运动周期最短的是各个化学键的振动，而这种运动对计算某些宏观性质并不产生影响。化学键的振动时间尺度为秒。很明显如此短的时间步长要模拟粒子的沉积行为(秒)，是非常困难的。因此时间步长就成为了分子动力学模拟的计算瓶颈。2.2 蒙特卡罗模型根据带球问题的变化规律，人为地够找一个合适的概率模型，依照该模型进行大量的统计试验，使它的某些统计参量正好是待求问题的解，这种求解方法就是Monte Carlo方法。当所要求解的问题是某种事件出现的概率，或者是某个随机变量的期望值时，它们可以通过对模型或者过程的观察或抽样试验，得到这种事件出现的频率，或者是某个随机变量的期望值时，它们可以通过模型或者过程的观察或者抽样试验，得到这种事件出现的频率，或者这个随机变数的平均值，并用它们作为问题的解。这就是蒙特卡洛的基本思想。蒙特卡洛方法通过数学方法来加以模拟，即进行一种模拟实验。它是以一个概率模型为基础，按照这个模型所描绘的过程，通过模拟实验的结果，作为问题的近似解。用Monte Carlo 方法处理问题可以分为俩类，一类是所求解问题本身具有概率和统计性，所遇随机性问题，通常采用直接模拟方法。首先必须根据实际问题的随机变量y=g（x1,x2,x3,xm）的分布。蒙特卡洛方法的计算过程实际上就是用数学方法模拟实际过程，它主要是在计算机上产生已知分布的随机变量样本，以代替昂贵的甚至难以实现的试验，所以它又被看作是用计算机来完成实验的一种方法。另一类是确定性问题，通常采用间接模拟方法。在解决确定性问题时，首先要建立一个有关的概率统计模型，是所求得解就是这个模型的概率分布或数学期望，然后对这个模型做随机抽样，最后用其算术平均值作为所求解的近似值。2.3 KMC方法动力学蒙特卡罗方法是将动力学和蒙特卡罗方法结合，它是一个综合模型，已成为原子尺度研究薄膜生长最有利工具。过渡态理论认为例子从一个位置迁移到另一个位置上时，中间一定要经过一个过渡态，并且形成这个过渡态必须吸收一定的能量。由概率论知道粒子迁移的速度可由粒子形成过渡态的概率决定。若给定的系统的态A，如果已知系统从态i转变到其他n个态的转变速率（ki-j），那么形同发生态转变所需要的平均时间为 粒子迁移的概率可有公式近似决定 P= v0exp(E/kT)其中v0约为粒子振动频率，单位为10-12；E为激活能；k为玻尔兹曼常数；T为温度上面来个方程组成了KMC方法的基础，KMC方法描述如下：（1）假定系统起始态为i，经过时间，系统按方程给出的几率随机地转化到态j，2）重新计算系统从态j转变到所有其它m个新态的数率，重复过程（1）.如果需要研究的时间尺度相对于时间爱你步长来说并不是足够长，时间步长取为，其中r是0，1之间的随机数。KMC方法的框架图如下：设置连续原子沉积时间t沉积原子数目Num是否下落新粒子更新每一个事件的扩散速率Ra计算所有可以执行事件的扩散速率Ra系统可执行事件的速率和物理时间进度沉积原子时间到tNmax仿真结果输出及数据统计图YesNoYesNoNo图2.1动力学蒙特卡罗仿真流程图3 KMC程序具体的实现本程序主要分为三部分，用户界面及主程序、data类、dimer类。主程序是用户控制程序；data类主要用于数据的读取与存储及搜索粒子信息；dimer类是dimer方法的实现，用于计算粒子的反应速率3。MainFormDataDimer更新读取图 3.1 仿真程序主框架模型实施具体步骤本课题所建模型实施步骤如下：1、通过控制界面设定沉积参数：沉积温度、沉积速率、入射Si含量、沉积粒子数；2、设定仿真参数：基底长和宽、确定相互粒子作用半径、影响粒子作用参数alpha；3、程序开始计算，首先产生基底粒子。实例化data类，开通粒子信息存储通道；4、随机产生新粒子；粒子x,y坐标由随机数产生，粒子z坐标由周围“最高”粒子决定，然后通过dimer方法计算找到粒子稳定位置；5、对新“下落”粒子周围进行驰豫。驰豫时认为粒子活动收敛精度为0.02。另外，对新粒子周围（半径为5的范围内）可能的迁移路径重新计算，更新其迁移概率；6、用马尔可夫链方法随机选取粒子及其迁移方向。7、计算迁移消耗时间，判断粒子跃迁可能性。消耗时间由系统迁移概率决定。8、若迁移时间大于沉积间隔时间，系统不跃迁；继续下落新粒子回到步骤4；9、如果迁移时间小于间隔时间，选中粒子在所选方向迁移。对迁出位置及迁入位置进行驰豫，并对其迁移概率重新计算。10、继续选择粒子，对其迁移时间计算，直至系统迁移时间大于沉积间隔时间；11、继续下落新粒子，直至粒子数大于设定沉积粒子数；12、提取过程中存储的数据，对粒子信息进行显示。其仿真程序流程图见图3.2图3.2 仿真程序框图4 仿真程序对模拟的验证粒子位置表4.1计算粒子位置点与格点位置比较typexyz格点x格点y格点z偏离113.7801159418.01997829.5412149913.7818.029.540.00122070409.5395480353.1810955139.5408255429.543.189.540.00144427911.0594050115.2999945549.5411507751.065.39.540.001295502075.2545004215.8998318115.9027167575.2615.915.90.006136318018.017658427.42425978218.0219057318.027.4218.020.005221163118.019259645.30284869818.0234527418.025.318.020.004537028020.139478615.30403421818.0237026420.145.318.020.005500573011.6604688739041905111.663.1815.90.004629165132.8536460120.142006918.0172472732.8620.1418.020.007209602平均粒子偏离格点位置0.0037194334，这和Ti/N的晶格常数2.12相比，是其 0.175倍；再加之粒子本身振动行为和不规整的沉积形貌原因，得到此数据是可信的。例子下落结构形貌 图4.1 温度为600K时 图4.2 温度为500K时 图4.3 温度为400K不同的温度形貌对比图可以看出，Ti/N随着温度的升高，其生长模式的变化。温度低时，薄膜表面较为粗糙。仿真过程中还发现Ti的迁移率大于N的迁移率，一般为3到4倍。9文献：1 范秋林. 宋力昕. 施尔畏. 胡行方. 纳米复合超硬薄膜的研究现状J. 无机材料学报. 2004; 19(1): 9162 D. Raabe. 计算材料学 北京：化学工业出版社.20023 孙士阳，Ti-Si-N表面生长过程的无格点KMC仿真及相应并行计算研究