vasp+phonopy-QHA計(jì)算材料熱膨脹系數(shù)與格林奈森常數(shù)等腳本
?腳本分享
通過使用QHA方法可以計(jì)算得到材料的熱膨脹系數(shù)，熱容和格林奈森參數(shù)等數(shù)據(jù)，這里使用的軟件為vasp和phonopy。通過腳本一次計(jì)算得到相關(guān)數(shù)據(jù)。
相關(guān)參考為：
http://phonopy.github.io/phonopy/examples.html
計(jì)算主要流程為:以聲子計(jì)算無虛頻的結(jié)構(gòu)文件為基礎(chǔ)，更改其縮放系數(shù)，在不同體積下計(jì)算其原子間受力、聲子結(jié)構(gòu)和熱性質(zhì)。以提取出的e-V的數(shù)據(jù)為基礎(chǔ)，使用phonopy-qha計(jì)算包括熱膨脹系數(shù)、熱容和格林艾森參數(shù)等在內(nèi)的熱學(xué)性質(zhì)。
首先將可計(jì)算得到無虛頻的vasp輸入文件放置在文件夾中，包括POSCAR、INCAR、KPOINTS和POTCAR，以及phonopy處理需要的設(shè)置文件，我這里使用vaspkit ?305功能產(chǎn)生的KPATH.phonopy文件當(dāng)作設(shè)置文件（請(qǐng)根據(jù)你的使用習(xí)慣更改腳本）。
KPATH.phonopy文件內(nèi)容如下：
NPOINTS = 501
DIM =? 1 1 1? ? ?
#擴(kuò)胞倍數(shù)隨具體情況而定
BAND = 0.000000 0.000000 0.000000 0.500000 0.000000 0.500000 0.625000 0.250000 0.625000, 0.375000 0.375000 0.750000 0.000000 0.000000 0.000000 0.500000 0.500000 0.500000 0.500000 0.250000 0.750000 0.500000 0.000000 0.500000
BAND_LABELS = $\Gamma$ X U K $\Gamma$ L W X
#聲子色散路徑自行修改
MP = 11 11 11
TETRAHEDRON = .TRUE.
#PDOS = 1, 2
BAND_CONNECTION = .TRUE.
FORCE_CONSTANTS = READ
ATOM_NAME = W
DIM = 1 1 1
TPROP=.TRUE.
TMAX=2000
TSTEP=5
FORCE_CONSTANTS = READ
# FORCE_SETS = READ
# IRREPS = 0? 0? 0
# SHOW_IRREPS = .TRUE.
# LITTLE_COGROUP = .TRUE.
計(jì)算腳本以及腳本注釋如下：
#!/bin/bash
mkdir QHA
#建立QHA計(jì)算文件夾
for i in $(seq 0.95 0.005 1.05)
#對(duì)i賦值0.95到1.05之間，每隔0.005取值。
do
mkdir $i
cp POSCAR? INCAR KPOINTS POTCAR KPATH.phonopy ./$i
cd ./$i
sed '2c '$i'' POSCAR > POS
mv POS POSCAR
#修改POSCAR中的縮放系數(shù)，并覆蓋POSCAR
pwd
mpirun -np 48 vasp_std |tee runlog
#執(zhí)行vasp計(jì)算
pwd
phonopy --fc? vasprun.xml?
#獲得力常數(shù)
phonopy -c POSCAR? -s? -p KPATH.phonopy
#獲得聲子色散和聲子態(tài)密度
phonopy -t KPATH.phonopy |tee? thermo.dat
#獲得熱學(xué)性質(zhì)
cp thermal_properties.yaml ../QHA/thermal_properties-$i.yaml
#轉(zhuǎn)移文件以便后續(xù)計(jì)算
d=$(awk 'NR==3{print $1}' CONTCAR)
V=$(echo `awk -v x=$i -v a=$d -v b=$d -v c=$d 'BEGIN{printf "%.14f\n",x*x*x*a*b*c}'`)?
#獲得當(dāng)前縮放系數(shù)下體積
E=$(grep "TOTEN" OUTCAR | tail -1| awk '{printf "%12.6f \n", $5}')
#獲得當(dāng)前縮放系數(shù)下能量
echo $V $E >> ../QHA/v-e.dat
cd ..
done
cd ./QHA/
phonopy-qha v-e.dat thermal_properties*.yaml |tee thermo.dat
#執(zhí)行QHA計(jì)算
注意：在筆者進(jìn)行測(cè)試時(shí)發(fā)現(xiàn)在部分情況下會(huì)出現(xiàn)體積計(jì)算值為0的結(jié)果，需要根據(jù)實(shí)際情況調(diào)整相關(guān)的計(jì)算方式。
計(jì)算結(jié)果展示在終端窗口中
phonopy案例圖圖如下
計(jì)算得到的其它熱學(xué)性質(zhì)存儲(chǔ)在以下文件中
可根據(jù)文件名查找。
另外，如果修改縮放系數(shù)聲子計(jì)算結(jié)果出現(xiàn)虛頻，phonopy-qha計(jì)算時(shí)會(huì)提醒
# Warning: thermal_properties-1.040.yaml has imaginary modes.
# Warning: thermal_properties-1.045.yaml has imaginary modes.
# Warning: thermal_properties-1.050.yaml has imaginary modes.
可根據(jù)實(shí)際情況修改縮放系數(shù)的區(qū)間。
愿有所成
IEchoQ
引喻失義? ?妄自菲薄
