版本:2018
QuantumATK 为电子传输计算提供了一种市场领先且可实现高度优化非平衡格林函数(NEGF)的方法。即便如此,实现 NEGF 计算快速和可靠收敛从而获得电子基态有时也可能具有挑战性。因此,本指南为收敛性较差的 NEGF 计算进行故障排除提供了行动要点。
我的计算不收敛是什么意思?
就像周期性块体的标准计算一样,NEGF 器件计算采用自洽场(SCF)方法迭代搜索电子基态。在 SCF 循环的每一步都解决了电子结构问题,并计算了总能量。当 SCF 迭代之间的能量变化降到某个公差以下,也就是找到自洽电子密度时,该计算达到收敛。
如果计算不收敛,则意味着在允许的 SCF 迭代次数内未找到电子结构问题(基态)的自洽解。在进一步分析之前,您的任务是修改 NEGF 计算,使其收敛。
如果计算未收敛,QuantumATK 将发出警告。只要计算没有达到合适的收敛时,QuantumATK 的计算器日志文件中将包含以下警告消息:
################################################################################ # # # Warning: The calculation did not converge to the requested tolerance! # # # ################################################################################
查看 QuantumATK 日志中的此类警告是一种基本习惯,因为如果计算未收敛,结果可能就没有用。取决于作业的执行方式,日志流将被发送到系统 stdout
和/或一个文件。从 O-2018.06 版本开始,如果计算未收敛,日志文件还将包含具有详细信息的 Convergence Report。
如上所述,器件计算收敛到基态很重要。如果在零偏压计算中没有发生这种情况,则通常是由于器件构型设置不正确或 NEGF 器件计算器中的设置不够理想。有限偏压计算有时可能会带来其他的收敛性问题,这也是我们在此要解决的问题。
以下列表为零偏压和有限偏压计算提供了故障排除指南。
收敛困难的一个非常普遍原因是器件的中央区域太短。具体来说,中心区域对于在其中的载流子而言太短,以至于不能从界面作用中屏蔽电极。换句话说,NEGF 形式主义要求静电势达到平衡,且在到达电极之前变平。如果不能满足要求,NEGF 算法将无法将其与块状电极匹配。这种情况下补救方法很简单,即延长中心区域。
有关屏蔽重要性的更多信息,请参见 NEGF: Device Calculators 的描述。
在下图中,我们显示了左侧含金属器件和右侧硅的 Hartree 电势差。垂直线表示两种材料之间的界面。
图190:沿着硅屏蔽区域长度的 Hartree 电势差图。请注意,80 Å 不足以使电势变平。
从上图我们可以看到,硅距离界面 80 Å 处电势并未变平坦,而是在 100 Å 处。金属那一侧,我们观察到仅在几埃内就足以屏蔽界面并按预期使电势变平坦。
通过增加掺杂或者如果已经掺杂的情况下通过增加掺杂水平,可以减小耗尽区的长度,即屏蔽发生处半导体的长度。在下图中,我们展示了硅耗尽长度与掺杂水平的关系。
图191:屏蔽(耗尽层)长度与掺杂浓度的关系图。例如,1019 cm-3 的掺杂对应于约 100 Å 的耗尽层长度。
收敛问题的另一个常见原因是电极的尺寸。电极必须足够长才能确保第二近邻晶胞之间在传输方向上没有直接相互作用,这一点非常重要。换一种说法就是电极扩展和电极的第一周期图像之间必须没有相互作用。一个好的经验法则是,电极在传输方向上的长度应是基组中最长基函数的两倍以上。基函数的长度可以在 New Calculator 小部件中绘制函数图查看。
更彻底的测试是使用 NanoLab 中的 Electrode Validator 工具。使用它之前首先要对电极进行自洽计算,然后在 LabFloor 上选择结果生成的 BulkConfiguration,然后点击右侧插件栏中的 Electrode Validator 工具。然后,NanoLab 将通知您电极是否有效,即是否足够长以满足上述标准。
列表中的下一项是确保计算使用具有正确边界条件的泊松求解器。如果要以这种方式对器件建模,则边界条件在 A 和 B 方向应具有周期性,如果在该方向上有栅,则边界条件应为 Neumann。最后,在 C 方向上应为 Dirichlet 边界条件。
接下来要检查的是在运输方向上的 k 点取样。为了匹配电极和中心区域之间的费米能级,电极 k 点网格沿 C 方向(传输方向)排布密集是很重要的。沿 C 方向的默认密度为 150(相比之下,否则为 2 至 7),有时甚至还有必要增大。注意,在大多数情况下,电极计算在总 NEGF 计算时间中是一个微不足道的部分,因此通常不必担心这部分的操作。
另外一个重要的点是检查两个横向上 k 点的数量合理性地收敛。如果它们在一个方向上太低会导致非物理性的结果,也会妨碍收敛。下图我们展示了一个计算的费米能级与横向 k 点数量的关系。在这种情况下,我们可以预料到 1、2 和 3 个 k 点会给出非常差的结果,并且可能会出现收敛问题,而从 5 个点开始,结果会合理性地收敛。
图192 费米能级与 k 点数量的关系图。在这里,我们显示的是沿非运输的 A 向和 B 向。
下一个选项是尝试将初始密度的计算更改为 EquivalentBulk
,替代默认的 NeutralAtom
。默认行为是创建一个将中性原子密度直接加和作为初始密度的猜测。EquivalentBulk
则是进行自洽计算,把中心区域视为块体构型,将所得密度作为 NEGF 计算的初始猜测。请注意这个计算可能要求很高,可能并不如 NEGF 计算本身那样可以有效地在许多计算核心上并行。
NeutralAtom
法在大多数情况下都适用,但对于某些系统就很有必要使用 EquivalentBulk
以从更好的初始猜测开始进行器件计算。但如果您的电极不相同和/或采用 Hückel 模型,则使用 EquivalentBulk
很可能不占优势。
特别是在金属或小间隙半导体中,如果没有充分地描述接近费米能级的状态,则会阻碍收敛。如果是这种情况,这些状态将在被占用和未被占用之间切换,从而导致收敛出现问题。这可以通过升高电子占有函数的温度(宽度)和/或尝试使用非默认占有函数之一来缓解。请注意,主页关于 Occupation Methods 的介绍中提到较高的温度可能会导致结果略有误差。
在 NanoLab 中设置计算器时,轮廓积分的下限会自动设置得很低,以使其包含构型中元素的所有价态。这个事实非常重要,接下来要做的是检查是否确实如此。如果您对此有疑问,可以对您的材料进行 DOS 计算,以确定最深的价态及其能量,并验证轮廓积分的下限在能量上低于这些状态。
有时,可以通过增加默认 SemiCircleContour
中半圆部分点的数量或切换到 OzakiContour
来提升收敛性,在 OzakiContour 中可以更容易地增加极点数(等效于轮廓点),因为只有一个参数可以调整。注意,每个轮廓点/极点都必须像 k 点一样单独计算,因此计算成本随轮廓点/极点数变化近似为线性。有关更多信息,请参见 NEGF 的手册页 NEGF: Device Calculators 和两种轮廓积分类型:SemiCircleContour
和 OzakiContour
。
图193 使用半圆轮廓时的轮廓积分的草图。有关更多详细信息,请参见手册页。
更改 Poisson 求解器也可以帮助收敛,例如从 MultigridSolver
到 ParallelConjugateGradientSolver
,因为它们使用了不同的数值方法。这意味着计算得出的静电中的数值噪声将有所不同,这反过来又会影响具有更复杂静电系统(例如栅)的收敛。还要注意,对于并行运行的大型计算,ParallelConjugateGradientSolver
通常会更快一些。
零偏压计算的最后一项检查是密度网格截断。与横向的 k 点一样,它应该足够高以提供正确的物理结果。确保就感兴趣的量进行适当收敛也将有助于收敛。下图我们展示了费米能级随网格截段的变化而变化,可以看出非常低的网格截断在确定费米能级时会带来不准确性。然而还要注意,此图和上文与 k 点有关的那个图在比例上的差异。
图194 计算的费米能级与网格截断的关系。
通常情况下,有限偏压的计算要比零偏压更难收敛,以上所有要点也适用于这种情况。但还有另一件事可以尝试,将在下一章节中进行说明,并且如果零偏压计算能够很好地收敛时,则应该首先尝试一下。
有限偏压计算应始终从较小偏压或零偏压的收敛计算中重新开始。如果不收敛,它通常将有助于减小偏压差异。例如您想要从零偏压开始进行 0.4 V 的计算,最好的方法是逐步计算,比如以 0.1 V 为步长。以下示例脚本可以完全做到这一点:↓ iv.py。
# Set bias # Positive: Forward # Negative: Reverse bias_list = [0.10, 0.20, 0.30, 0.40]*Volt # Read DeviceConfiguration zero_bias_file = 'device_zero_bias.hdf5' device_configuration = nlread(zero_bias_file, DeviceConfiguration)[-1] for bias in bias_list: if processIsMaster(): print("Bias is now: ", bias ) # Get the calculator calculator = device_configuration.calculator() # Set the bias voltage calculator=calculator(electrode_voltages=(bias/2, -bias/2)) # Attach the calculator and use the old initial state device_configuration.setCalculator( calculator(), initial_state=device_configuration) device_configuration.update() nlsave('device_bias_%.2f.hdf5' % bias.inUnitsOf(Volt), device_configuration)
这种方法对于较高偏压会变得更重要,因为可能需要更短的增量才能实现收敛。因此在上述情况下,也许 0.1 V 和 0.2 V 时收敛,但 0.3 V 时不会。建议尝试 0.25 V,从 0.2 V 开始以 0.05 V 递进。当然,这种方法是否可行取决于所研究的系统以及您希望实现的最终目标。
如果您已经验证了 NEGF 计算中使用的物理模型是合理且适当的(请参阅上文),但计算仍然无法收敛,那么最后一个选择就是调整控制 SCF 循环的参数。
SCF 参数在 ATK Python 的 IterationControlParameters
类中定义——导航至 New Calculator 程序后在 Iteration control parameters 条目调整设置。
即使是稳健且表现良好的 NEGF 计算,有时朝着基态收敛的速度也会比预期的慢,因此可能需要额外的 SCF 步骤才能正确收敛。
在这种情况下,只需增加参数 Maximum steps。
然而,首先要从失败的计算中查阅 QuantumATK 的日志文件,确认 SCF 循环实际上是否运行良好,即是否收敛到了合理的基态。如果没有,那么额外的 SCF 步骤就不可能解决问题——可以考虑调整 Pulay 混合。
首先,能带能量和哈密顿矩阵元素都应稳定地收敛于恒定值。检查 ATK 日志文件:在每次 SCF 迭代之后查找报告能带能量 E
和与上一次 SCF 迭代的能量差(dE
)和汉密尔顿矩阵元素(dH
) 的所在行。在终端,可以很方便地使用 grep 命令实现:
user@machine $ grep dE atk.log | 0 E = -70.6148 dE = 5.719151e-01 dH = 7.266251e-01 | | 1 E = -56.8671 dE = 1.374766e+01 dH = 5.675282e-01 | | 2 E = -58.5695 dE = 1.702319e+00 dH = 1.737177e-01 | | 3 E = -59.6745 dE = 1.105067e+00 dH = 9.875918e-02 | | 4 E = -59.3529 dE = 3.216787e-01 dH = 4.048304e-02 | | 5 E = -59.0695 dE = 2.833359e-01 dH = 2.447425e-02 | | 6 E = -59.2082 dE = 1.386839e-01 dH = 2.978777e-03 | | 7 E = -59.0397 dE = 1.684759e-01 dH = 1.492602e-02 | | 8 E = -59.078 dE = 3.826926e-02 dH = 1.691435e-03 | | 9 E = -59.05 dE = 2.800614e-02 dH = 2.041091e-03 | | 10 E = -59.0559 dE = 5.931576e-03 dH = 1.624175e-03 | | 11 E = -59.0553 dE = 6.427625e-04 dH = 2.825693e-04 | | 12 E = -59.0543 dE = 1.008929e-03 dH = 1.880459e-04 | | 13 E = -59.0542 dE = 1.964723e-05 dH = 6.034241e-05 |
dE
和 dH
应当朝着指定的公差(默认为10-4)持续稳定地减小。
其次 Mulliken populations(马利肯布居数)应在物理上合理。检查 ATK 日志文件中最近几次 SCF 迭代的 Density Matrix Report 。针对每个原子,报告了 DM
和 DD
的值。前者是器件密度矩阵(马利肯布居数)的迹,而后者是 DM
与原子上的中性原子电荷之间的差(由赝势给出)。
在下面的示例中,我们看到如何满足这些要求:
+------------------------------------------------------------------------------+ | Density Matrix Report DM DD | +------------------------------------------------------------------------------+ | 0 Si [ 1.907 , 1.100 , 1.564 ] 3.89764 -0.10236 | | 1 Ni [ 0.001 , -0.000 , 2.342 ] 18.21499 0.21499 | | 2 Si [ 0.001 , 2.201 , 3.120 ] 3.89379 -0.10621 | | 3 Si [ 0.001 , -0.000 , 4.675 ] 3.89556 -0.10444 | | 4 Ni [ 0.001 , 2.201 , 5.453 ] 18.21180 0.21180 | | 5 Si [ 1.907 , 1.100 , 6.231 ] 3.89446 -0.10554 | | 6 Si [ 0.001 , 2.201 , 7.787 ] 3.89564 -0.10436 | | 7 Ni [ 1.907 , 1.100 , 8.565 ] 18.21115 0.21115 |
DD
值都为负,而 Ni 原子的相应正值是是它的两倍。
SCF 循环使用 PulayMixer
在每个 SCF 步骤中通过将之前一些步数中的密度混合到当前步骤的猜测中来生成电子密度。Damping factor 参数控制在混合中使用前几个密度的百分数,History steps 控制要使用前面步骤中密度的数量。
通常默认的 Pulay mixing 参数应该可以实现快速、稳定的 SCF 收敛,但有时可以进行一些调整:
SCF 循环将会一直持续至达到精度公差或最大步数为止。该公差是对总能带能量和哈密顿矩阵元素都进行测量的——在两种情况下,一个 SCF 步骤与另一个 SCF 步骤之间的绝对差都必须低于指定的公差,以使计算被认为是准确和收敛的。
如果任何尝试都失败了,请在在线 QuantumATK 论坛(https://forum.quantumatk.com)上发帖或联系 QuantumATK 支持(quantumatk-support@synopsys.com)。在使用这两种方式联系时,请附上您的 QuantumATK 脚本和日志文件。