计算材料学与第一性原理、分子动力学、蒙特卡洛计算方法_第一性原理分子动力学计算模拟经典教程

【引言】计算材料学：“定做”材料的高级理论阶段

计算材料学(Computational Materials Science)是近年来飞速发展的一门新兴交叉学科。它综合了凝聚物理、材料物理学、理论化学、材料力学和工程力学、计算机算法等多个相关学科。学科旨在利用现代高速计算机，模拟材料的各种物理化学性质，深入理解材料从微观到宏观多个尺度的各类现象与特征，并对于材料的结构和物性进行理论预言，从而达到设计新材料的目的。与传统的物理、化学、材料等学科相比，计算材料学仍是一门正处在迅猛发展时期的新兴学科。

计算材料学主要包括两个方面的内容：一方面是计算模拟，即从实验数据出发，通过建立数学模型及数值计算，模拟实际过程；另一方面是材料的计算机设计，即直接通过理论模型和计算，预测或设计材料结构与性能。前者使材料研究不是停留在实验结果和定性的讨论上，而是使特定材料体系的实验结果上升为一般的、定量的理论，后者则使材料的研究与开发更具方向性、前瞻性，有助于原始性创新，可以大大提高研究效率。因此，计算材料学是连接材料学理论与实验的桥梁。

目前常用的计算方法包括第一性原理从头计算法，分子动力学方法，蒙特卡洛方法，元胞自动机方法、相场法、几何拓扑模型方法、有限元分析等。

以GPU 为代表的高性能计算技术有效提高了计算机的模拟能力，结合算法以及理论，在计算机虚拟环境下从纳观、微观、介观、宏观尺度对材料进行多层次研究，也可以模拟超高温、超高压等极端环境下的材料服役性能，模拟材料在服役条件下的性能演变规律、失效机理，进而实现材料服役性能的改善和材料设计，有效减少了在优化材料和设计新工艺方面所必须进行的大量试验。材料模拟和材料制备工艺大幅进步，极大地促进了新产品的优化和开发。

下面就对三大经典的、常见材料计算方法：第一性原理、分子动力学和蒙特卡洛方法进行简单介绍。

【人类群星闪耀时】第一性原理：从头算，算到尾，算到宇宙尽头

第一性原理（First Principle），是从量子力学理论出发的计算方法，它“号称”自己仅需要原子精细结构常数、电子质量及带电量、原子核质量及电量、普朗克常量和光速这几个已知的参数，便可根据原子核和电子相互作用的原理及其基本运动规律，经过多个近似处理后直接求解薛定谔方程，进而得到材料（几乎）所有的基态性质。狭义的第一性原理计算，是指基于Hatree-Fock自洽场计算方法的“从头算”（ab initio），广义的第一性原理计算在此基础上还包含了密度泛函理论（DFT）计算。但也有人提出，希望“从头算”专指HF方法，而第一性原理专指DFT方法。关于HF方法和DFT方法更具体的简介，可查看本专栏的第一篇文章。

在材料计算中，第一性原理方法常见于超晶胞、表面、界面、团簇等晶体结构模型，用于优化其几何结构、计算体系能量，得到如能带结构、态密度、电荷分布等重要的信息，进而对体系进行理论量化的分析。理论上说，在有第一性原理要求的几个参数之后，的确可以通过量子力学“从头算到尾”，但在实际操作中，除了第一性原理要求的已知参数外，我们还常常加入一些“经验参数”，它们通常来自于第一性原理计算中前人得出的、已得到大量实例验证的规律，或者来自于实验领域最直接的结果，这样的做法通常称之为“半经验的”，它能有效地减少计算资源的损耗，进而保证计算工作在最优化的条件下进行。

当然，“第一性原理”这个词现在远远不仅用于计算领域，著名的大牛马一龙非常吹捧“第一性原理”的思维方式，还曾有斯坦福的教授在清华大学的讲座中提出要用“第一性原理”进行投资。Anyway，这种方法的思想是“自上而下”的，具有人类对自然认知最原始的艺术美感，对日常生活也具有一定的启发性。

【黄金时代的余晖】分子动力学：微观世界的经典力学方法

