在这篇教程中，你会学习到如何利用“刚性”（rigid body）结构限制和简单的优化步骤，来快速、可靠的弛豫一个器件结构的内坐标。

器件结构内坐标的弛豫，是一个冗长费时间的过程。主要有两个原因：

双电极体系由三个区域组成：左右两个电极，以及一个由电极扩展层和中心散射区组成的中心区域。在第一性原理计算中，每一部分的几何结构需要单独优化，但是同时中心区域的两电极扩展层必须和相应的电极匹配。

1. 优化电极。当优化完电极后，需要确定中心区域的电极扩展层是否匹配优化后的电极。
2. 优化中心区域。中心区域的最佳长度通常并不是提前已知的。对于周期性结构，需要使应力最小化。但是，应力场并不单单由器件结构确定，因此需要其他的方法确定。

在第二步中，我们必须在确保电极扩展层内坐标不变的前提下最小化器件总能量。在ATK里至少有两种方法实现。

首先是完全弛豫电极结构，并在此基础上构建双电极器件结构。之后，我们需要取出并弛豫中心区域的晶胞，但是必须对电极扩展层的原子进行“刚性”结构限制。之后可以恢复添加两端电极，重新得到器件（或是界面）结构。这个方法计算效率比较高，通常可以得到大部分的弛豫结构。我们称这种方法为“Bulk Rigid Relaxation (BRR)”。 BRR方法不计算整个器件的总能量，因此重新得到的器件结构可能不是能量最低的几何构型。部分电子态性质可能会对这种差异敏感。

更有力的方法则必须使器件处于能量最小化状态（相对于内坐标和中心区域长度），即在改变中心区域长度的同时进行整个双电极器件的计算。这是一个一维原子链最小化的问题，在每一步都针对原子位置作弛豫。因此，我们称之为1DMIN。为此，我们发展了用ATK进行这种最小化的方案，使之尽可能有效地处理这个问题，但是这种方法的计算量远大于BRR方法。因此，我们建议在对中心区域使用BRR弛豫之后，再进行1DMIN计算。

这篇教程介绍了利用BRR和1DMIN方法对电极和中心区域的弛豫。教程中的实例是用经典力场方法对银-金接触表面进行弛豫。利用经典势可以使计算更快捷，同样的方法可以适用于DFT计算。通常，我们会建议利用DFT计算来检验ATK-Classical计算的结果。

生成一个“device_relaxation”的新项目。我们从一个未弛豫的双电极结构开始。其结构可以参照此处教程 Building an interface between Ag(100) and Au(111)，或者从此处下载结构脚本如何弛豫器件体系的几何结构。 在项目中保存py文件，并将其拖拽至Builder中。 两电极结构如下图所示。它由左电极Ag(100)、中心区域、右电极Au(111)构成。利用插件，将器件分割成左电极、右电极和中心区域。

电极的几何结构是由Ag(100)和Au(111)构成（参考Building an interface between Ag(100) and Au(111)）。其左右晶胞的结构采用了实验上的几何数据。我们假设实验结构也是理论上能量最低的结构。（实际理论研究中需要弛豫晶胞以验证这一点）

在接触器件构建中，我们选择Strain second surface（对第二个表面施加应力）处理表面晶胞。这意味着在形成界面的过程中，Ag(100)结构将会保持实验上的结构，而Au(111)表面将在x、y方向上将被拉伸。当晶胞在x、y方向上拉伸之后，它将按照材料的泊松比在z方向上弛豫。

接下来，我们将把这种弛豫应用于Au(111)晶胞，方法是在保持x、y方向固定的情况下对z方向上的应力进行弛豫。获得的Au(111)结构将会应用到中心区域的金原子中。

打开Builder并重命名右电极结构为Au111，双击激活结构，并发送至Scripter，添加：

New Calculator

OptimizeGeometry

打开New Calculator，选择ATK-Classical计算方法，并选择“EAM_Zhou_2004”势。点击OK。 下一步，点开 OptimizeGeometry ：

1. 设置force tolerance为0.01
2. 设置stress tolerance为0.001
3. 取消z方向的晶胞固定（Constrain cell）

