CP2K专题

CP2K 分子红外与拉曼光谱

H2O与CO2的真实结构、完整输入、有限位移、模式与IR/Raman活动数据。

更新于 2026-10-28约 11 分钟阅读按任务与环境选择适用配置
阅读本篇目录

计算分子红外和拉曼,需要先把结构优化到驻点,再计算位移前后的力、偶极矩和极化率。只设置 INTENSITIES TRUE 不会自动替你打开这些输出。本篇用H₂O与CO₂分开说明非线性和线性分子。

下载完整输入、原始输出、频率表和绘图脚本,解压后进入 run.sh 所在目录。本文在鸢算科技P100上用CP2K2026.2双精度 cp2k.psmp 实算;构建采用GNU13.3、MPICH5.0.1、CUDA12.9、OpenBLAS0.3.33、ScaLAPACK2.2.3。分析环境为ASE3.29.0。

1 提交与查看结果#

sbatch submit.slurm water

模板申请四张P100、四个MPI进程,每进程一个CPU核和一个OpenMP线程。NPROC_REP 1让不同进程计算不同位移,脚本按MPI全局编号分配GPU。在已分配的计算资源里也可以执行 bash run.sh water。套餐用户用 yskit 生成对应分区的提交脚本,再保留这条运行命令;不要在登录节点直接运行。

脚本依次执行“优化结构 → 提取最终XYZ → 有限位移力、偶极与极化率 → 振动模式 → 数据提取 → 画图”。inputs为可编辑输入,reference为本文结果,work为本次新结果。每个条件请新解压一份,脚本不会覆盖已有work。

2 先优化,再算振动#

采用PBE、DZVP-MOLOPT-SR-GTH、GTH-PBE赝势、600/50 Ry网格和12 Å非周期MT盒。振动对数值噪声敏感,电子收敛EPS_SCF为1E−8;优化最大力1E−5 Hartree/Bohr,约0.000514 eV/Å,比一般结构优化更紧。

完整输入:water-opt#

&GLOBAL
 PROJECT water
 RUN_TYPE GEO_OPT
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME ../../DATA/BASIS_MOLOPT
  POTENTIAL_FILE_NAME ../../DATA/GTH_POTENTIALS
  &QS
   EPS_DEFAULT 1.0E-12
  &END QS
  &MGRID
   CUTOFF 600
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-8
   MAX_SCF 100
   SCF_GUESS ATOMIC
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
 &END DFT
 &SUBSYS
  &CELL
   ABC 12 12 12
   PERIODIC NONE
  &END CELL
  &COORD
     O         5.9999998739        6.0000000010        5.9992568887
     H         5.9999999865        6.7657755697        6.5949200315
     H         5.9999998417        5.2342244900        6.5949201085
  &END COORD
  &KIND O
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q6
  &END KIND
  &KIND H
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q1
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
 &END PRINT
&END FORCE_EVAL
&MOTION
 &GEO_OPT
  OPTIMIZER BFGS
  MAX_ITER 100
  MAX_FORCE 1.0E-5
  RMS_FORCE 5.0E-6
  MAX_DR 1.0E-4
  RMS_DR 5.0E-5
 &END GEO_OPT
&END MOTION

脚本取 water-pos-1.xyz 的最后一帧写成optimized.xyz。振动输入使用该文件,而不是重新读取最初的未优化坐标。

完整输入:water-ir-raman#

