这里给出一个CP2K计算shell脚本：

主要功能是建立一个计算文件夹，在其下实现QS下GPW方法的晶格优化，脚本中可选择结构优化，或者MD，以及对角化或者OT方法，示例材料是六角结构单层InP3的计算，可是CP2K并不能得到和平面波方法一致的晶格结构，所以仅供脚本参考：

#!/bin/bash

#1文件夹名称，并创建此文件夹，同时文件夹名也用于后期计算文本文件输入和输出等标识
file_pre="cp2k_InP3"
file_comp="ubuntu"
file_date="210312"
workdir="${file_pre}_${file_comp}_${file_date}_cellOpt_diagNosmear"
if test ! -d ${workdir}
then
echo "${workdir} not exist and make!!!"
mkdir ${workdir}
else
echo "${workdir} already exist and not make!!!"
fi
cd ${workdir}

#2创建一个文本文件，用于计算结束后提取计算输出文件中的能量等信息
etot_file="${workdir}_etot"
if test -r "${etot_file}"
then
rm "${etot_file}"
fi
cat >> ${etot_file} << EOF2
#../${workdir}/${etot_file}
EOF2
grid_header=true

#3CP2K的运行环境变量，根据自己安装情况选择
ubuntu_set="/安装路径/cp2k/tools/toolchain/install/setup"
cp2k_set_name=${ubuntu_set}
source ${cp2k_set_name}

#openmpi等并行运算的前缀
pre_cp2k="mpirun -np 16  "

#cp2k运行名称：

cp2k_name="cp2k.popt"

#是否是计算中断后restart运算，是则t，否则f
res_j="f"

#4 cp2k的QS的GPW方法重要的收敛参数设置
cutoff_sh=350

#用于参数收敛测试循环
for rel_cutoff_sh in 60
do

file_count="cutoff${cutoff_sh}_relcut${rel_cutoff_sh}"
pre_fix="${file_pre}_${file_count}"
cp2k_filename_in="${file_pre}_${file_comp}_${file_date}_${file_count}.inp"
cp2k_filename_out="${file_pre}_${file_comp}_${file_date}_${file_count}.oup"
#${pre_cp2k} ${cp2k_name} -i ./${restart_cp2k_filename_in} -o ./${restart_cp2k_filename_out}
restart_cp2k_filename_in="${pre_fix}-1.restart"
restart_cp2k_filename_out="${pre_fix}_restart.out"

cat > ${cp2k_filename_in} << EOF0
&GLOBAL
 PROJECT ${pre_fix}
 RUN_TYPE CELL_OPT
 #ENERGY ENERGY_FORCE CELL_OPT GEO_OPT MD
 PRINT_LEVEL LOW 
 !#MEDIUM LOW
&END GLOBAL
&FORCE_EVAL
 METHOD Quickstep
 &SUBSYS
  &KIND Sn
   ELEMENT  Sn
   BASIS_SET DZVP-MOLOPT-SR-GTH
   !#DZVP-MOLOPT-SR-GTH SZV-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE
   !# GTH-PBE GTH-PADE  
  &END KIND
  &KIND Ge
   ELEMENT   Ge
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE
   !# GTH-PBE GTH-PADE-q4
  &END KIND
  &KIND Ga
   ELEMENT   Ga
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE
   !# (GTH-PBE-q13 GTH-PBE) GTH-PBE-q3  GTH-PADE-q3 (GTH-PADE-q13 GTH-LDA-q13 GTH-PADE GTH-LDA)
  &END KIND
  &KIND In
   ELEMENT   In
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE
   !# (GTH-PBE-q13 GTH-PBE) GTH-PBE-q3  GTH-PADE-q3 (GTH-PADE-q13 GTH-LDA-q13 GTH-PADE GTH-LDA)
  &END KIND
  &KIND P
   ELEMENT P   
   BASIS_SET DZVP-MOLOPT-SR-GTH
   !#BASIS_MOLOPT:DZVP-MOLOPT-SR-GTH TZVP-MOLOPT-GTH TZV2P-MOLOPT-GTH
   !#BASIS_pob:pob-TZVP  plus-pob-TZVP
   !#BASIS_MOLOPT_UCL:
   POTENTIAL GTH-PBE
   !#GTH-PADE  GTH-PBE 
  &END KIND
  &KIND O
   ELEMENT   O
   BASIS_SET DZVP-MOLOPT-GTH
   !#BASIS_MOLOPT:DZVP-MOLOPT-GTH TZVP-MOLOPT-GTH TZV2P-MOLOPT-GTH TZV2PX-MOLOPT-GTH
   !#BASIS_pob:
   !#BASIS_MOLOPT_UCL:
   POTENTIAL GTH-PBE
   !#GTH-PADE  GTH-PBE 
  &END KIND

  &CELL
   !#InP3_HEX60_QEdeg0.02relax
   ABC [angstrom] 7.539600 7.539600 20
   ALPHA_BETA_GAMMA 90 90 60
   PERIODIC XYZ
   SYMMETRY HEXAGONAL
   !# CUBIC HEXAGONAL ORTHORHOMBIC TETRAGONAL_AB
   MULTIPLE_UNIT_CELL 1 1 1
  &END CELL
  &COORD