发送脚本至 Job Manager 并进行计算，计算耗时不到 1分钟，当执行完之后，检查log文件。 我们将会在 LabFloor看到一个“Au111.nc”的nc文件，它包含了弛豫前后的两个结构的结构数据。选中两个结构，并拖拽至Viewer。在Viewer中，鼠标选中晶胞，我们会发现，弛豫前的结构，C方向矢量为7.064 Å ，弛豫后为7.011 Å。这个结果同样可以在“Au111.log”文件中获得。 弛豫后的结果压缩了0.75%。为了应用到中心区域，打开Builder，双击打开中心区域结构。

1. 在stash面板重命名为“central”
2. 鼠标选中所有金原子
3. 打开右侧Bulk Tools ‣ Stretch Cell插件
4. 在C轴中输入0.9925
5. 使用“Selected atoms are” stretched.
6. 点击Apply。

现在，我们仅需要进行中心器件弛豫。当进行几何优化时，我们希望电极扩展层可以自由移动，因此我们需要在两端添加真空层： 鼠标点击绘制窗口的背景区域，确保没有原子被选中。点开右侧的Bulk Tools ‣ Stretch Cell ，利用C-vector slider，延长晶胞30%的长度，或者简单将因子设置为1.3。然后点击Apply应用，晶胞就会相应的增大。 点开Coordinate Tools ‣ Center，取消勾选A、B，点击Apply，可将中间区域制动设置到中心。

我们要固定电极扩展层的原子，但为了在后面便于使用“Device from bulk”工具添加电极时确定电极的长度，我们需要比常规的电极扩展层多固定一层。为此，选中最外层的5层银原子，点开Selection Tools ‣ Tags。在“Click here to write a new tag…”输入框中，输入Left。

接下来，同样选择最外层的4个金原子， 标记为Right.

发送结构至Script Generator，设置如下参数：

1. 添加New Calculator和OptimizeGeometry项目。
2. 选择ATK-Classical方法，和EAM_Zhou_2004势做计算。
3. 在几何优化中设置如下参数：
4. 设置force tolerance为0.01 eV/Å.
5. 点击Save trajectory，设置文件名为central.nc。
6. 单击Add Constraints，对两电极扩展层的设置“Rigid”结构限制.
7. 检查Left和Right两个标签要标记的原子是否正确。在两个标签的下拉菜单中选择“Rigid”。
8. 点击OK，并返回Scripter。

发送脚本至Job Manager，并运行计算。计算非常快，一旦计算结束，我们将会看到central.nc文件出现在LabFloor上。

可以从Viewer中看到优化的轨迹动画。

优化后的中心区域结构保存在central.nc文件中的gID0002中。这是我们对中心区域结构的最初猜测。将其拖拽至Builder中。

打开Device Tools ‣ Device From Bulk，将bulk结构转换为器件。保留默认参数即可。

结构如下图所示，右键或F2重命名为“device”。

我们现在可以搜索器件全局最小能量结构了。通过Save as,将器件结构保存为NetCDF文件，并命名为“device_in.nc”。

现在我们将用SciPy minimization程序搜索可以使能量最小化的中心区域长度。我们将从device_in.nc的器件结构开始。 我们需要两个脚本： device_relaxation/device.py和 device_relaxation/optimize.py。将两个脚本和device_in.nc放在同一个目录下。第一个脚本是读取结构并设置计算参数，第二个脚本包含了计算所需的类和函数的定义。

两个脚本都可以在VNL的Project Files中看到。拖拽device.py至Jobs，并运行脚本。

一旦计算完成，将会出现一个数据图窗口。蓝色点是弛豫结构的总能量。晶胞的最初长度是28.51 Å，其他的蓝点是算法模拟出来的结果。其中位于28.58 Å的红色点为中心区域能量最低的点。优化程序只简单地绘制出所示的势能曲线，并找到它的最小值。

手动关闭弹出的窗口可以结束脚本运行。

弛豫器件结构保存在输出文件“device_out.nc”，并包含总能量。红色点处的结构（最终找到的能连个最低点）是最后一个器件结构（id最大的那一个）。

在Viewer中，能量最小点的结构如图所示：

这个例子里，1DMIN弛豫后结构与BRR弛豫的结构相比较，变化十分细微。这可以看出，通常情况BRR方法的结果已经足够精确，并不需要再进行1DMIN的步骤。1DMIN仅用于极其精确的计算。

下载pdf, Infographics。点击顶部按钮查看图片。

• 英文教程：英文。

本文翻译：闫强