Quantum Espresso + Phonopy 计算声子过程

以Ammonia (N4H12) 为例。

 

第一步： 松弛结构/ Geometry Optimization

 

Input:

 &control
    calculation = 'relax',
    nstep = 300,
    etot_conv_thr = 1.0d-4,
    forc_conv_thr = 2.0d-3,
    restart_mode = 'from_scratch',
    prefix = 'ammonia',
    outdir = 'Waves',
    pseudo_dir = '/home/z8j/svn_personal/qe/pp'
 /
 &system
    ibrav = 0,
    ntyp = 2,
    nat = 16,
    occupations = 'fixed',
    ecutwfc = 80
 /
 &electrons
    mixing_beta = 0.7,
    conv_thr = 1.d-8,
    electron_maxstep = 1000,
    startingwfc = 'atomic+random'
 /
 &ions
    ion_positions = "default"
 /

ATOMIC_SPECIES
 H  1.000  H.pbe-rrkjus.UPF
 N  14.000 N.pbe-rrkjus.UPF

CELL_PARAMETERS angstrom
 5.1305 0.0000 0.0000
 0.0000 5.1305 0.0000
 0.0000 0.0000 5.1305

K_POINTS automatic
  4 4 4 1 1 1

ATOMIC_POSITIONS angstrom
N       1.0789441   1.0789441   1.0789441
N       1.4863058   4.0515558   3.6441941
N       3.6441941   1.4863058   4.0515558
N       4.0515558   3.6441941   1.4863058
H       1.9095721   1.3477823   0.5710246
H       1.9942253   3.2209279   3.9130323
H       3.1362746   0.6556779   3.7827176
H       4.5594753   4.4748221   1.2174676
H       1.3477823   0.5710246   1.9095721
H       1.2174676   4.5594753   4.4748221
H       3.2209279   3.9130323   1.9942253
H       0.5710246   1.9095721   1.3477823
H       0.6556779   3.7827176   3.1362746
H       4.4748221   1.2174676   4.5594753
H       3.9130323   1.9942253   3.2209279
H       3.7827176   3.1362746   0.6556779

然后Run pw.x:

PBS file:

#!/bin/bash

#PBS -A sns
#PBS -q long
#PBS -m ea
#PBS -M 123123@foxmail.com
#PBS -j oe
#PBS -l qos=condo
#PBS -W group_list=cades-virtues
#PBS -l walltime=50:30:00
#PBS -l nodes=4:ppn=32
#PBS -N ammonia-rlx-qe@/lustre/or-hydra/cades-virtues/z8j/run/rmg_test_ORNL/Ammonia/qe

module load /software/tools/modules/env/cades-virtues
module load /software/tools/modules/compilers/intel/2016.1
module load /software/tools/modules/utils/intel/mkl/2016.1
module load /software/tools/modules/mpi/openmpi/intel/1.10.2
module load /software/user_tools/current/modules/cades-virtues/env/intel
module load /software/user_tools/current/modules/cades-virtues/espresso/5.4.0
module list

export OMP_NUM_THREADS=1

cd /lustre/or-hydra/cades-virtues/z8j/run/rmg_test_ORNL/Ammonia/qe

date

mpirun -np 128 --bind-to none --map-by ppr:32:node:pe=1 /software/user_tools/current/cades-virtues/apps/quantum-espresso/intel/espresso/5.4.0/pw.x < qe.in > qe.00.out

date

 

第二步： 写好用来构建超胞用的单胞/Prepare unit cell file for supercell building

 

这一步比较简单，直接从第一步的output file里读取出现成的原子位置即可，然后修改下必要的参数，如名字 Kpoint mesh。有一点需要注意，phonopy只支持crystal 的原子位置类型，所以需要自己转换一下。

 &control
    calculation = 'scf',
    nstep = 300,
    etot_conv_thr = 1.0d-4,
    forc_conv_thr = 2.0d-3,
    restart_mode = 'from_scratch',
    prefix = 'ammonia',
    outdir = 'Waves',
    pseudo_dir = '/home/z8j/svn_personal/qe/pp'
 /
 &system
    ibrav = 0
    ntyp = 2,
    nat = 16,
    occupations = 'fixed',
    ecutwfc = 80
 /
 &electrons
    mixing_beta = 0.7,
    conv_thr = 1.d-8,
    electron_maxstep = 1000,
    startingwfc = 'atomic+random'
 /

ATOMIC_SPECIES
 H  1.000  H.pbe-rrkjus.UPF
 N  14.000 N.pbe-rrkjus.UPF

CELL_PARAMETERS angstrom
 5.1305 0.0000 0.0000
 0.0000 5.1305 0.0000
 0.0000 0.0000 5.1305

K_POINTS automatic
  4 4 4 1 1 1

ATOMIC_POSITIONS (crystal)
N       0.208255064    0.208255064    0.208255064
N       0.291744916    0.791744916    0.708255064
N       0.708255064    0.291744916    0.791744916
N       0.791744916    0.708255064    0.291744916
H       0.373650671    0.263677961    0.110622992
H       0.389376988    0.626349329    0.763677961
H       0.610622992    0.126349329    0.736322020
H       0.889376988    0.873650671    0.236322020
H       0.263677961    0.110622992    0.373650671
H       0.236322020    0.889376988    0.873650671
H       0.626349329    0.763677961    0.389376988
H       0.110622992    0.373650671    0.263677961
H       0.126349329    0.736322020    0.610622992
H       0.873650671    0.236322020    0.889376988
H       0.763677961    0.389376988    0.626349329
H       0.736322020    0.610622992    0.126349329

这一步不用运行，只是为下一步做准备。

 

第三步： 用Phonopy构建超胞/ Build supercells with displaced atoms by Phonopy

将上一步的input 拷出来并命名为qe.in，然后运行Phonoy（如果想省略掉-c选项的话，重命名为unitcell.in）：

python /home/zjyx/softwares/phonopy-1.10.10/bin/phonopy --pwscf -d --dim="4 4 4" -c qe.in

这样当前目录下会出现若干个supercell-???.in文件。这些文件都是带有微扰的超胞，需要在下一步中用来做电子自洽计算。另外这些文件都只含有结构信息，缺少设置信息，需要把下面的信息添加为每个supercel.???.in的开头部分：

 &control
    calculation = 'scf',
    nstep = 300,
    etot_conv_thr = 1.0d-4,
    forc_conv_thr = 2.0d-3,
    tprnfor = .true.
    tstress = .true.
    restart_mode = 'from_scratch',
    prefix = 'ammonia',
    outdir = 'Waves',
    pseudo_dir = '/u/sciteam/zhang7/svn_personal/qe/pp'
 /
 &system
    ibrav = 0
    ntyp = 2,
    nat = 1024,
    occupations = 'fixed',
    ecutwfc = 80
 /
 &electrons
    diagonalization = 'david',
    mixing_beta = 0.7,
    conv_thr = 1.d-8,
    electron_maxstep = 1000,
    startingwfc = 'atomic+random'
 /

K_POINTS automatic
  1 1 1 0 0 0

与前面的输入文件相比，这儿多了这两个选项：tprnfor = .true. 和 tstress = .true.。这样的话才会在电子自洽计算完成后输出每个原子的受力情况，用以构建后来的Force constant matrix.

 

第四步: 超胞的电子自洽计算/Electron SCF calculation of supercells

依次运行pw.x < supercell-???.in即可。

 

第五步： 用Phonopy计算