分子动力学（Molecular Dynamics，简称MD），是从经典物理的统计力学出发的计算方法，它通过对分子间相互作用势函数及运动方程的求解，分析其分子运动的行为规律，模拟体系的动力学演化过程，给出微观量（如：分子的坐标与速度等）与宏观可观测量（如：体系的温度、压强、热容等）之间的关系，从而研究复合体系的平衡态性质和力学性质，是研究材料内部流体行为、通道运输等现象有效的研究手段。

分子动力学的基本思想是：根据玻恩-奥本海默近似，可以将原子核与电子分离开，再将原子核假设为组成体系的经典粒子进行研究。按照经典力学牛顿定律，体系中的经典粒子受力为：

其中，F为粒子所受的力，a为粒子的加速度，而m是粒子的质量。我们知道，力可以表示为势能函数对坐标的一阶导数，加速度可以表示为速度对坐标的一阶导数，我们将这组方程改写为以下两组方程：

其中，v是粒子的速度，U为相应的势能函数。显然，解这两组方程是非常繁杂的，我们通常会采用数值解进行求解，因为严格的解析解只可能出现在简单的势函数作用的体系之下，这在实际的研究中是几乎不可能出现的。有理论力学基础的同学应该对这个解法不陌生，它与解相空间的哈密顿方程有点相似。这样一个随时间演化的体系, 我们可以根据坐标描述其运动轨迹，根据速度描述其运动轨迹上的变化，从而完全描述其演化过程，这是所谓的经典力学的“确定性”。通过给定势函数，赋予体系初始的坐标和速度，我们会得到一系列包含了整个分子动力学过程的坐标与速度，再通过对坐标与速度的统计，我们可以得到需要的体系相应的热力学与动力学性质，分子动力学流程草图如下：

通过方程可以看出，势函数（也称之为“力场”）的选取对分子动力学的结果非常重要，加之分子动力学对计算资源的消耗非常大，一个常见的分子动力学程序可能往往要跑一周以上，因此合理选取势函数、正确使用算法也成了分子动力学实际操作中的重点。对于势函数的选取，通常选用一些经典势（如：伦纳德-琼斯势，学过固体物理的同学肯定不会陌生）。当然，势函数也可以通过第一性原理计算，通过电荷密度得到Hellmann-Feynman力得到，这种方法又被称之为“第一性原理分子动力学”，目前应用也非常广泛，此处不再展开。

【掷骰子吧！】蒙特卡洛：科学研究不可以“随意”，但可以“随机”

蒙特卡洛方法（Monte Carlo，简称MC），是以概率论和数理统计为理论基础、使用随机数（或更常见的伪随机数）解决实际问题的一种随机抽样统计方法，它常用于求解一些带有“随机”性质的实际生活问题和研究一些带有现条件下难以观测的物理量的实验。根据大数定律，要使得随机的独立事件呈现出具有一定规律的统计结果，需要大量地进行重复实验，例如十九世纪末蒲丰就曾提出用“投针实验”的统计结果来估算圆周率，这个实验大概需要40万次重复，才能使得其置信水平达到0.95以上。因此，在计算机上运用蒙特卡洛方法进行模拟实验，本身就具有天然的优势。实际上，在冯诺依曼等人研究核武器时，就曾运用蒙特卡洛方法来分析中子在核反应堆中的传输过程。由于量子力学不确定性原理的制约，科学家只能通过随机抽样、模拟超大数量的中子行为，并使用其统计结果作为设计核武器的依据。相传，“蒙特卡洛方法”这个名字正是冯诺依曼起的，他觉得用这座世界闻名的赌城来命名这个方法，更加具有一丝神秘的气息。

蒙特卡洛方法与前文说到的分子动力学方法既有相似之处，也有显著的不同。首先，读者应该可以体会到，这两种方法都蕴含着统计的思想，通过体系的统计结果得到相应的体系性质。不同之处在于，分子动力学方法中每一个粒子的运动，取决于它所处的势函数，严格按照牛顿定律在相空间进行演化；而蒙特卡洛方法中每一个粒子的运动，则是取决于抽样所给定的概率分布，例如：给定一个麦克斯韦-玻尔兹曼分布，每个粒子的运动状态就完全随机从中抽取得到。也就是说，分子动力学中的统计结果来自于“确定的”经典力学结果，而蒙特卡洛方法的统计结果来自于“随机的”概率统计结果。此外，蒙特卡洛方法相比前两者，还具有程序相对简单、消耗计算资源较小等优点，有时候，也许还会有一些常规思维意想不到的“意外发现”。