离子弛豫
Quantum ESPRESSO · Si 固定胞弛豫：末态力与内部坐标
固定面内条件前，先学会读取原子弛豫
离子弛豫利用原子力更新位置，同时保持晶胞不变。本例从 Si SCF 的两原子结构出发，将第二个 Si 沿 x 方向移开约 0.108 Å，再观察 BFGS 优化怎样减小残余力。第一个原子固定，用来去除整体平移自由度。截断能与 k 网格沿用入门算例，参数比较见收敛测试。
界面堆叠比较和应变扫描常常规定面内晶格，只允许层间距和内部坐标调整。Si 的这次固定晶胞位移测试先把这一操作缩到两个原子：如果直接使用对称位置，原子力可能已经为零，便看不到优化如何修复内部坐标。本例有意加入位移，检查给定晶胞中能否回到小残余力的结构，为后续固定结构性质建立起点。Giannozzi 等的 QE 方法论文第 4.1 节区分原子坐标与晶胞自由度；本页只开放前者，末态不能回答平衡体积是多少。
pw.x 输入说明 · PWscf 用户手册 · QE 的 Si 结构示例
完整算例包保留真实输入、OUT、XML、轨迹数据及原绘图脚本。解包后在 si-pbe 读取下面的 relax 和 relax-check；赝势沿用 SCF 页的官方来源。
下载包保留输入、输出、XML 与作图数据，未打包 tmp/si.save 中的电荷密度和波函数。阅读输出、重新作图可直接使用包内文件；重新计算时，按本页完整的 relax 输入，从有位移的初始结构生成电子态与优化轨迹。
修改坐标与离子优化设置
原始对称位置是 (0.25, 0.25, 0.25)，输入改为 (0.27, 0.25, 0.25)，单位为 alat。本例晶格常数为 5.397607551 Å，这个位移约 0.108 Å。下面是实际运行的完整文件。从零准备时，先在 si-pbe 下用 mkdir relax 建立计算目录；下载包中的 relax 已有原始输入输出，复算应另建目录，并保持 ../pseudo 指向赝势目录。
[preston@preston-System-Product-Name si-pbe]$ cp scf/scf.in relax/relax.in
[preston@preston-System-Product-Name si-pbe]$ vi relax/relax.in
[preston@preston-System-Product-Name si-pbe]$ cat relax/relax.in
&CONTROL
calculation = 'relax'
etot_conv_thr = 1.0d-7
forc_conv_thr = 1.0d-4
nstep = 40
prefix = 'si'
outdir = './tmp'
pseudo_dir = '../pseudo'
tprnfor = .true.
tstress = .true.
/
&SYSTEM
ibrav = 2
A = 5.397607551
nat = 2
ntyp = 1
ecutwfc = 60
ecutrho = 640
occupations = 'fixed'
/
&ELECTRONS
conv_thr = 1.0d-10
/
&IONS
ion_dynamics = 'bfgs'
/
ATOMIC_SPECIES
Si 28.085 Si.pbe-n-rrkjus_psl.1.0.0.UPF
ATOMIC_POSITIONS alat
Si 0.00 0.00 0.00 0 0 0
Si 0.27 0.25 0.25 1 1 1
K_POINTS automatic
8 8 8 0 0 0
[preston@preston-System-Product-Name si-pbe]$
相对于 SCF，关键变化在 calculation='relax'、&IONS 和坐标。ion_dynamics='bfgs' 指定离子优化算法；etot_conv_thr=1.0d-7 的单位为 Ry，比较相邻离子步的整胞能量变化，forc_conv_thr=1.0d-4 的单位为 Ry/Bohr，检查允许移动分量上的力。两类条件都要满足；OUT 中的 Total force 只是帮助看趋势，不能代替逐分量检查。这里的力阈值约为 0.00257 eV/Å，后面还会用独立静态计算核对。
nstep=40 给出允许的离子步数上限，达到上限不等于优化完成。&ELECTRONS 中的 conv_thr=1.0d-10 则决定每个离子位置上的电子自洽精度。若电子误差造成的力变化已经接近离子阈值，只继续收紧力阈值会让优化难以稳定结束；因此下面保留了收紧电子阈值后的实际力对照。
原子每移动一次，电子密度所对应的外势就改变，需要在新坐标上重新自洽再求力。离子步外面的 BFGS 更新与里面的电子迭代因此嵌套：增加电子步数解决当前坐标的电子求解，增加离子步数允许更多位置更新。若某个结构的电子循环未求完，后面的力不能仅凭比上一点小就接受；末态也要用自己的坐标复核，不能沿用初始结构的 SCF 密度代替。
第一个 Si 后面的 0 0 0 固定三个分量，第二个 Si 的 1 1 1 允许三个分量移动。输入里没有 &CELL；晶格由 ibrav=2 和 A 固定。这次计算只能回答“给定这个晶胞，原子是否回到力较小的位置”，不能回答平衡晶格常数是多少。要让晶胞参与优化，接晶格优化。
提交优化并读取每一步的力
[preston@preston-System-Product-Name si-pbe]$ cat relax/run.sh
#!/bin/bash
#SBATCH --job-name=atlas-si-relax
#SBATCH --nodes=1
#SBATCH --ntasks=4
#SBATCH --cpus-per-task=1
#SBATCH --time=00:20:00
#SBATCH --output=_out.%j.log
#SBATCH --error=_err.%j.log
unset DISPLAY XAUTHORITY
ulimit -s unlimited
ulimit -c 0
export OMP_NUM_THREADS=1
export OPENBLAS_NUM_THREADS=1
export MKL_NUM_THREADS=1
cd "$SLURM_SUBMIT_DIR"
/usr/bin/mpirun --bind-to none -np 4 <qe_bin>/pw.x -in relax.in > relax.out 2> relax.err
[preston@preston-System-Product-Name si-pbe]$
输入和脚本可分别下载：relax.in、run.sh。这里使用 4 个 MPI 进程。提交后先看队列，再追踪 relax.out；空的队列只说明任务不在运行，仍要回到输出中查结束原因。
[preston@preston-System-Product-Name si-pbe]$ cd relax
[preston@preston-System-Product-Name relax]$ sbatch run.sh
Submitted batch job 777
[preston@preston-System-Product-Name relax]$ cd ..
先从 OUT 头部核对这次读到了哪一份输入，以及实际运行版本和进程数。不要把上一次留下的文件当成刚提交任务的输出。
[preston@preston-System-Product-Name si-pbe]$ head -n 30 relax/relax.out
Program PWSCF v.7.5 starts on 22Sep2026 at 21:35:40
This program is part of the open-source Quantum ESPRESSO suite
for quantum simulation of materials; please cite
"P. Giannozzi et al., J. Phys.:Condens. Matter 21 395502 (2009);
"P. Giannozzi et al., J. Phys.:Condens. Matter 29 465901 (2017);
"P. Giannozzi et al., J. Chem. Phys. 152 154105 (2020);
URL http://www.quantum-espresso.org",
in publications or presentations arising from this work. More details at
http://www.quantum-espresso.org/quote
Parallel version (MPI), running on 4 processors
MPI processes distributed on 1 nodes
7463 MiB available memory on the printing compute node when the environment starts
Reading input from relax.in
Current dimensions of program PWSCF are:
Max number of different atomic species (ntypx) = 10
Max number of k-points (npk) = 40000
Max angular momentum in pseudopotentials (lmaxx) = 4
R & G space division: proc/nbgrp/npool/nimage = 4
Subspace diagonalization in iterative solution of the eigenvalue problem:
a serial algorithm will be used
Parallelization info
[preston@preston-System-Product-Name si-pbe]$
relax.out 的主体是“电子自洽 → 力 → 更新原子位置 → 下一次电子自洽”反复出现。下面把这次的力和 BFGS 结束信息读出来。初始总力为 0.058861 Ry/Bohr；第一步后降到 0.042098，再到 0.023403。第四个电子周期已经接近极小值，最后一个电子周期才出现 BFGS 收敛信息。
[preston@preston-System-Product-Name si-pbe]$ grep -E 'Total force|bfgs converged|Final energy' relax/relax.out
Total force = 0.058861 Total SCF correction = 0.000001
Total force = 0.042098 Total SCF correction = 0.000003
Total force = 0.023403 Total SCF correction = 0.000001
Total force = 0.000193 Total SCF correction = 0.000001
Total force = 0.000000 Total SCF correction = 0.000000
bfgs converged in 5 scf cycles and 4 bfgs steps
Final energy = -22.8385922964 Ry
[preston@preston-System-Product-Name si-pbe]$
这里明确写了 5 个 SCF 周期、4 步 BFGS。一个 SCF 周期求解当前坐标的电子态；一个 BFGS 步用这些力更新坐标，电子迭代号会在下一个结构重新计数。初始结构的电子循环第 8 轮曾出现三条 c_bands: 3 eigenvalues not converged，第 9、10 轮未再出现，随后给出电子收敛行；后续四个结构的电子循环也分别收敛。检查这类警告时，要连着读它所在的电子轮次和之后的收敛记录，具体说明见 QE 故障排查。BFGS 收敛行直接说明了离子优化的停止原因。若输出写的是到达最大步数，或 bfgs 没有收敛，即使程序已经结束，也不能把那份结构当成通过优化。
读取最终坐标，复核末态力
继续读最终坐标。第二个 Si 的 x 分量回到 0.2500000911，y、z 保持在 0.25；第一个原子保持固定。由于这是 relax，最终坐标段没有重新优化出的晶胞。
[preston@preston-System-Product-Name si-pbe]$ grep -A8 'Begin final coordinates' relax/relax.out
Begin final coordinates
ATOMIC_POSITIONS (alat)
Si 0.0000000000 0.0000000000 0.0000000000 0 0 0
Si 0.2500000911 0.2500000000 0.2500000000
End final coordinates
[preston@preston-System-Product-Name si-pbe]$
不要把前面打印的 Total force = 0.000000 理解成数学上的零。逐原子力保留了更多有用信息：
[preston@preston-System-Product-Name si-pbe]$ grep -A8 'Forces acting' relax/relax.out | tail -n 9
Forces acting on atoms (cartesian axes, Ry/au):
atom 1 type 1 force = 0.00000029 0.00000000 0.00000000
atom 2 type 1 force = -0.00000029 0.00000000 0.00000000
Total force = 0.000000 Total SCF correction = 0.000000
SCF correction compared to forces is large: reduce conv_thr to get better values
[preston@preston-System-Product-Name si-pbe]$
第二个原子的残余 x 力约为 −2.9×10⁻⁷ Ry/Bohr，小于本次输入的 1.0×10⁻⁴ Ry/Bohr 条件。后面的提示也保留下来：SCF 修正相对已经很小的残余力仍可能显得大，因此不能用这些末位数字讨论极高精度的力。
为核对这个提示，本次另外保留了 relax-check，使用最终坐标做固定结构 SCF，并把 conv_thr 收紧为 1.0d-12。这是一份新输入、新输出，不覆盖刚才的优化。实时查看时使用了下面的命令，看到电子迭代从第 7 次继续向后推进；按 Ctrl-C 只退出 watch，不会取消提交给 Slurm 的作业。
[preston@preston-System-Product-Name si-pbe]$ watch -n 2 "tail -n 8 relax-check/scf.out"
收紧电子阈值后的力如下。读这一段时，应比较它是否仍远小于离子收敛条件，而不是要求它和上一次输出的最后一位完全相同。
[preston@preston-System-Product-Name si-pbe]$ grep -A8 'Forces acting' relax-check/scf.out
Forces acting on atoms (cartesian axes, Ry/au):
atom 1 type 1 force = -0.00000000 0.00000000 -0.00000000
atom 2 type 1 force = 0.00000000 -0.00000000 0.00000000
Total force = 0.000000 Total SCF correction = 0.000000
Computing stress (Cartesian axis) and pressure
[preston@preston-System-Product-Name si-pbe]$
最后再看原优化输出的末尾。原生程序用时约 1 分 42 秒；这次前段和其他短作业同时运行，WALL 时间还包含了资源竞争。不能只用这个时间推算复杂材料的优化成本。
[preston@preston-System-Product-Name si-pbe]$ tail -n 12 relax/relax.out
interpolate : 0.17s CPU 0.55s WALL ( 39 calls)
Parallel routines
PWSCF : 51.81s CPU 1m41.71s WALL
This run was terminated on: 21:37:22 22Sep2026
=------------------------------------------------------------------------------=
JOB DONE.
=------------------------------------------------------------------------------=
[preston@preston-System-Product-Name si-pbe]$
保留优化轨迹，进入固定结构计算
完整 relax.out保留每一步坐标、能量与力。逐步数据表中的第一点是初始结构，5 个电子周期对应 4 次 BFGS 位置更新：
| 电子周期 | 能量 (Ry/原胞) | Total force (Ry/Bohr) |
|---|---|---|
| 1 | -22.83255668 | 0.058861 |
| 2 | -22.83552923 | 0.042098 |
| 3 | -22.83765081 | 0.023403 |
| 4 | -22.83859223 | 0.000193 |
| 5 | -22.83859230 | 0.000000 |
末轮 Total force=0.000000 是打印精度下的值；逐原子 x 分力为约 ±2.9×10⁻⁷ Ry/Bohr，小于本次 1×10⁻⁴ Ry/Bohr 的力阈值。上面的紧电子阈值静态复核用于检查末态力对电子求解精度的敏感性。
第 4→5 个电子周期的能量只变化 0.00000007 Ry，表中总力却从 0.000193 降到打印为零。能量和力反映极小值附近的不同变化，不能因为能量几乎重合就提前停止。Total force 汇总了力信息，正式比较仍回到允许移动的分量。紧电子阈值复核中逐原子力也落在打印精度内，支持末态力仍低于这次的离子阈值；它没有测定无限精度下的零力。固定晶胞下压力仍为 38.45 kbar，说明这次修复的是内部坐标，求平衡体积应转到变胞问题。
对于规定面内晶格的异质结，应保留该晶格和真空方向，检查所有允许移动原子的末态力，并从最终坐标读取层间距。对同一应变下的各个堆叠使用一致约束，才能把能量差解释为该条件下的构型比较；若开放晶胞，比较的问题就改变了。
Ba₂N 原文 Fig. 1(a)，第 165101-2 页同时给出俯视、侧视、晶胞边界和二维布里渊区。俯视图识别面内排列，侧视图把 Ba–N–Ba 的三层高度分开，便能读清结构图中的层间距与计算晶胞的真空是两回事。将本页末态按同样方式展示时，从最终 ATOMIC_POSITIONS 和固定晶格导出自己的结构文件，在 VESTA 显示晶胞与两个 Si，分别保存俯视和侧视；先把 alat 坐标换算为实际晶格下的位置，不按 Ba₂N 图把 Si 画成三原子层。本页 ±2.9×10⁻⁷ Ry/Bohr 的末力检查解释的是这份几何怎样停止，界面建模页再把这类视图接到层间距和共同晶格；结构图本身不代替力、结合能或声子的计算。