P        0.644488351   0.211022298   0.093732097
P        0.644397199   0.644449544   0.093704383
P        0.211153257   0.644449544   0.093704383
In       0.833294676   0.833412648   0.126252974
In       0.166686801   0.166628397   0.120714080
P        0.355577966   0.355515472   0.153470304
P        0.788905562   0.355515472   0.153470304
P        0.355496188   0.789007625   0.153450475
  #UNIT ANGSTROM
  SCALED 
  &END COORD
  &PRINT
   &ATOMIC_COORDINATES
   &END ATOMIC_COORDINATES
  &END PRINT
 &END SUBSYS
  
 &DFT
  BASIS_SET_FILE_NAME BASIS_MOLOPT_UCL
  BASIS_SET_FILE_NAME BASIS_MOLOPT
  BASIS_SET_FILE_NAME BASIS_MOLOPT_LnPP1
  #! BASIS_MOLOPT BASIS_def2_QZVP_RI_ALL   BASIS_pob
  POTENTIAL_FILE_NAME  POTENTIAL
  &QS
   METHOD GPW
   EPS_DEFAULT 1.0E-10
   !#Default value: 1.00000000E-010
   EXTRAPOLATION ASPC
   !# ps  ASPC
   #EXTRAPOLATION_ORDER 3
  &END QS
  &MGRID
   NGRIDS 5
   CUTOFF ${cutoff_sh}
   REL_CUTOFF ${rel_cutoff_sh}
  &END MGRID
  &XC
   &XC_FUNCTIONAL PBE
   !# PBE PADE
   &END XC_FUNCTIONAL
  &END XC
  &SCF
   SCF_GUESS ATOMIC
   EPS_SCF 1.0E-6
   !#Default value: 1.00000000E-005
   MAX_SCF 200
  &DIAGONALIZATION .TRUE.
    ALGORITHM STANDARD
  &END DIAGONALIZATION
   
   &OT .FALSE.
    PRECONDITIONER FULL_ALL
    !#FULL_KINETIC(use for very large systems-DEFAULT) 
    !#FULL_ALL (Most effective state selective preconditioner based on diagonalization, requires ENERGY_GAP to be an underestimate of the HOMO-LUMO gap. This preconditioner is recommended for almost all systems, except very large systems)
    ENERGY_GAP 0.002
    MINIMIZER DIIS
    !#CG(DEFAULTmost reliable, use for difficult systems) DIIS(less reliable than CG, but sometimes about 50% faster)  SD(not recommended)  BROYDEN
   &END OT

   &MIXING T
    METHOD BROYDEN_MIXING
    !#BROYDEN_MIXING BROYDEN_MIXING_NEW DIRECT_P_MIXING PULAY_MIXING
    ALPHA 0.4
    NBROYDEN 8
   &END MIXING
   
#   ADDED_MOS 40
   &SMEAR .FALSE.
    METHOD FERMI_DIRAC
    !# ENERGY_WINDOW:WINDOW_SIZE FERMI_DIRAC:ELECTRONIC_TEMPERATURE  LIST
    ELECTRONIC_TEMPERATURE [K] 300
   &END SMEAR
   
   &PRINT
    &RESTART OFF
    &END
   &END PRINT
  &END SCF

  &POISSON
   PERIODIC XYZ
   #POISSON_SOLVER PERIODIC
  &END POISSON
  #PLUS_U_METHOD MULLIKEN
  !# MULLIKEN LOWDIN MULLIKEN_CHARGES
 &END DFT
 #!force_eval
 &PRINT
  &FORCES HIGH
  &END FORCES
  &STRESS_TENSOR
  &END STRESS_TENSOR
 &END PRINT

