py4vasp 空间屏蔽与诱导电流

用水分子实际NICS与磁场诱导电流网格，画完整分子平面、局部放大及矢量图。先用实际数据完成提取和绘图，再按相同顺序准备自己的计算。

1. 先用现成数据走通流程#

下载真实空间数据与终端实操包，解压进入资料目录：

unzip py4vasp-realspace-practical-wave289-*.zip
cd py4vasp-realspace-practical-wave289
bash preview.sh

先在 preview/figures/ 看5张静态图，再执行：

bash run.sh
ls results/data
ls tools-output

run.sh依次准备读取载体、提取数值、调用py4vasp原生绘图、生成静态PNG。四步也可分开运行 bash prepare.sh、bash extract.sh、bash tools.sh、bash plot.sh。脚本保护已存在的输出，重复操作请换一个新解压目录。后处理使用已分配资源内1个CPU；不会重新运行VASP，也不需要GPU。

后处理模块为 py4vasp/0.11.3（Python3.11独立环境）；静态图用 pymatgen/2026.9.24 环境的NumPy2.3.5、Matplotlib3.11.2。HTML用浏览器打开，通过CDN加载Plotly；离线查看可直接用PNG。

包内 carriers/*.npz 是本次实际计算导出的自有数值。prepare.sh在work/重建只包含指定数组的最小HDF读取载体，并写明来源；它不是原始完整VASP输出，也不能用来读取未提供的物性。包中不含PAW、波函数或VASP程序。reference/保留独立提取的数据，inputs/提供各例完整INCAR、POSCAR、KPOINTS。

2. 生成空间磁响应网格#

本例使用3原子水分子、正交真空晶胞14×15.526478×14.596309 ų，550 eV、Gamma点、EDIFF=1E−6。VASP6.6.1 CPU/HDF5（oneAPI2024.2、MPI、HDF51.14.6）以1个MPI进程计算，基态电子收敛23步；随后完成磁响应。原生网格为 216×240×224，包含9分量空间屏蔽及外场Bx/By/Bz各自的3分量诱导电流。

LCHIMAG = .TRUE.
LNMRSHIELD = .TRUE.
NUCIND = .TRUE.
WRT_NMRCUR = 4
LLRAUG = .TRUE.
ISYM = 0

NUCIND和WRT_NMRCUR的这些网格输出需要对应新版本支持。CPU响应运行前增大程序栈；此HDF5构建使用MPI，设置OMP_NUM_THREADS不会启用OpenMP。提交脚本中使用：

ulimit -s unlimited
export OMP_STACKSIZE=256M
export OMP_NUM_THREADS=1
# 在自己的Slurm资源分配内执行，使用CPU/HDF5版本
mpirun -np 1 vasp_std > vasp.log

本包已有真实响应数据，bash run.sh只后处理，不重新执行这项响应。若要另算体系，应自行准备PAW并提交CPU作业；作业线程数和资源申请须一致。

3. 读取张量，明确NICS正负号#

module load py4vasp/0.11.3
python - <<'PYCODE'
import numpy as np
import py4vasp
c = py4vasp.Calculation.from_path('work/Water-current-NICS')
sigma = c.nics.to_numpy()
sigma_iso = c.nics.to_numpy('isotropic')
print(sigma.shape, sigma_iso.shape)
print('trace/3差:', np.max(np.abs(
    sigma_iso-np.trace(sigma,axis1=-2,axis2=-1)/3)))
chemical_nics = -sigma_iso
np.save('sigma-isotropic.npy', sigma_iso)
np.save('chemical-NICS.npy', chemical_nics)
PYCODE

原生字段是屏蔽张量σ，单位ppm，本图画 σiso=trace(σ)/3。若使用化学文献中NICS=−σ的约定，需要显式取负，不能把σ图改个标题就当NICS符号。这个水分子示例用于演示空间响应操作，不能据此研究芳香性；网格也不能直接替代核位置上含完整PAW贡献的核屏蔽结果。

4. 诱导电流是矢量，不是一个分量#

module load py4vasp/0.11.3
python - <<'PYCODE'
import py4vasp
c = py4vasp.Calculation.from_path('work/Water-current-NICS')
j = c.current_density.read('nmr')
for name in ['current_bx', 'current_by', 'current_bz']:
    print(name, j[name].shape)
c.current_density.to_contour('NMR(bx)', a=13/27).to_plotly().write_html(
    'current-Bx.html', include_plotlyjs='cdn')
c.current_density.to_quiver('NMR(bx)', a=13/27).to_plotly().write_html(
    'current-vectors.html', include_plotlyjs='cdn')
PYCODE

0.11.3的选择写 NMR(bx)，不是其他版本示例中的NMR(x)。current_bx指外场沿Bx时的(jx,jy,jz)整条矢量，末轴长度3，并非只有jx。单位为hartree/bohr²的磁响应约定；它不是外加偏压输运电流，不能直接标成nA。空间网格为平面波部分，不包含完整一中心电流。

水分子位于原始a=0.5的yz平面。图中底色为完整三分量电流模，箭头为该面内jy、jz，箭头长度表示相对大小。右图明确放大分子区域，左图保留整个晶胞。

5. 大网格怎样交付与重画#

为减少下载和原生交互图的开销，包内三维读取载体精确保留每8个原始网格点，无插值，对应27×30×28，交互工具示例在a=13/27截面。这个截面不能冒充原始分子面a=0.5。

reference/Water-molecular-plane.npz另行保存原始完整分辨率的a=0.5分子平面及原子位置，本文两张静态图均由该完整平面重画。不要用抽样三维载体重新计算核屏蔽、积分峰值或宣称网格收敛。需要自己的完整体数据时，在自己的原生输出目录调用同样接口。

bash tools.sh
bash plot.sh
ls tools-output/Water-*.html
ls plots/figures/Water-*.png

tools.sh若已运行过，请在新解压目录操作；若只需看已有图，可直接打开随包的figures/，或先运行preview.sh。

完整计算输入#

下列为对应原生计算的输入，PAW和前置电荷/波函数请按正文说明自行准备。后处理包默认不运行这些计算。

SYSTEM=Water-NMR-current-NICS
ENCUT=550
PREC=Accurate
EDIFF=1E-6
EDIFFG=-0.02
IBRION=-1
NSW=0
ISIF=2
ISMEAR=0
SIGMA=0.05
LREAL=.FALSE.
LASPH=.TRUE.
LWAVE=.TRUE.
LCHARG=.TRUE.
NELM=160
ALGO=Normal
ISYM=0
LCHIMAG=.TRUE.
LNMRSHIELD=.TRUE.
NUCIND=.TRUE.
WRT_NMRCUR=4
LLRAUG=.TRUE.

H  O 
 1.0000000000000000
    14.0000000000000000    0.0000000000000000    0.0000000000000000
     0.0000000000000000   15.5264780000000009    0.0000000000000000
     0.0000000000000000    0.0000000000000000   14.5963089999999998
 H   O  
   2   1
Direct
  0.5000000008650107  0.5494280959188510  0.4795763630523225
  0.5000000008711729  0.4505719043605292  0.4795763636714860
  0.4999999982638164  0.4999999997206199  0.5204205681413534
Cartesian
  0.0000000000000000  0.0000000000000000  0.0000000000000000
  0.0000000000000000  0.0000000000000000  0.0000000000000000
  0.0000000000000000  0.0000000000000000  0.0000000000000000

Gamma mesh
0
Gamma
1 1 1
0 0 0

继续阅读#

• py4vasp 势场分量与截面图
• py4vasp 选带电荷密度
• py4vasp 表面功函数
• 全部quantity与文件前提