Efficientstochasticgenerationofspecialquasirandoms

特殊准随机结构的高效随机模拟方法摘要我们提出一种新的算法来生成特殊准随机结构(SQS)，即对处在真无序状态的每个给定原子数量的超胞进行最好的周期性超胞近似。这种方法基于蒙特卡洛模拟法模拟退火循环的过程，获得目标函数，寻求与最大值完美匹配的关联函数（而不是仅仅为了缩小SQS关联函数和无序态关联函数的差距而预设一个关联函数）。该方法优化超晶胞的形状的同时优化了原子的站位，因此能保证获得详尽的空间构型并且没有预设的超胞形状的影响。这种方法已经被应用于AlloyTheoreticAutomatedToolkit(ATAT)软件的“mcsqs”代码中，是计算多组元多亚晶格系统的一般框架，在某种程度上减少了用户需要指定输入信息的数量，提高了并形化效率。关键词：特殊准随机结构、蒙特卡洛、合金原理、第一原理计算方法、无序态、固溶体1.简介因为许多具有可调谐性的技术重要材料以无序合金的形式存在，对无序晶体合金准确和有效的建模在材料科学应用广泛。一个流行、简易、通俗的方法是仅仅构建巨大的超胞和随机化每个晶格点阵的站位。然而，这种方法不是最有效率的，由于同元素原子在会占据临近的晶格点阵，影响局部关联函数，在一个有限的超胞内的随机生成结构有相当大的可能性违背了“完美”的随机性。运用SQS的概念可以系统的改进这种方法。SQS代表了在真无序状态的每个给定原子数量的超胞进行最好的周期性超胞近似。SQS是最符合标准的因为在SQS中相邻点阵站位的特定集合的关联函数可以匹配相应的真正完全无序的状态的关联函数。SQS已经非常成功地用于获得无序材料的电学和热力学性质（见[2-8]）。SQS通常通过完全枚举法，由所有可能的原子站位的超晶胞所构成，这是一个复杂的指数阶的算法。这是唯一已知的保证获得最佳的SQS的方法，但是这只适用于计算易处理的小晶胞（当前最多可以处理25个原子，即使用效率最高的结构枚举算法[9]）。然而，电子结构方法现在可以正常处理相当大的晶胞（超过50个原子），这使得完全枚举法的使用不切实际。所以，现在一般通过随机数法来替代生成SQS[5,10,11]。我们的成果主要有四方面。一、设计了一种新的SQS生成算法，寻找与最大值完美匹配的关联函数的SQS（而不是仅仅为了缩小SQS关联函数和无序态关联函数的差距而预设一个关联函数）。我们归纳了一个通常用于随机SQS研究的被称为“rewards”的目标函数，在误差范围内完美相关。二、我们运用了一种计算多组元多亚晶格系统的一般框架的方法[12]来研究大类技术相关合金（大多数在不同的亚晶上展现出无序性）。三、这种方法同时优化超晶胞的形状和原子的站位，因此可以获得详尽的空间结构并没有预设的超胞形状的影响。最后，我们在ATAT上实施了这个算法[12-14]，在某种程度上减少了用户需要指定输入信息的数量，提高了并行化效率。2.聚簇展开形式让我们先简单介绍一下多组元多亚晶格系统的聚簇展开形式[12]。σ代表变量的矢量，σi指占据晶格上i点的矢量，σi=0，…，Mi=1,如果Mi代表的元素种类只能占据点阵i。之后，我们需要一种方法来量化原子占据多重点阵如i1，i2之间的相关性。在计算中，有着特定相关性的点阵的集合被称为聚簇，用α表示。在二元系统中（每个点自由度都为一），定义哪个阵点在聚簇中至关重要（阵点在聚簇中时αi=1，不在时αi=0）。在多组元系统之中，聚簇中的每个阵点i都有一个αi值，范围从0到Mi-1，这表示了哪个组元的自由度被考虑了（αi=0时阵点不在集团之中）。关联函数与聚簇之间的关系是…α表示所有与α聚簇镜面对称的等价聚簇的平均值，Γ𝛼′(𝜎)是一个聚簇函数，定义为其中𝛾аi，Mi(𝜎𝑖)满足𝛾0，Mi(𝜎𝑖)=1，已经以下正交关系（二元合金之中，通常选择𝛾0，2(0)=1，𝛾0，2(1)=1，𝛾1，2(0)=−1，𝛾1，2(1)=+1，这样能直接关联到著名的伊辛模型。）多组元系统在在ATAT中通常的选择方法也已经叙述过了[12]。尽管在原理上来说所有的晶体阵点都可以算出结果，但是通常选择𝛾0，Mi(𝜎𝑖)=1，因为这样只输出在聚簇α上的阵点的结果。3.运算法则我们的运算法则能这样解释：𝜌𝛼（𝜎𝑟𝑛𝑑）表示完全无序状态的给定组元的关联函数，𝜌𝛼（σ）表示候选的SQS相σ的关联函数。Δ𝜌𝛼(σ)=𝜌𝛼(σ)−𝜌𝛼（𝜎𝑟𝑛𝑑）。因为在无序状态下，阵点的站位都是独立的，所以𝜌𝛼（𝜎𝑟𝑛𝑑）可以很容易的算出：𝛾𝑎′𝑖，Mi(𝜎𝑖)𝛼可以通过平均组元的阵点i计算得出。我们的目标函数（为获得最小值）为：其中，L表示l的最大值，A表示用户自定义的聚簇，表示用户自定义的重量。目标函数（3）中的-ωL是我们的方法中一个重要的特色。它反映了SQS的质量通常是用来衡量完全无序状态的关联函数的值，并能够匹配完全正确。特别的，当逐步扩大扩展的超晶匹配关联性的范围直到收敛时，可以发现在关联性范围较小时相关的性能可以更好的匹配。通常，SQS的建立是通过一对关联函数，虽然若要更加精确，还需要考虑多体的关联函数。方程（3）中的第二部分考虑了剩余关联函数的绝对误差，用于指导随机搜索的正确方向以扩展范围，达到完美匹配。通过以下步骤达到目标函数（3）的最小化：1.每个超胞都有指定的原子数量n并且进行计算。（尽管晶胞计算（n）是复杂的多项式，在计算上运用穷举法还是容易处理的）。2.对于每个超胞，根据用户自定义的每个亚点阵的组元，晶格位置的占据也是完全随机的。该算法需要获得精确的nsxts。此时原子的类型为t，亚点阵为s，ns表示超胞上占据亚点阵s阵点的原子个数，xts是用户自定义的组元成分。这些可以通过随机改变某一个ns来获得。通过随机排列s亚点阵中的ns可能的结构，可以置乱原子的站位。（如果nsxts不是整数的话代码终止，这表示给定的超胞尺寸与设计的组元成分无法匹配）。3.满足方程（3）的标准的各个超胞，就可以定义为最佳的SQS。（剩余的其他超胞可以留作他用）。4.目标函数（3）之后用来模拟退火循环，配置采样与exp（-Q/T）成正比，其中T表示用户自定义的虚拟温度，通过获得。在计算步骤中，算法会尝试随机改变超晶胞的形状和在超晶胞中交换两个原子在同一亚晶格的排列。Metropolis算法会判断接受或者拒绝这种改变步骤。5.新获得的超晶胞和结构会与到目前为止获得的最佳SQS进行比较,后者如果没有前者优化，就会被更新。每次获得更好的SQS，结果都会被输出到硬盘。上述算法可以很容易的并行化：只需一个简单的运行多个独立的实例代码（确保随机数生成器有着不同的种子），让它们分别分析空间结构的不同部分并保证过程中最佳SQS的不断生成。虽然不太可能，上述方法获得的SQS也有可能表示最小的晶胞。这是由于(看不懂。。。)。ATAT工具中的一个命令（cellcvrt-sbestsqs.out）允许用户检查是否存在最小晶胞。4.讨论几句关于使用SQS的提醒。一、通过SQS计算所获得的特定体系的物理性能的精确性并不仅仅取决于完全无序状态的关联函数与体系的关联函数是否相配。关联函数对不同的性能有着不同的影响，且这种关系本身就独立于系统，这在之前的研究中也有提到[1,10,15]。因此，我们不能事先得知（通过纯几何原理）为了获得精确值所取的关联函数应该是多少。所以，我们扩大范围，调查一系列的SQS，获得完美相配的关联函数收敛性以研究所需的性能。综上所述，基础理论预示着这样的收敛会快速发生在性能计算中，如总能量。正如此领域著名的“orderN”方法所述，许多物理参数都能导致局部总能量的扩张[16]。这种扩张通常表现为分簇层次结构[13]：小的聚簇（在空间拓展和阵点数量的角度）起着更主要的作用。然而，快速收敛也并不是与所有物理性能都有必然的联系。如带隙和电导率这些性能就可能受长波长组元的强烈影响，使得SQS的预测不准确。在预测张量的值，如弹性常数时，SQS也存在一个潜在问题。能获得优秀的标量性质的SQS可能获得一些不自然的各向异性张量性质[17]。这是由于耦合关联函数与性能的关系时，标量与张量运用了不同的等效对称规则[18]。幸运的是我们可以通过对称化技术来修正这一问题。5.常用代码5.1输入文件在下列描述中，方括号”[]”内的数据可以用适当值取代，没有在括号内的条目都是其字面的意思。Mcsqs代码需要2个输入文件。第一，一个定义随机状态的文件（缺省值为rndstr.in，会被-1=[filename]选项覆盖）。文件格式如下（与ATAT中的lat.in文件相类似，用来运行maps和corrdump代码，旦包含了局部原子站位信息）：1.建立坐标系，定义方向矢量，主要两种方法2.将每个点阵的矢量u⃑，v⃑，w⃑⃑⃑，在坐标系中表示出来，如下3.最后，给定每个晶粒在点阵中的站位xi（与点阵的矢量在同一坐标系上），同时给定站位于点阵i的可能的原子的元素类型t1i，t2i,…和相应的站位该原子的浓度C1i，C2i，…。需要注意每个点的浓度值之和应该为1。表1是这个文件的一个示例。所有对称等效的阵点必须站位情况相同。为满足此要求，两个亚晶格之间的“虚拟”相（不占据阵点）也必须计算在内。第二个输入文件（缺省值为clusters.out，会被-cf=[filename]选项覆盖）表示在匹配无序状态的关联函数时，哪些聚簇需要考虑在内。这个文件可以通过corrdump代码的功能来生成，命令行如下：Corrdump-1=rndstr.in–ro–noe–nop–clus-2=…-3=…-ro命令允许corrdump调用与mcsqs（命令中的rndstr.in）相同的晶格输出文件。-noe和-nop命令用来跳过mcsqs没有使用的空的或者特定的聚簇。-clus命令使代码只生成聚簇。-2=…-3=…命令表示选取的范围。若要获得更多信息，可以输入corrdump–h命令。下面是一些可选的输入选项：创建一个sqsparam.in文件（或使用-pf=[filename]命令）可以提供一个方法控制以下两个变量对过程的优化：方程（3）中的重量ω及退火循环过程中的模拟温度T。这个文件仅包含两个数据[ω][T](7)当代码运行时，这两个数据会根据表现实时调整。如果输入-rt=[rt]命令，这个文件会在每个rtMonteCarlo步骤时进行重读。创建一个stopsqs文件以随时停止代码。如果设置了-rt选项，sqscell.out文件就会选择特定的超胞进行输入。通常（不设置-rt），sqscell.out文件会通过穷举法来随机生成超胞。此文件先生成ns个超胞，然后选择其中满足方程（6）的超胞。-rt选项在两种情况下有用：一是当代码需要重算多次（或者平行计算）时，可以节省生成超胞的时间，二是当暴力算法花费的时间太长时，可以事先给定一些特殊的超胞。5.2命令列选项以下的命令可以从mcsqs[options]中调出，每个命令可能有多个选项。-n=[int]指定需要的原子或晶胞的个数（必须是rndstr.in中的阵点个数的整数倍）。这是唯一需要的命令行参数。-l=[string]输入文件指定随机结构（缺省值：rndstr.in）-cf=[string]输入文件指定聚簇（缺省值：clusters.out）-tol=[real]允许的SQS与理想随机状态的差值（缺省值：1e-3）。方程（3）中的|∆𝜌𝛼(σ)|项小于这个值时被看做0。-wr=[real]为获得关联函数与目标函数的完美匹配，质量（方程（3）中）的分配范围（缺省值：1）-T=[real]蒙特卡洛模拟退火循环时的温度（缺省值：1）-pf=[string]指定输入文件的最优化参数的名称（缺省值：sqsparam.in）。此文件（方程（7）中有介绍）一直会被读取（不论-pf选项是否被修改），如果存在，将会重置-wr和-T中的给定值。-rc从sqscell.out文件中读取超胞。如果未设定，代码会自动生成超胞并写入sqscell.out。-ip=[int]