&GLOBAL
 PROJECT water
 RUN_TYPE VIBRATIONAL_ANALYSIS
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME ../DATA/BASIS_MOLOPT
  POTENTIAL_FILE_NAME ../DATA/GTH_POTENTIALS
  &QS
   EPS_DEFAULT 1.0E-12
  &END QS
  &MGRID
   CUTOFF 600
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-8
   MAX_SCF 100
   SCF_GUESS ATOMIC
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
  &PRINT
   &MOMENTS ON
    PERIODIC .FALSE.
    REFERENCE COM
   &END MOMENTS
  &END PRINT
 &END DFT
 &PROPERTIES
  &LINRES
   EPS 1.0E-8
   MAX_ITER 100
   PRECONDITIONER FULL_SINGLE_INVERSE
   &POLAR ON
    DO_RAMAN .TRUE.
    PERIODIC_DIPOLE_OPERATOR .FALSE.
   &END POLAR
  &END LINRES
 &END PROPERTIES
 &SUBSYS
  &CELL
   ABC 12 12 12
   PERIODIC NONE
  &END CELL
  &TOPOLOGY
   COORD_FILE_NAME optimized.xyz
   COORD_FILE_FORMAT XYZ
  &END TOPOLOGY
  &KIND O
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q6
  &END KIND
  &KIND H
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q1
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
 &END PRINT
&END FORCE_EVAL
&VIBRATIONAL_ANALYSIS
 NPROC_REP 1
 DX 0.005
 FULLY_PERIODIC .FALSE.
 INTENSITIES .TRUE.
 THERMOCHEMISTRY .TRUE.
 TC_TEMPERATURE 298.15
 TC_PRESSURE 101325
 &PRINT
  &MOLDEN_VIB ON
  &END MOLDEN_VIB
 &END PRINT
&END VIBRATIONAL_ANALYSIS

DX 0.005单位Bohr,约0.002646 Å。三个原子的中心有限位移共18个正负方向,分配给四个MPI进程;它不是18个独立材料体系。FULLY_PERIODIC FALSE配合分子边界处理刚体自由度。IR需要 DFT/PRINT/MOMENTS,Raman需要 PROPERTIES/LINRES/POLAR,非周期分子采用局域偶极算符。

3 查看原生频率和活性#

grep -E 'VIB\|Frequency|VIB\|IR int|VIB\|Raman' work/cases/water-ir-raman/output.log
cat work/analysis/water-ir-raman.csv

实际优化水分子、红外与非共振拉曼活性谱

4 换成CO₂:保留线性几何,释放振动自由度#

新解压一份,执行 sbatch submit.slurm co2。CO₂优化时只允许三个原子的z坐标调整,保持x=y=6 Å的严格直线;在本例正交固定盒中,XY约束与笛卡尔方向一致。优化结束后,振动输入完全移除该约束,让所有原子三个方向都参与有限位移。若把约束一直带进振动,可能会漏掉弯曲模式。

完整输入:co2-opt#

&GLOBAL
 PROJECT water
 RUN_TYPE GEO_OPT
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME ../../DATA/BASIS_MOLOPT
  POTENTIAL_FILE_NAME ../../DATA/GTH_POTENTIALS
  &QS
   EPS_DEFAULT 1.0E-12
  &END QS
  &MGRID
   CUTOFF 600
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-8
   MAX_SCF 100
   SCF_GUESS ATOMIC
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
 &END DFT
 &SUBSYS
  &CELL
   ABC 12 12 12
   PERIODIC NONE
  &END CELL
  &COORD
 C 6 6 6
 O 6 6 4.82
 O 6 6 7.18
  &END COORD
  &KIND O
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q6
  &END KIND
  &KIND C
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q4
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
 &END PRINT
&END FORCE_EVAL
&MOTION
 &GEO_OPT
  OPTIMIZER BFGS
  MAX_ITER 100
  MAX_FORCE 1.0E-5
  RMS_FORCE 5.0E-6
  MAX_DR 1.0E-4
  RMS_DR 5.0E-5
 &END GEO_OPT
 &CONSTRAINT
  &FIXED_ATOMS
   LIST 1 2 3
   COMPONENTS_TO_FIX XY
  &END FIXED_ATOMS
 &END CONSTRAINT
&END MOTION

完整输入:co2-ir-raman#