STRESS_TENSOR ANALYTICAL
!#ANALYTICAL DIAGONAL_ANALYTICAL DIAGONAL_NUMERICAL NUMERICAL 
&END FORCE_EVAL

&MOTION
 &MD
 ENSEMBLE NVT
 !# NVE NVT  
 TEMPERATURE [K] 300
 TIMESTEP [fs] 1
 STEPS 6000
 &THERMOSTAT
  REGION GLOBAL
  #TYPE NOSE
 &END THERMOSTAT
 &END MD

 &GEO_OPT
 TYPE MINIMIZATION
 MAX_DR    3.0E-03
 MAX_FORCE 4.5E-04
 RMS_DR    1.5E-03
 RMS_FORCE 3.0E-04
 #MAX_ITER 200
 OPTIMIZER BFGS
  !#CG
  !#Good choice for 'small' systems (use LBFGS for large systems) 
&END GEO_OPT

&CELL_OPT
 TYPE DIRECT_CELL_OPT
  #!DIRECT_CELL_OPT(default)  GEO_OPT MD
 CONSTRAINT Z
 !#YZ Z
 EXTERNAL_PRESSURE  0 0 0  0 0 0  0 0 0
 KEEP_SYMMETRY T
 KEEP_ANGLES T
 MAX_DR 3.0E-003
 MAX_FORCE 4.5E-004
 #RMS_DR    1.5E-03
 #RMS_FORCE 3.0E-04
 #MAX_ITER 200
 OPTIMIZER BFGS
  !#CG LBFGS BFGS(default)
 PRESSURE_TOLERANCE 100
  !#100 bar
&END CELL_OPT

#&CONSTRAINT
# &FIXED_ATOMS
#  COMPONENTS_TO_FIX XYZ
#  LIST 1
# &END FIXED_ATOMS
#&END CONSTRAINT

!# MOTION PRINT
&PRINT
&TRAJECTORY
 &EACH
  MD 2
  !# 1
 &END EACH
&END TRAJECTORY
&VELOCITIES 
 &EACH
  MD 2
  !# 1
 &END EACH
&END VELOCITIES
&FORCES 
 &EACH
  MD 2
  !# 1
 &END EACH
&END FORCES
&RESTART_HISTORY
 &EACH
  MD 500
  !# 500
  CELL_OPT 50
  !# 1
  GEO_OPT  50
  !# 1
 &END EACH
&END RESTART_HISTORY
&RESTART
 BACKUP_COPIES 3
 &EACH
  MD 200
  !# 20
  CELL_OPT 20
  !# 1
  GEO_OPT  20
  !# 1
 &END EACH
&END RESTART
&END PRINT
&END MOTION
EOF0

if [ "${res_j}" = "t" ]
then
echo  "restart the cp2k (${pre_cp2k} ${cp2k_name} -i ./${restart_cp2k_filename_in} -o ./${restart_cp2k_filename_out} ) calculation..."
${pre_cp2k} ${cp2k_name} -i ./${restart_cp2k_filename_in} -o ./${restart_cp2k_filename_out}
echo "done!!!"
else
echo  "running the cp2k ( as "${pre_cp2k} ${cp2k_name} -i ./${cp2k_filename_in} -o ./${cp2k_filename_out}" )"
${pre_cp2k} ${cp2k_name} -i ./${cp2k_filename_in} -o ./${cp2k_filename_out}
echo "done!!!"
fi

total_energy=$(grep -e '^[ \t]*Total energy' ${cp2k_filename_out} | awk '{print $3}')
ngrids=$(grep -e '^[ \t]*QS| Number of grid levels:' ${cp2k_filename_out} | awk '{print $6}')

if ${grid_header} ; then
        printf "# Cutoff (Ry) | Total Energy (Ha)" >> ${etot_file}
        for ((igrid=1; igrid <= ngrids; igrid++)) ; do
            printf "| NG on grid %d" $igrid >> ${etot_file}
        done
        printf "\n" >> ${etot_file}
grid_header=false
fi
printf "%10.2f %10.2f %15.10f" ${cutoff_sh} ${rel_cutoff_sh} $total_energy >> ${etot_file}
for ((igrid=1; igrid <= ngrids; igrid++)) ; do
 grid=$(grep -e '^[ \t]*count for grid' ${cp2k_filename_out} | \
               awk -v igrid=$igrid '(NR == igrid){print $5}')
 printf "  %6d" $grid >> ${etot_file}
done
printf "\n" >> ${etot_file}

done
 

更多推荐

CP2K计算脚本