VASPsol++溶剂计算入门

本专题用鸢算科技集群实际计算的H₂O、H₃O⁺与Cu(111)小模型，示范VASPsol++的非线性、非局域电解质和恒电势操作。使用 vasp/6.3.2-8173m-vtst-sol++ 的 vasp_std，原构建采用Intel ifort / oneAPI 2023.2和AVX-512，当前模块提供oneAPI 2024.2运行环境。它是CPU扩展构建，不能把其他版本VASP GPU的模块直接替换进本教程。

下载完整输入、数据、9张图和提取脚本，保存为 vaspsol-plus.zip：

python3 -m zipfile -e vaspsol-plus.zip .
cd vaspsol-plus-practical-wave283
bash preview.sh
ls preview/figures

examples/保存48个实际条件的INCAR、POSCAR和KPOINTS，reference/是相应数值表；预览只读取这些数值，不重新计算。POTCAR和VASP程序需使用您自己的授权资源。新的DFT任务通过Slurm提交。

1 先准备气相参照#

流程是：先准备结构与POTCAR → 优化气相原子位置 → 将同一几何用于各溶剂单点 → 提取能量和场分布。这里气相优化只移动原子，14 Å晶胞不变；电子收敛为1E-6 eV，最大原子力小于0.02 eV/Å。

H₂O的O/H价电子合计8；H₃O⁺的中性价电子合计9，设 NELECT=8 表示少一个电子。不能仅凭目录名中的“离子”就假定软件自动改变电子数。

cp -r examples/water-gas-relax water-gas
cp submit.sh water-gas/
cd water-gas
# 按POSCAR元素顺序准备O、H的有授权POTCAR
sbatch submit.sh

submit.sh沿用CPU分区56个MPI进程、数学库单线程的正式入口；提交前用 sinfo -p cpu-8173m 确认资源可用。本文不提供该CPU节点的性能数据。不要在登录节点直接运行生产DFT。

气相H₂O的完整INCAR：

SYSTEM=water-gas-reference
ENCUT=450
PREC=Accurate
EDIFF=1E-6
EDIFFG=-0.02
NELECT=8
NBANDS=12
ISPIN=1
ISMEAR=0
SIGMA=0.05
ALGO=Normal
NELM=180
IBRION=2
POTIM=0.5
NSW=80
ISIF=2
ISYM=0
LREAL=.FALSE.
LASPH=.TRUE.
LWAVE=.FALSE.
LCHARG=.FALSE.
LSOL=.FALSE.

本例气相H₂O共3次力计算，最后最大力0.010877 eV/Å，E0为−14.23103574 eV；H₃O⁺共8次，最后最大力0.004183 eV/Å，E0为−10.63810511 eV。带电周期盒含背景与有限尺寸效应，这不是已收敛的绝对溶剂化自由能标定。

2 在同一结构上做溶剂单点#

以0.1 mol/L电解质中的H₃O⁺为例，下面是实际使用的完整输入。

INCAR：

SYSTEM=hydronium-nonlinear-01M
ENCUT=450
PREC=Accurate
EDIFF=1E-6
NELECT=8
NBANDS=12
ISPIN=1
ISMEAR=0
SIGMA=0.05
ALGO=Normal
NELM=180
IBRION=-1
NSW=0
ISYM=0
LREAL=.FALSE.
LASPH=.TRUE.
LWAVE=.TRUE.
LCHARG=.FALSE.
LSOL=.TRUE.
LSOL_SCF=.TRUE.
ISOL=2
C_MOLAR=0.1
R_ION=4
LVHAR=.TRUE.

POSCAR：

O  H                                    
   1.00000000000000     
    14.0000000000000000    0.0000000000000000    0.0000000000000000
     0.0000000000000000   14.0000000000000000    0.0000000000000000
     0.0000000000000000    0.0000000000000000   14.0000000000000000
   O    H 
     1     3
Direct
  0.5000032075302921  0.5000000009152141  0.4975238326463285
  0.5677832618111174  0.4999999998158029  0.4779673459451298
  0.4661067660176011  0.5587041054579234  0.4779686969402537
  0.4661067646409967  0.4412958938110597  0.4779686958968637
 
  0.00000000E+00  0.00000000E+00  0.00000000E+00
  0.00000000E+00  0.00000000E+00  0.00000000E+00
  0.00000000E+00  0.00000000E+00  0.00000000E+00
  0.00000000E+00  0.00000000E+00  0.00000000E+00

KPOINTS：

Molecular reference
0
Gamma
1 1 1
0 0 0

LSOL打开溶剂，ISOL=2选择非线性、非局域模型，C_MOLAR指定浓度，R_ION=4设置离子空腔半径。默认溶剂参数对应水；把介电常数换成另一溶剂的数值，并不等于已为那种溶剂重新参数化全部模型。

cd ..
cp -r examples/screening/hydronium-nonlinear-01M ion-01M
cp submit.sh ion-01M/
cd ion-01M
# 按O、H顺序准备POTCAR
sbatch submit.sh
tail OSZICAR
grep -A35 'Solvation model:' OUTCAR

此例 NSW=0，只计算电子与溶剂响应，不优化离子位置；力可能大于0.02 eV/Å，不能把它叫作“溶剂中优化完成的结构”。如需溶剂中的平衡几何，应另行设置原子优化，并检验最后的力。

3 比较实际结果#

这些是固定几何的电子能量差，未包含标准态、振动零点能、熵等完整热力学修正。ISOL=1与2的默认 NC_K、TAU也不同，图中二者的差异不能全部归因于“线性/非线性”一项。

下一步看屏蔽与场分布或恒电势操作。