&GLOBAL
 PROJECT water
 RUN_TYPE VIBRATIONAL_ANALYSIS
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME ../DATA/BASIS_MOLOPT
  POTENTIAL_FILE_NAME ../DATA/GTH_POTENTIALS
  &QS
   EPS_DEFAULT 1.0E-12
  &END QS
  &MGRID
   CUTOFF 600
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-8
   MAX_SCF 100
   SCF_GUESS ATOMIC
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
  &PRINT
   &MOMENTS ON
    PERIODIC .FALSE.
    REFERENCE COM
   &END MOMENTS
  &END PRINT
 &END DFT
 &PROPERTIES
  &LINRES
   EPS 1.0E-8
   MAX_ITER 100
   PRECONDITIONER FULL_SINGLE_INVERSE
   &POLAR ON
    DO_RAMAN .TRUE.
    PERIODIC_DIPOLE_OPERATOR .FALSE.
   &END POLAR
  &END LINRES
 &END PROPERTIES
 &SUBSYS
  &CELL
   ABC 12 12 12
   PERIODIC NONE
  &END CELL
  &TOPOLOGY
   COORD_FILE_NAME optimized.xyz
   COORD_FILE_FORMAT XYZ
  &END TOPOLOGY
  &KIND O
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q6
  &END KIND
  &KIND C
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q4
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
 &END PRINT
&END FORCE_EVAL
&VIBRATIONAL_ANALYSIS
 NPROC_REP 1
 DX 0.005
 FULLY_PERIODIC .FALSE.
 INTENSITIES .TRUE.
 THERMOCHEMISTRY .TRUE.
 TC_TEMPERATURE 298.15
 TC_PRESSURE 101325
 &PRINT
  &MOLDEN_VIB ON
  &END MOLDEN_VIB
 &END PRINT
&END VIBRATIONAL_ANALYSIS

非线性H₂O应有3N−6=3个内部模式,线性CO₂应有3N−5=4个。CO₂两条弯曲模式几乎简并,不能因为频率重合就只保留一条。本文CO₂参考采用八个MPI/八张P100,公开脚本采用四个MPI/四张P100;改变的是位移的并行分配,仍计算相同18个位移。

CO2实际线性结构、简并弯曲与IR和Raman选择性

5 本例得到什么#

分子/模式 频率/cm⁻¹ IR强度/km·mol⁻¹ Raman活性/Å⁴·amu⁻¹
H₂O 模式1 1612.1373 69.6924 6.04528
H₂O 模式2 3718.1052 19.546 48.6739
H₂O 模式3 3822.7933 95.6427 14.0597
CO₂ 模式1 670.1025 25.1708 1.95968e-10
CO₂ 模式2 670.1033 25.1707 1.41208e-10
CO₂ 模式3 1294.7412 1.60565e-05 15.0776
CO₂ 模式4 2354.8074 515.215 3.39661e-07

CO₂对称伸缩在IR中很弱,反对称伸缩在非共振Raman中很弱;水的三个内部模式均为正频率。结果来自有限基组和网格下的谐近似,不乘经验缩放系数。Raman列是活性,不是加入激光频率、温度因子和仪器响应后的实验散射强度;两种图的纵轴分别归一化,不能互相比较峰高绝对值。

提取与重画#

module load ase/3.29.0
python analyze.py work
python plot_results.py work
ls work/analysis
ls work/figures

CSV给出模式编号、频率、IR强度和Raman活性;JSON保留质量、最终结构和热修正。原生 water-VIBRATIONS-1.mol 是Molden振动文件,可以用支持该格式的分子查看器检查模式位移。项目名前缀water是本包统一输出名,不代表CO₂输入成了水。

只重画本文结果,把上述命令中的work换成reference。脚本用实际模式的线位置与活性叠加高斯曲线,展示宽度σ=10 cm⁻¹;这是绘图展宽,不是计算得到的声子寿命。

设置依据:CP2K振动手册、极化率与Raman。继续看同位素与位移步长和气相热化学。

需要根据你的体系选择资源?联系讨论计算需求。实际可用软件与权限,以当前账户和运行环境为准。

查看配图

放大后可滚动查看细节,按 Esc 关闭。

100%查看原图 ↗