目录
动力学矩阵和声子的计算
软件版本：QuantumATK-2023.12
动力学矩阵（DynamicalMatrix）是一个用于计算许多与晶格动力学相关特性的对象，例如迁移率（Mobility）、拉曼光谱（RamanSpectrum）、形变势（DeformationPotential）、介电张量（DielectricTensor）等。在本教程中，您将学习如何计算动力学矩阵，并使用它来执行声子能带结构（PhononBandstructure）计算。特别是，您将学习如何设置计算，以及哪些参数至关重要。您还将学习如何通过最佳并行化策略优化此类计算的性能。
- 您将首先使用LDA交换关联函数执行DFT-LCAO计算。您将了解声子计算中使用的所有相关参数。
- 最后，您可使用力场计算（ForceFieldCalculator）模块中的经典势计算，并与前面计算结果进行比较。
预备知识
本教程介绍QuantumATK中的动力学矩阵功能，只需要基本熟悉QuantumATK图形界面即可。这可以通过以下基本教程和指南实现：
- 如何使用QuantumATK进行快速简便的原子尺度建模
- Workflow简介
- 介绍
- 在Nanolab数据视图中组织数据
- 几何优化：CO/Pd（100）
- 本地执行QuantumATK脚本的作业管理（Job Manager）
创建工作流
转到Builder并从内置数据库导入Silicon（alpha）结构。然后将结构发送到工作流生成器，并添加以下脚本块：
- 原子轨道线性组合方法计算（LCAOCalculator）
- 几何优化（OptimizeGeometry）
- 动力学矩阵（DynamicalMatrix）
- 声子能带结构（PhononBandstructure）
生成的工作流应如下图所示：
现在我们将浏览每个模块的设置并编辑重要参数。
LCAOCalculator设置
为了获得动力学矩阵的良好结果，必须设置好几个关于收敛的关键参数：
注意
- 对于一般的计算模拟，我们建议使用默认的 PseudoDojo-Medium 基组，或提高基组精度到 High。
- 还请注意，所需的 Density mesh cutoff 是特定于材料和赝势的，因此，虽然对于硅（Silicon）来说，50 Hartree足够了，但对于其他材料来说，太低了。一般情况下，应该从默认值开始，然后增加，直到达到收敛。
在Workflow Builder中执行以下操作：
- 打开LCAOCalculator模块。
- 保持Pseudopotential为Pseudo Dojo不变，在 Basis Set 下选择High。
- 在 Main 部分的 Numerical Accuracy 中，将Broadening设置为300 Kelvin，将 Density mesh cutoff 增加到 50 Hartree，并将k点密度设置为4 Angstrom。
- 在Iteration Control区域中，将Damping factor设置为0.5。同时将Tolerance降低到0.000001。
提示
- 通常需要从DynamicalMatrix计算的默认值增加网格截断能和k点密度，以确保声子能带结构的收敛。
- 对于一般的非金属，我们建议使用300开尔文扩展。
- 对于DynamicalMatrix计算，通常建议使用比默认参数更严格的公差。
- 通过调整阻尼因子和使用的历史步数，可以减少迭代步数。
晶格优化设置
打开OptimizeGeometry模块，并如下图所示修改参数，以执行比默认收敛参数更严格的力和应力最小化。
提示
- 我们通常建议在计算动力学矩阵时减少目标应力和力。
- 为了获得计算电子性质的良好结构，默认设置通常就足够了。
动力学矩阵设置
动力学矩阵是能量的二阶导数，它对应于力的一阶导数。力的一阶导数是通过有限差分方案计算得出的，其中心晶胞中每一个对称唯一原子都会沿着所有笛卡尔方向发生微小位移，这一过程也被称为冻结声子计算。
对于具有小体积单晶胞的晶体，通常需要周期性重复单晶胞的超胞，以准确考虑动力学矩阵中的长程相互作用。QuantumATK 可以自动处理这一需求。默认情况下，单晶胞会被重复扩展，以便将单晶胞中原子周围、依据预设且依元素而定的距离内的所有原子都包含进重复的超胞中。这对于大多数体系而言已经足够。同时，也允许直接设置在每个方向上重复的次数。更多详细信息，您可以在DynamicalMatrix手册条目的“注释”部分查看。
- 在动力学矩阵对象中，我们将采用默认设置。
- 打开PhononBandstructure模块，并选择每段400个点，以及沿Brillouin区路径 G, X, K, G, L 进行计算。
注意：
计算的总模拟时间（TAT）高度依赖于重复系统的次数，原因很简单，即重复后的系统可能会变得非常庞大。对于结构复杂或较大的系统，如果自动检测到的重复次数导致总模拟时间非常高，有时可以减少重复次数。然而，对于实际研究计算，重要的是，要始终检查与诸如声子能带结构等因素的收敛性，以确保结果的准确性。
为了测试使用不同重复次数的效果，我们再次进行计算，将重复次数设置为3x3x3和5x5x5。完整的计算流程可在此处下载：Silicon_phonon_bandstructure.hdf5。若要更改重复次数，需将Repetitions设置为Custom，然后选择所需的数值。
运行计算
您现在可以运行DFT声子能带结构计算了。这是通过使用Send To图标将脚本发送到Job Manager来完成的。
但在开始作业之前，有几项细节需要考虑：
并行化
在QuantumATK中，动力学矩阵的计算可以通过针对位移进行并行化来实现。所需位移的数量取决于体系，可以通过打开Dynamical Matrix对象并按下“Calculate……”按钮，在Number of Displacements下计算得出。
每个位移所分配的进程数量可以通过Processes per displacement参数进行控制。一般来说，我们推荐以下两种并行策略之一：
- 默认设置ProcessesPerNode会将位移任务分配到各个节点上。当在单个节点上运行时，这意味着所有位移都将按顺序计算。
- 根据考虑的体系和硬件配置的不同，有时将所有位移尽可能并行地在多进程中运行更为有利。对于这种策略，应取总进程数减一（其中一个将作为调度进程预留），然后将剩余进程均匀分配到位移数量上。
本教程中，我们将使用单个20核节点上的默认ProcessesPerNode。
体系规模
本教程中运行的三个体系，分别重复3x3x3、5x5x5和7x7x7单元胞，分别包含54、250和686个硅原子。除了运行时间外，你还应考虑所需的内存量。
计时
下表报告了在上述参数设置下，运行完整声子能带结构工作流程所需的总时间。这些计算使用20个MPI进程并行运行： Repetition 3x3x3 5x5x5 7x7x7 Total time 1m29s 22m26s 1h39m44s
分析结果
当任务完成后，您将在Data中找到PhononBandStructure。请将其高亮显示并点击Open以查看。
下图展示了之前定义的3x3x3、5x5x5和7x7x7系统的分析结果。
您会立即注意到，3x3x3重复单元的体系太小，无法给出哪怕是定性良好的结果，而中间尺寸的体系已经是一个相当不错的近似。然而，您也可以发现，在5x5x5和7x7x7体系接近X点的声学模式之间存在微小的差异。
现在，动力学矩阵与BulkConfiguration一起存储。这意味着您可以读取动力学矩阵，并在短短几秒钟内执行另一个PhononBandstructure计算（例如，在不同的布里渊区路径上）或PhononDensityOfStates计算。
使用力场加速计算
QuantumATK提供了ForceFieldCalculator中的大量经典势函数。这些势函数非常适合用于研究材料的振动特性，比如本教程中计算的硅的声子能带结构。
注意
经典势仅对有限数量的系统可用。随着您的系统变得越来越复杂，比如包含界面、缺陷或多元素时，您经常会发现没有相关的经典势存在。在这种情况下，QuantumATK提供了使用力矩张量势来训练您自己的机器学习势的功能。如果您想尝试这项功能但尚未获得访问权限，请联系我们。
将硅块体发送到WorkflowBuilder。添加一个ForceFieldCalculator。现在可以按照上述说明设置应力优化和声子能带结构计算。
在ForceFieldCalculator中，选择适用于硅的最恰当势，例如Tersoff势或Stillinger-Weber势。
每个势都附带有文献参考，您应检查这些参考以确定势是否适合您的体系及您感兴趣的性质。
完全按照LCAOCalculator案例设置OptimizeGeometry模块和PhononBandstructure分析对象，然后运行计算。此计算在普通台式机上只需几秒钟。
使用Tersoff_si_2005势[1]计算的硅的声子能带结构如下面的图所示。
参考文献
[1] Erhart and K. Albe, “Analytical potential for atomistic simulations of silicon, carbon, and silicon carbide”, Phys. Rev. B, 71, p. 035211, (2005).