1. 处理收敛问题
在电荷自洽计算和结构弛豫过程中，有时pw.x计算会算到最大迭代次数或报错退出，这是计算不收敛的问题。广义地说，任何计算所关心的物理量，如果随着截断能、k点密度等收敛量的变化而有明显变化，也属于不收敛的问题。
笔者还不知道任何数学方法可以保证对Kohn-Sham方程求解的收敛，实际上，不收敛的情况确实存在，好消息是，收敛的情况并不少见——任何值得发表的、可以作出确定性结论的结果都应该建立在计算收敛的基础之上。本文试图全面概括QE计算不收敛的原则和解决方案，难免挂一漏万，烦请告知。如果有体系确实不复杂而又难以收敛的算例，也欢迎提供。
收敛问题一方面是物理的原因，也就是计算不符合实际，包括：
（1）建模型不符合实际：输入的原子坐标、元素种类出错；对称性的计算设置。
（2）使用的近似不符合实际：性能良好的赝势是实现计算的必要条件；对于金属的k空间取样密度和展宽方法；具有磁性，强的自旋轨道耦合；局域d、f轨道的过渡金属元素on-site库仑排斥（Hubbard U或Hybrid functional）；分子晶体等具有的范德瓦尔斯力；表面slab模型的偶极电场等。这些不收敛的情况也加深了我们对体系的认识，可能需要相应的方法才能收敛以及正确地刻画物理体系。
另一方面是数值计算的原因，包括矩阵对角化算法、自洽迭代算法、结构优化算法不收敛，即算法出错或迭代次数过多。结构不符合实际是一部分数值计算出错的根源：实际中不存在的结构，在计算时会产生与原子初始电荷偏离较多的电荷分布，增加了收敛的难度。当然，不收敛的原因也包括数值计算算法本身具有的不稳定性，计算精度和计算机性能的矛盾，计算机内存、硬盘空间、硬盘文件数和CPU资源不足等。软件版本的升级迭代可能会影响收敛的情况。
电荷收敛的有利条件是初始电荷（自旋）密度分布接近于电子-离子体系能量极小值状态的电荷（自旋）密度分布。实现收敛的原则是建立的模型、使用的近似符合实际，同时适当减少冗余。如果计算中出现电荷不收敛，首先要检查的是结构是否合理（例如，两个位置重合的原子，晶格常数的单位，甚至错误的元素种类）。
有时在scf或nscf计算会出现类似以下情况：
iteration # 1 ecut= 70.00 Ry beta=0.80
Davidson diagonalization with overlap
c_bands: 1 eigenvalues not converged
c_bands: 3 eigenvalues not converged
c_bands: 1 eigenvalues not converged
这是一个warning，是指对角化迭代过程中，最高的1-4个本征值计算没有收敛（阈值由diago_thr_init控制），超过5个本征值不收敛（或nbnd/4）则会报错退出。对总能和能级有一定影响，如果不是最后一步，通常可以忽略，用本页所述方法可能会消除这个warning。
1.1 实现自洽计算收敛的原则与方法
自洽计算收敛是在电荷自洽计算过程中，总能不再明显变化。判据是pw.x输出中的estimated scf accuracy小于给定值conv_thr，默认值(1.D-6)相对宽松，根据计算的物理量有时还需进一步提高收敛标准。
实现收敛的方法：
方法一：
修改电荷混合参数。(1)降低 mixing_beta。默认值是0.7，将mixing_beta调低至0.3 ~ 0.1 ( even smaller 0.05 for very long slab)。
mixing_beta的值越低，混合新的电荷就越少，混合后的电荷和原电荷越接近，有利于保证收敛。
(2)增大mixing_ndim，也就是上溯几个步骤的电荷会考虑到mixing中, 默认是8，同时计算所需内存会增大。
mixing_ndim = 12
方法二：
提高ecutwfc和ecutrho，以提高平面波个数和电荷及势能计算网格的密度（nr1, nr2, nr3）。对于超软赝势（US PP）ecutrho默认是ecutwfc的4倍（模守恒赝势4倍即可），增加到10倍
ecutrho=10*ecutwfc
方法三：
对于金属体系（或窄带隙，半金属等）：(1)加一些空带: 增大nbnd，具体取值根据体系的电子总数，默认最少4个空带，增加到足够空带，让最上面的空带占据几率趋向于零；(2)增加k点网格密度，以消除半满带对总能收敛的影响；(3)同时逐步地增大展宽，直到总能收敛。
occupations = 'smearing',
smearing='marzari-vanderbilt',
degauss=0.01
对于绝缘体、半导体，
occupations = 'fixed',
或使用极小的degauss
occupations = 'smearing',
smearing='gauss',
degauss=1.0d-9,
方法四：
对于对角化出错的体系，在检查结构是否合理（例如，两个位置重合的原子）的基础上，尝试使用diagonalization='cg'，默认是’david’，同时设置startingwfc='random'。减小初始迭代时的对角化阈值diago_thr_init=1.0d-6（scf默认1.0d-2，这一选项也是nscf计算阈值）。
方法五：
真的需要那么大的真空尺寸吗，在平板模型中（slab），真空和平板越厚，收敛遇到的问题就越多。检查过表面能的收敛了吗（而不只是总能的收敛）？请从薄的平板做起，逐渐增厚。
1.2 实现结构驰豫收敛的原则与方法
结构驰豫收敛的判据是相邻两个离子步总能变化小于etot_conv_thr(1.0D-4)，以及各个离子受力小于forc_conv_thr(1.0D-3)，对于vc-relax还要求cell受压强小于press_conv_thr (0.5)。默认只进行nstep个离子步（默认50），在使用bfgs优化结构时，离子步scf收敛阈值conv_thr还会逐渐减小直到conv_thr/upscale，upscale默认为100.0。这些收敛默认值相对宽松，有必要时还需进一步提高收敛标准。建议使用默认的ion_dynamics='bfgs', cell_dynamics='bfgs'，而'damp'系列方法收敛很慢，只适用初始很接近优化结构的情况。
如果relax计算的前几步电荷正常收敛，而进行到某一步报错，即结构优化不收敛，参考以下方法。
方法一：
对于原胞做vc-relax得到理论晶格常数。用理论晶格常数建超胞，而对于超胞只优化ion，即calculation='relax'，不做vc-relax。
方法二：
如果relax进行了若干步后停止，这种情况有时是因为体系优化后带隙消失，查看带隙可以用这里的脚本(gappw.sh relax.out)。这时重新按照金属进行自洽计算
occupations = 'smearing',
smearing='gauss',
degauss=1.0d-2,
方法三：
如果relax达到最大迭代步数nstep后停止，这对于原子个数较多（约64个）的体系也属于正常情况，relax收敛标准不是过高的时候，可以用最后一步的结构作为初始结构继续优化，并提高nstep（默认50）。
方法四：
磨刀不误砍柴工，用较高的精度计算力。适当增大截断能ecutwfc, ecutrho，减小conv_thr（1.0d-8至1.0d-9），即提高scf的收敛精度，以计算出较准确的力。vc-relax计算时降低press_conv_thr至0.1。
方法五：
计算初始化时会自动寻找体系对称性，并通过找到的对称性减轻计算量，如果模型中的对称性是符合实际的，应该保留，但是也存在部分情况，初始对称性过高，体系对称性自发破缺（如铁电、简并基态的Jahn-Teller效应），需要关掉对称性nosym=.true.，使体系可以弛豫到对称性低的结构。
方法六：
尽量不要对大体系中所有原子进行晶格优化，例如在表面slab模型中，对于slab内部或基底原子固定不优化。
另二篇
存在没有bug的程序吗？
大师如是说：
“任何一个程序，无论它多么小，总存在着错误。”
初学者不相信大师的话，
“如果一个程序小得只执行一个简单的功能，那么会怎样？”他问。
“这样一个程序将没有意义，”大师说，“但假设这样一个程序存在的话，操作系统最后将失效。产生一个错误。”
但初学者不满足。
“如果操作系统不失效，那么会怎样？”他问。
“没有不失效的操作系统，” 大师说，“但假设这样一个操作系统存在的话，硬件最后将失效，产生一个错误。”
初学者仍不满足。
“如果硬件不失效，那么会怎样？”他问。
大师长叹一声。
“没有不失效的硬件，”他说，“但假设这样的硬件存在的话，用户就会想让这个程序做一件不同的事，这件事也是一个错误！”
没有错误的程序是一则谬论，世间难寻。假设存在着一个没有任何错误的程序，那么这个世界将会不复存在。
—-《编程之禅》第四篇（金）第二节
存在一种方法能让自洽计算总是收敛吗？
大师如是说：
“任何自洽场计算程序，无论如何，总存在不收敛的情况。因为当体系复杂，原子个数多时将难以收敛。”
初学者不相信大师的话，
“如果自洽场计算程序只计算一个简单的体系，那么会怎样？”他问。
“这样一个自洽场计算程序将没有意义，”大师说，“但假设这样一个自洽场程序存在的话，计算会由于截断能过高、k点过多、收敛阈值过高造成不收敛。”
但初学者不满足。
“如果收敛相关要求设置不高，那么会怎么样？”他问。
“没有一个自洽场计算程序会限制收敛量的设置，” 大师说，“但假设这样一个自洽场程序存在的话，计算会由于使用的近似不足，造成误差，从而不符合实际，因此你得重新计算，这也是一种不收敛。”
初学者仍不满足。
“如果误差在可以容忍范围，那么会怎样？”他问。
大师长叹一声。
“这样一个自洽场计算结果将缺少物理意义，”他说，“但假设这样一个自洽场计算程序存在的话，客户就会让这个程序计算一个不同的性质，计算一个不同的模型，这也是一种不收敛！”
总是收敛的自洽场计算程序是一则谬论，世间难寻。假设存在着一个总是收敛的自洽场程序，那么这个世界将会不复存在。
本文总阅读量次