mirror of https://gitlab.com/QEF/q-e.git
387 lines
8.5 KiB
Bash
Executable File
387 lines
8.5 KiB
Bash
Executable File
#!/bin/sh
|
|
|
|
# run from directory where this script is
|
|
cd `echo $0 | sed 's/\(.*\)\/.*/\1/'` # extract pathname
|
|
EXAMPLE_DIR=`pwd`
|
|
|
|
# check whether echo has the -e option
|
|
if test "`echo -e`" = "-e" ; then ECHO=echo ; else ECHO="echo -e" ; fi
|
|
|
|
$ECHO
|
|
$ECHO "$EXAMPLE_DIR : starting"
|
|
$ECHO
|
|
$ECHO "This example will calculate the water dipole and calculate the work"
|
|
$ECHO "function on a Ni slab with a CO molecule adsorbed using the dipole"
|
|
$ECHO "correction"
|
|
|
|
# set the needed environment variables
|
|
. ../../../environment_variables
|
|
|
|
# required executables and pseudopotentials
|
|
BIN_LIST="pw.x pp.x average.x"
|
|
PSEUDO_LIST="C.pbe-rrkjus.UPF O.pbe-rrkjus.UPF Ni.pbe-nd-rrkjus.UPF HUSPBE.RRKJ3"
|
|
|
|
$ECHO
|
|
$ECHO " executables directory: $BIN_DIR"
|
|
$ECHO " pseudo directory: $PSEUDO_DIR"
|
|
$ECHO " temporary directory: $TMP_DIR"
|
|
$ECHO
|
|
$ECHO " checking that needed directories and files exist...\c"
|
|
|
|
# check for directories
|
|
for DIR in "$BIN_DIR" "$PSEUDO_DIR" ; do
|
|
if test ! -d $DIR ; then
|
|
$ECHO
|
|
$ECHO "ERROR: $DIR not existent or not a directory"
|
|
$ECHO "Aborting"
|
|
exit 1
|
|
fi
|
|
done
|
|
for DIR in "$TMP_DIR" "$EXAMPLE_DIR/results" ; do
|
|
if test ! -d $DIR ; then
|
|
mkdir $DIR
|
|
fi
|
|
done
|
|
cd $EXAMPLE_DIR/results
|
|
|
|
# check for executables
|
|
for FILE in $BIN_LIST ; do
|
|
if test ! -x $BIN_DIR/$FILE ; then
|
|
$ECHO
|
|
$ECHO "ERROR: $BIN_DIR/$FILE not existent or not executable"
|
|
$ECHO "Aborting"
|
|
exit 1
|
|
fi
|
|
done
|
|
|
|
# check for gnuplot
|
|
GP_COMMAND=`which gnuplot 2>/dev/null`
|
|
if [ "$GP_COMMAND" = "" ]; then
|
|
$ECHO
|
|
$ECHO "gnuplot not in PATH"
|
|
$ECHO "Results will not be plotted"
|
|
fi
|
|
|
|
# check for pseudopotentials
|
|
for FILE in $PSEUDO_LIST ; do
|
|
if test ! -r $PSEUDO_DIR/$FILE ; then
|
|
$ECHO
|
|
$ECHO "Downloading $FILE to $PSEUDO_DIR...\c"
|
|
$WGET $PSEUDO_DIR/$FILE $NETWORK_PSEUDO/$FILE 2> /dev/null
|
|
fi
|
|
if test $? != 0; then
|
|
$ECHO
|
|
$ECHO "ERROR: $PSEUDO_DIR/$FILE not existent or not readable"
|
|
$ECHO "Aborting"
|
|
exit 1
|
|
fi
|
|
done
|
|
$ECHO " done"
|
|
|
|
# how to run executables
|
|
PW_COMMAND="$PARA_PREFIX $BIN_DIR/pw.x $PARA_POSTFIX"
|
|
PP_COMMAND="$PARA_PREFIX $BIN_DIR/pp.x $PARA_POSTFIX"
|
|
AVERAGE_COMMAND="$BIN_DIR/average.x"
|
|
$ECHO
|
|
$ECHO " running pw.x as: $PW_COMMAND"
|
|
$ECHO " running pp.x as: $PP_COMMAND"
|
|
$ECHO " running average.x as: $AVERAGE_COMMAND"
|
|
$ECHO " running gnuplot as: $GP_COMMAND"
|
|
$ECHO
|
|
|
|
#
|
|
# Run the calculation for the Ni+CO slab with the dipole correction
|
|
# emaxpos, the starting of the potential inversion, has to be carefully
|
|
# placed in a position where there's almost no charge.
|
|
#
|
|
cat > ni+co.scf.in << EOF
|
|
&control
|
|
calculation='scf',
|
|
restart_mode='from_scratch',
|
|
prefix='ni+co',
|
|
pseudo_dir = '$PSEUDO_DIR',
|
|
outdir='$TMP_DIR'
|
|
|
|
tefield=.true.,
|
|
dipfield=.true.,
|
|
/
|
|
&system
|
|
|
|
nat=5, ntyp=3,
|
|
ibrav=0, celldm(1)=4.70366666,
|
|
|
|
ecutwfc = 30.0
|
|
occupations='smearing', smearing='m-v', degauss=0.03
|
|
|
|
edir=3,
|
|
emaxpos=0.55,
|
|
eopreg=0.06,
|
|
eamp=0,
|
|
|
|
/
|
|
&electrons
|
|
mixing_beta = 0.3
|
|
conv_thr = 1.0d-6
|
|
/
|
|
|
|
ATOMIC_SPECIES
|
|
C 1.0 C.pbe-rrkjus.UPF
|
|
O 1.0 O.pbe-rrkjus.UPF
|
|
Ni 1.0 Ni.pbe-nd-rrkjus.UPF
|
|
|
|
CELL_PARAMETERS alat
|
|
1.00000000 0.00000000 0.00000000
|
|
0.00000000 1.41421356 0.00000000
|
|
0.00000000 0.00000000 9.10000001
|
|
|
|
ATOMIC_POSITIONS (alat)
|
|
|
|
C -0.00364039 0.02119538 1.54673745
|
|
O -0.00634860 0.04192428 2.02021975
|
|
Ni 0.48527378 0.00197332 0.97713547
|
|
Ni -0.00049546 0.70236680 0.45417840
|
|
Ni 0.50000000 0.00000000 0.00000000
|
|
|
|
K_POINTS {gamma}
|
|
EOF
|
|
$ECHO
|
|
$ECHO " running pw.x for Ni+CO slab...\c"
|
|
$PW_COMMAND < ni+co.scf.in > ni+co.scf.out
|
|
check_failure $?
|
|
$ECHO " done"
|
|
|
|
#
|
|
# Extract the potential from the SCF calculation
|
|
#
|
|
cat > ni+co.pp.in << EOF
|
|
&inputpp
|
|
prefix='ni+co',
|
|
outdir='$TMP_DIR',
|
|
filplot = 'ni+co.vpot'
|
|
plot_num= 11
|
|
/
|
|
EOF
|
|
$ECHO
|
|
$ECHO " running pp.x to extract ni+co potential...\c"
|
|
$PP_COMMAND < ni+co.pp.in > ni+co.pp.out
|
|
check_failure $?
|
|
$ECHO " done"
|
|
|
|
#
|
|
# Average the potential to obtain the planar average along edir
|
|
#
|
|
cat > ni+co.avg.in << EOF
|
|
1
|
|
ni+co.vpot
|
|
1.D0
|
|
150
|
|
3
|
|
3.000000
|
|
EOF
|
|
$ECHO
|
|
$ECHO " running average.x to obtain the potential along Z...\c"
|
|
$AVERAGE_COMMAND < ni+co.avg.in > ni+co.vpot-z
|
|
check_failure $?
|
|
$ECHO " done"
|
|
mv avg.dat ni+co.vpot-z
|
|
|
|
#
|
|
# If gnuplot is present plot the potential. In the graph it's clear where
|
|
# the potential reaches a constant value thanks to the dipole correction.
|
|
#
|
|
|
|
if [ "$GP_COMMAND" = "" ]; then
|
|
break
|
|
else
|
|
cat > ni+co.gnuplot <<EOF
|
|
set terminal postscript enhanced color
|
|
set output 'ni+co.eps'
|
|
|
|
set multiplot
|
|
|
|
set size 1,1
|
|
set origin 0,0
|
|
|
|
set style line 1 lt 1
|
|
set style line 2 lt 2
|
|
set style line 3 lt 3
|
|
|
|
set arrow 1 from 8,-2.2976 to 8,4.75 ls 3
|
|
set arrow 2 from 8,4.75 to 8,-2.2976 ls 3
|
|
#set arrow 3 from 18,4.488 to 16.5,-4 nohead ls 3
|
|
|
|
#set object 1 rect from 16.5,4.488 to 19.05048,4.896 fc ls 3 fs empty
|
|
|
|
set xlabel "Length (Angstroms)"
|
|
set ylabel "Energy (eV)"
|
|
|
|
set label 'W_{funct} = 7.04 eV' at 8.2,0
|
|
|
|
plot '< tail -150 ni+co.vpot-z' u (\$1*0.52918):(\$2*13.6) w l ls 1 title 'Pot.' , '< tail -150 ni+co.vpot-z' u (\$1*0.52918):(-2.2976) w l ls 2 title 'E_{f}'
|
|
|
|
|
|
|
|
set size 0.4,0.4
|
|
set origin 0.45,0.1
|
|
set xrange [8:18]
|
|
set yrange[4.74:4.76]
|
|
|
|
set xtics 8.2,2.8,18
|
|
set ytics 4.742, 0.005, 4.758
|
|
|
|
set style line 1 lt 1
|
|
set style line 2 lt 2
|
|
set style line 3 lt 3
|
|
|
|
set arrow 1 from 10,4.7566 to 10,4.7529 ls 3
|
|
set arrow 2 from 10,4.7529 to 10,4.7566 ls 3
|
|
#unset arrow 3
|
|
|
|
set label "Dipole corr." at 10.8,4.748
|
|
|
|
unset xlabel
|
|
unset ylabel
|
|
|
|
plot [8:18][4.74:4.76] "< tail -150 ni+co.vpot-z" u (\$1*0.52918):(\$3*13.6) w l ls 1 notitle, "< tail -150 ni+co.vpot-z" u (\$1*0.52918):(4.7566) w l ls 2 notitle, "< tail -150 ni+co.vpot-z" u (\$1*0.52918):(4.7529) w l ls 2 notitle
|
|
|
|
unset multiplot
|
|
set output
|
|
reset
|
|
EOF
|
|
$ECHO " Plotting averages and work function results ...\c"
|
|
$GP_COMMAND < ni+co.gnuplot
|
|
$ECHO " done"
|
|
rm ni+co.gnuplot
|
|
fi
|
|
|
|
|
|
###############################
|
|
# WATER
|
|
###############################
|
|
|
|
|
|
#
|
|
# Run the calculation for the water molecule with the dipole correction
|
|
#
|
|
cat > water.scf.in << EOF
|
|
&control
|
|
calculation='scf',
|
|
restart_mode='from_scratch',
|
|
prefix='water',
|
|
pseudo_dir = '$PSEUDO_DIR',
|
|
outdir='$TMP_DIR'
|
|
|
|
tefield=.true.,
|
|
dipfield=.true.
|
|
/
|
|
&SYSTEM
|
|
ibrav=1, celldm(1) = 15
|
|
nat=3, ntyp=2
|
|
|
|
ecutwfc=30.0
|
|
occupations='smearing', degauss=0.01
|
|
|
|
edir=3
|
|
eamp=0.D0
|
|
eopreg=0.1
|
|
emaxpos=0.6
|
|
/
|
|
&ELECTRONS
|
|
mixing_beta = 0.7
|
|
conv_thr = 1.0d-8
|
|
/
|
|
ATOMIC_SPECIES
|
|
O 15.9994 O.pbe-rrkjus.UPF
|
|
H 1.00794 HUSPBE.RRKJ3
|
|
ATOMIC_POSITIONS { Angstrom }
|
|
O 0.0 0.0 0.0
|
|
H 0.77 0.0 0.62
|
|
H -0.77 0.0 0.62
|
|
K_POINTS { gamma }
|
|
EOF
|
|
$ECHO
|
|
$ECHO " running pw.x for water molecule ...\c"
|
|
$PW_COMMAND < water.scf.in > water.scf.out
|
|
check_failure $?
|
|
$ECHO " done"
|
|
|
|
#
|
|
# Extract the potential from the SCF calculation
|
|
#
|
|
cat > water.pp.in << EOF
|
|
&inputpp
|
|
prefix='water',
|
|
outdir='$TMP_DIR',
|
|
filplot = 'water.vpot'
|
|
plot_num= 11
|
|
/
|
|
EOF
|
|
$ECHO
|
|
$ECHO " running pp.x to extract water potential...\c"
|
|
$PP_COMMAND < water.pp.in > water.pp.out
|
|
check_failure $?
|
|
$ECHO " done"
|
|
|
|
#
|
|
# If gnuplot is present plot the potential. In the graph it's clear how the
|
|
# dipole influences the potential slope
|
|
#
|
|
cat > water.avg.in << EOF
|
|
1
|
|
water.vpot
|
|
1.D0
|
|
90
|
|
3
|
|
3.000000
|
|
EOF
|
|
$ECHO
|
|
$ECHO " running average.x to obtain the potential along Z...\c"
|
|
$AVERAGE_COMMAND < water.avg.in > water.vpot-z
|
|
check_failure $?
|
|
$ECHO " done"
|
|
mv avg.dat water.vpot-z
|
|
|
|
|
|
if [ "$GP_COMMAND" = "" ]; then
|
|
break
|
|
else
|
|
cat > water.gnuplot <<EOF
|
|
set terminal postscript enhanced color
|
|
set output 'water.eps'
|
|
|
|
set style line 1 lt 1
|
|
set style line 2 lt 2
|
|
set style line 3 lt 3
|
|
|
|
set arrow 1 from 2.5,-0.15 to 2.5,0.93 ls 3
|
|
set arrow 2 from 2.5,0.93 to 2.5,-0.15 ls 3
|
|
|
|
set label "Dipole correction" at 3,0.5
|
|
|
|
set yrange [-1.2:1.2]
|
|
|
|
set xlabel "Length (Angstroms)"
|
|
set ylabel "Energy (eV)"
|
|
|
|
plot "< tail -90 water.vpot-z" u (\$1*0.52918):(\$2*13.6) w l ls 1 title "Pot.", "< tail -90 water.vpot-z" u (\$1*0.52918):(-0.15) w l ls 2 notitle, "< tail -90 water.vpot-z" u (\$1*0.52918):(0.93) w l ls 2 notitle
|
|
|
|
set output
|
|
reset
|
|
|
|
EOF
|
|
$ECHO " Plotting averages and work function results ...\c"
|
|
$GP_COMMAND < water.gnuplot
|
|
$ECHO " done"
|
|
rm water.gnuplot
|
|
fi
|
|
|
|
$ECHO " cleaning $TMP_DIR...\c"
|
|
rm -rf $TMP_DIR/ni+co.* $TMP_DIR/water.*
|
|
|
|
$ECHO "Calculations done."
|
|
$ECHO ""
|
|
$ECHO "----------------------------------------------------"
|
|
$ECHO " results/ni+co.eps -> Graph for the slab calculation"
|
|
$ECHO " results/water.eps -> Graph for the water molecule"
|
|
$ECHO "----------------------------------------------------"
|