How to run with the same pseudos in S and A

How to run with the same pseudos
in SIESTA and ABINIT
-215.78
Siesta (DZP)
Abinit (13.005 Ha)
Free Energy (eV)
-215.79
-215.80
-215.81
-215.82
5.30
5.35
5.40
5.45
Lattice constant (Ang)
5.50
Objectives
Run examples with the same pseudos (same decomposition in
local part and Kleinman-Bylander projectors) in SIESTA and ABINIT.
Compare total energies
Download the last versions of both codes,
SIESTA and ABINIT
Regarding the ABINIT code, you can download the required version from:
http://personales.unican.es/junqueraj/Abinit.tar.gz
But the merge of the relevant subroutines into the main trunk
will be done soon
Regarding SIESTA, the code is available at the usual web site:
http://www.icmab.es/siesta
A few modifications to be done before running:
SIESTA
Edit the file atom.F in the Src directory
I.
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC FIELD.
1. Replace “nrval”1by “nrwf” in the call to schro_eq inside the subroutine KBgen
(file atom.F) CiteI.. ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC
I. ENERGY FUNCTIONAL
FOR A DIELECTRIC INSIDE AN ELECTRIC FIEL
if(erefkb(ikb,l).ge.1.0d3)
then
Cite1 .
nnodes=ikb
Cite1 .
nprin=l+1
call schro_eq(Zval,rofi,vps(1,l),ve,s,drdi,
if(erefkb(ikb,l).ge.1.0d3)
then
if(erefkb(ikb,l).ge.1.0d3)
then
.
nrwf,l,a,b,nnodes,nprin,
nnodes=ikb
nnodes=ikb erefkb(ikb,l),rphi(1,ikb))
.
nprin=l+1
nprin=l+1
call schro_eq(Zval,rofi,vps(1,l),ve,s,drdi,
call
schro_eq(Zval,rofi,vps(1,l),ve,s,drdi,
if(nkbl(l).eq.1)
then
.
nrwf,l,a,b,nnodes,nprin,
.
nrwf,l,a,b,nnodes,nprin,
call in
ghost(Zval,rofi,vps(:,l),vlocal,
2. Replace “nrval”
by. “nrwf”
the call
toerefkb(ikb,l),rphi(1,ikb))
ghost inside the subroutine KBgen (file
.
ve,s,drdi,nrwf,l,a,b,nrwf,
.
erefkb(ikb,l),rphi(1,ikb))
atom.F)
.
erefkb(ikb,l),rphi(:,ikb),ighost)
if(nkbl(l).eq.1) then
else
if(nkbl(l).eq.1)
then
call ghost(Zval,rofi,vps(:,l),vlocal,
.call ghost(Zval,rofi,vps(:,l),vlocal,
ve,s,drdi,nrwf,l,a,b,nrwf,
real(dp),
parameter :: Rmax_kb_default = 60.0d0
.
erefkb(ikb,l),rphi(:,ikb),ighost)
.
ve,s,drdi,nrwf,l,a,b,nrwf,
else
.
erefkb(ikb,l),rphi(:,ikb),ighost)
else
parameter :: to
Rmax_kb_default
= 60.0d0
3. Increase the default
ofreal(dp),
Rmax_kb_default
60.0 bohrs (file
atom.F)
1
2
3
::J. Rmax_kb_default
= 60.0d0
G. real(dp),
Lee, Y.-H. Shin,parameter
and J. Y. Son,
Am. Ceram. Soc., 95, 2773
(2012).
E. Wiesendanger, Ferroelectrics, 6, 263 (1974).
H. P. Sun, W. Tian, X. Q. Pan, J. H. Haeni, and D. G. Schlom, Appl. Phys. Lett., 84, 3298 (2004).
Definition of the Kleinman-Bylander projectors
X. Gonze et al., Phys. Rev. B 44, 8503 (1991)
The normalized Kleinman-Bylander projectors are given by
where
and
Note that these
are the eigenstates of the semilocal pseudopotential
(screened by the pseudovalence charge density).
are the radial part of the wave function multiplied by
,
For the sake of simplicity, we assume here only one projector per angular
momentum shell. If more than one is used, they must be orthogonalized
Definition of the Kleinman-Bylander projectors
(old choice in SIESTA for the atomic eigenstates)
In the standard version of SIESTA, the
Schrödinger equation for the isolated
atom while generating the KB
projectors is solved inside a box
whose size is determined by nrval.
This is usually a very large radius (of
the order of 120 bohrs)
s - shell
p - shell
d - shell
f - shell
ul
PS
(r)
1,0
0,5
0,0
0
10
Rmax_KB
20
30
radius (bohr)
40
50
Then, this wave functions is
normalized inside a sphere of a much
smaller radius, determined by
Rmax_KB (default value = 6.0 bohr)
Definition of the Kleinman-Bylander projectors
(new choice in SIESTA for the atomic eigenstates)
0,8
s - shell
p - shell
d - shell
f - shell
0,6
0,4
ul
PS
(r)
In the new version of SIESTA, everything
is consistent, and the Schrödinger
equation and the normalization are
solved with respect the same boundary
conditions
0,2
0,0
0
10
20
30
40
radius (bohr)
50
60
Rmax_KB
Almost no change in total energies observed, but the Kleinman-Bylander
energies might be very different, specially for unbounded orbitals
A few modifications to be done before running:
ABINIT
Edit the file src/65_psp/psp5nl.F90
1. Uncomment the last two lines
for the sake
of comparing
I. ENERGY
FUNCTIONAL
FOR Kleinman-Bylander
A DIELECTRIC INSIDE
energies and cosines with the ones obtained with SIESTA
Cite1 .
!DEBUG
write(std_out,*) ’EKB=’,(ekb(iq),iq=1,4)
write(std_out,*) ’COSKB=’,(ckb(iq),iq=1,4)
!ENDDEBUG
if(nkbl(l).eq.1) then
call
ghost(Zval,rofi,vps(:,l),vlocal,
2. Remember to compile the code
enabling
the FOX library
.
ve,s,drdi,nrwf,l,a,b,nrwf,
.
erefkb(ikb,l),rphi(:,ikb),ighost)
else
./configure --with-trio-flavor=netcdf+etsf_io+fox
real(dp), parameter :: Rmax_kb_default = 60.0d0
1
AN ELE
Examples to run SIESTA and ABINIT
with the same pseudos
1. Visit the web page:
http://personales.unican.es/junqueraj
and follow the links:
Teaching
Métodos Computacionales en Estructura de la Materia
Hand-on sessions
Pseudos
2. Download the Pseudos and input files for both codes
3. Untar the ball file
$ tar –xvf Siesta-Abinit.tar
This will generate a directory called Comparison-Siesta-Abinit with 4 directories:
$ cd Comparison-Siesta-Abinit
$ ls -ltr
$
Si
(example for a covalent semiconductor, LDA)
$
Al
(example for a sp-metal, LDA)
$
Au
(example for a noble metal, includes d-orbital, LDA)
$
Fe
(example for a transition metal, includes NLCC, GGA)
Examples to run SIESTA and ABINIT
with the same pseudos
In every subdirectory it can be found:
$ cd Si
$ ls –ltr
$
Pseudo
(files to generate and test the pseudopotential)
$
Optimized-Basis
(files to optimize the basis set)
$
Runsiesta
(files to run Siesta)
$
Runabinit
(files to run Abinit)
How to generate and test
a norm-conserving pseudopotential
Generate the pseudopotential using the ATM code as usual,
following the notes in the Tutorial
“How to generate a norm conserving pseudopotential”
Copy the input file
in the corresponding
atom/Tutorial/PS_Generation
directory
I. ENERGY
FUNCTIONAL FOR
A DIELECTRIC INSIDE AN ELECTRIC
FIELD.
and run
Cite1 .
$ ../../Utils/pg.sh Si.tm2.inp
==> Output data in directory Si.tm2
==> Pseudopotential in Si.tm2.vps and Si.tm2.psf (and maybe in Si.tm2.xml)
if(nkbl(l).eq.1) then
call ghost(Zval,rofi,vps(:,l),vlocal,
The pseudopotentials
will be on the same parent directory:
.
ve,s,drdi,nrwf,l,a,b,nrwf,
.vps
(required to test the pseudopotential)
. (unformatted)
erefkb(ikb,l),rphi(:,ikb),ighost)
else
.psf (formatted)
.xml (inreal(dp),
XML format)
to run Abinit)
parameter(required
:: Rmax_kb_default
= 60.0d0
Remember to test the pseudopotential using the ATM code as usual,
following the notes in the Tutorial
“How
to test the transferability of a norm conserving pseudopotential”
1
G. Lee, Y.-H. Shin, and J. Y. Son, J. Am. Ceram. Soc., 95, 2773 (2012).
2
3
E. Wiesendanger, Ferroelectrics, 6, 263 (1974).
H. P. Sun, W. Tian, X. Q. Pan, J. H. Haeni, and D. G. Schlom, Appl. Phys. Lett., 84, 3298 (2004).
Running the energy versus lattice
constant curve in SIESTA
Run the energy versus lattice constant curve in SIESTA as usual.
You can use both the .psf or the .xml pseudopotential.
Follow the rules given in the tutorial
“Lattice constant, bulk modulus, and equilibrium energy of solids”
The input file has been prepared for you (file Si.fdf).
Since we are interested in compare the performance of the basis set, it is
important to converge all the rest of approximations (Mesh Cutoff, k-point grid,
etc.) as much as possible
2
I. ENERGY
FUNCTIONAL
FOR to
A DIELECTRIC
INSIDE
AN called
ELECTRIC
FIELD.
At the end,
we would
be able
write a file
(here
Si.siesta.latcon.dat)
that looks like this:
Cite1 .
# Lattice Constant (Ang)
5.30
5.32
5.34
5.36
5.38
5.40
5.42
5.44
5.46
5.48
Free Energy (eV)
-215.790381
-215.802337
-215.811082
-215.816643
-215.819095
-215.818374
-215.814586
-215.808022
-215.798710
-215.786727
KBgen: Kleinman-Bylander projectors:
l= 0
rc= 1.936440
el= -0.796617
Ekb=
These data have been obtained with a
double-zeta plus polarization basis
set, optimized at the theoretical lattice
constant with a pressure of 0.05 GPa
4.661340
kbcos=
0.299549
Running the energy versus cutoff energy in ABINIT
To check the equivalent cutoff energy in Abinit:
1.  We run the same system (same lattice vectors and internal coordinates) at the
same level of approximations (same exchange and correlation functional,
Monkhorst-Pack mesh etc.) at a given lattice constant.
Here it has been written for you (file Si.input.convergence)
3
#Number of atoms, chemical species and atom types
natom
2
# Number of atoms in the unit cell
ntypat 1
# Number of types of atoms
typat
1 1
# Type of atoms
znucl
14.0
# Gives nuclear charge for each type of
#Coordinates and cell variables
rprim
0.0
0.5
0.5
0.5
0.0
0.5
0.5
0.5
0.0
acell
5.38
5.38
5.38
Angstrom
xred -0.125 -0.125 -0.125
0.125 0.125 0.125
#PlaneWave cutoff and k-grid mesh integration
ecutsm 0.5
# Energy cutoff smearing (Ha)
nband
10
# Number of bands
kptopt 1
# Kpoints option
# 0 = read directly nkpt, kpt, kptnrm and wtk
# 1 = rely on ngkpt or kptrlatt, as well as
#
on nshiftk and shiftk to set up
#
the k points.
#
Full symmetry taken into account.
# 2 = 1, but only time reversal symmetry
#
is taken into account.
# 3 = 1, but do not take into account
#
any symmetry
# A negative value = rely on kptbounds,
#
and ndivk to set up a band structure
#
calculation along different lines
ngkpt
3 3 3
# Number of grid points for k points
#
generation
# This is a 3x3x3 FCC grid,
#
based on the primitive vectors
#
of the reciprocal space.
#
For a FCC real space lattice,
#
like the present one,
#
it actually corresponds to the
#
so-called 6x6x6 Monkhorst-Pack grid,
#
if the following shifts are used :
nshiftk 4
#
shiftk 0.5 0.5 0.5
# Shift for k points
0.5 0.0 0.0
0.0 0.5 0.0
0.0 0.0 0.5
#Exchange-correlation
ixc
2
# Integer for exchange-correlation choice
# 0 = No xc
# 1 = LDA or LSD, Teter Pade parametrization
# 2 = LDA, Perdew-Zunger-Ceperley-Alder
# 3 = LDA, old Teter rational polynomial
#
parametrization, fit to Ceperley-Alder
#
data (no spin-polarization: no sp)
# 4 = LDA, Wigner functional (no sp)
# 5 = LDA, Hedin-Lundqvist functional (no sp)
# 6 = LDA, "X-alpha" functional (no sp)
# 7 = LDA or LSD, Perdew-Wang 92 functional
# 8 = LDA or LSD, x-only part of the PW 92
# 9 = LDA or LSD, x- and RPA part of the PW92
# 10= GGA, Perdew-Burke-Ernzerhof
Diamond structure at the lattice
constant that minimizes the energy
in SIESTA
6 × 6 × 6 Monkhorst-Pack mesh
Ceperley-Alder (LDA) functional
Running the energy versus cutoff energy in ABINIT
To check the equivalent cutoff energy in Abinit:
1.  We run the same system (same lattice vectors and internal coordinates) at the
same level of approximations (same exchange and correlation functional,
Monkhorst-Pack mesh etc.) at a given lattice constant.
Here it has been written for you (file Si.input.convergence)
ndtset
2. Change the cutoff energy for the plane waves
$ more Si.files
$ more Si.files Si.input.convergence
Si.input.convergence
Si.out
3. Edit the .files file and
select
Sii
Si.out
Sio
the input file and theSii
pseudo
t1x
file (in XML format) Sio
Si.xml
t1x
Si.xml
$ abinit < Si.files > Si.log
4. Run the code
$ abinit < Si.files > Si.log &
1
2
3
ecut1
ecut2
ecut3
ecut4
ecut5
ecut6
ecut7
ecut8
ecut9
ecut10
ecut11
ecut12
ecut13
ecut14
&ecut15
ecut16
ecut17
ecut18
ecut19
ecut20
ecut21
21
4.00
5.00
6.00
7.00
8.00
9.00
10.00
11.00
12.00
13.00
14.00
15.00
16.00
17.00
18.00
19.00
20.00
25.00
30.00
35.00
40.00
G. Lee, Y.-H. Shin, and J. Y. Son, J. Am. Ceram. Soc., 95, 2773 (2012).
E. Wiesendanger, Ferroelectrics, 6, 263 (1974).
H. P. Sun, W. Tian, X. Q. Pan, J. H. Haeni, and D. G. Schlom, Appl. Phys. Lett., 84, 329
Running the energy versus cutoff energy in ABINIT:
bulk Si (covalent semiconductor)
Write the total energy as a function of the cutoff energy and edit the
corresponding file that should look like this
3
$ grep Total Si.out > Si.abinit.convergence.dat
$ more Si.abinit.convergence.dat
# Cutoff energy (Ha)
Total energy (eV)
4.0
-2.14194449515052E+02
5.0
-2.14976484971340E+02
6.0
-2.15411139789484E+02
7.0
-2.15628684537067E+02
8.0
-2.15720857649398E+02
9.0
-2.15750170913495E+02
10.0
-2.15763328854496E+02
11.0
-2.15779811817589E+02
12.0
-2.15799418159205E+02
13.0
-2.15818945848621E+02
14.0
-2.15839165600763E+02
15.0
-2.15857713873448E+02
16.0
-2.15872186183511E+02
17.0
-2.15882236230850E+02
18.0
-2.15889235327402E+02
19.0
-2.15893738919208E+02
20.0
-2.15896300628071E+02
25.0
-2.15898209214678E+02
30.0
-2.15900079984311E+02
35.0
-2.15901734504665E+02
40.0
-2.15902017447268E+02
$ gnuplot
gnuplot> plot "Si.abinit.convergence.dat" u 1:2 w l
gnuplot> set terminal postscript
Terminal type set to ’postscript’
Options are ’landscape noenhanced defaultplex \
leveldefault monochrome colortext \
dashed dashlength 1.0 linewidth 1.0 butt \
palfuncparam 2000,0.003 \
"Helvetica" 14 ’
gnuplot> set output "Si.abinit.convergence.ps"
gnuplot> replot
-214
"Si.abinit.convergence.dat" u 1:2
-214.2
-214.4
-214.6
Equivalent PW
cutoff a DZP basis
set at 5.38 Å
-214.8
-215
-215.2
-215.4
-215.6
-215.8
-216
0
5
10
15
20
25
30
35
40
3
$ $ more Si.abinit.latcon.dat
# Lattice constant (Ang)
Total Energy (eV) (Cutoff = 13.005 Ha)
$ $ more5.30
Si.abinit.latcon.dat-2.15791847116547E+02
# Lattice
constant (Ang)
Total Energy (eV) (Cutoff = 13.005 Ha)
5.32
-2.15803559248711E+02
5.30
-2.15791847116547E+02
5.34
-2.15811826556959E+02
5.32
-2.15803559248711E+02
5.36
-2.15816906510895E+02
5.38
-2.15819045481337E+02
5.34
-2.15811826556959E+02
5.40
-2.15818274869825E+02
5.36
-2.15816906510895E+02
5.42
-2.15814668146456E+02
5.38
-2.15819045481337E+02
5.44
-2.15808351816212E+02
5.40
-2.15818274869825E+02
5.46
-2.15799487297689E+02
5.42
-2.15814668146456E+02
5.48
-2.15788012963537E+02
Running the energy versus lattice
constant curve in ABINIT
1. Same input as before but…
5.44
-2.15808351816212E+02
5.46
-2.15799487297689E+02
#input for
bulk Si in the diamond structure.
ecut
5.48
-2.15788012963537E+02
13.005
# Energy cutoff (Ha)
#input for bulk Si in the diamond structure.
ndtset
10
ecut
13.005
acell1
acell2
ndtset
acell3
acell4
acell1
acell5
acell6
acell2
acell7
acell3
acell8
acell4
acell9
acell5
acell10
5.30
5.32
10
5.34
5.36
5.30
5.38
5.40
5.32
5.42
5.34
5.44
5.36
5.46
5.38
5.48
acell6
5.40
acell7
5.42
$ more Si.files
acell8
5.44
Si.input.latcon
Si.out
acell9
5.46
Sii
acell10 5.48
5.30
5.32
5.34
5.36
5.30
5.38
5.40
5.32
5.42
5.34
5.44
5.36
5.46
5.38
5.48
5.40
5.42
5.44
5.46
5.48
… setting the plane wave cutoff to the
equivalent one to a DZP basis set
# Energy cutoff (Ha)
5.30
Angstrom
5.32
Angstrom
5.34
Angstrom
5.36
Angstrom
5.30 Angstrom
Angstrom
5.38
5.40
5.32 Angstrom
Angstrom
5.42
5.34 Angstrom
Angstrom
5.44
5.36 Angstrom
Angstrom
5.46
Angstrom
5.38
Angstrom
5.48
Angstrom
5.40
5.42
5.44
5.46
5.48
… and changing the lattice constant
embracing the minimum
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
2. Change in the .files the name of the input file
$ more Si.files
Si.input.convergence
Si.out
Si.out
$ gnuplot
Sii
Sii
gnuplot> plot "Si.abinit.convergence.dat" u 1:2 w l
Sio
gnuplot> set terminal postscriptSio
t1x
Terminal type set to ’postscript’
t1x
Options are ’landscape noenhanced defaultplex \
Si.xml
Si.xml
leveldefault monochrome colortext
\
Sio
$t1x
more Si.files
Si.xml
Si.input.latcon
dashed dashlength 1.0 linewidth 1.0 butt \
$ gnuplot
palfuncparam
2000,0.003 \
gnuplot>
plot "Si.abinit.convergence.dat"
1:2 w l
$ abinit < uSi.files
"Helvetica" 14 ’
gnuplot> set terminal postscript
gnuplot> set output "Si.abinit.convergence.ps"
Terminal
type set to ’postscript’
gnuplot> replot
Options are ’landscape noenhanced defaultplex \
leveldefault monochrome colortext \
3. Run the code
> Si.log &
Running the energy versus lattice
constant curve in ABINIT
Write the total energy as a function of the lattice constant and
edit the corresponding file that should look like this
3
$ grep "Total" Si.out > Si.abinit.latcon.dat
$ more Si.abinit.latcon.dat
# Lattice constant (Ang)
Total Energy (eV) (Cutoff = 13.005 Ha)
5.30
-2.15791847116547E+02
5.32
-2.15803559248711E+02
5.34
-2.15811826556959E+02
5.36
-2.15816906510895E+02
5.38
-2.15819045481337E+02
5.40
-2.15818274869825E+02
5.42
-2.15814668146456E+02
5.44
-2.15808351816212E+02
5.46
-2.15799487297689E+02
5.48
-2.15788012963537E+02
#input for bulk Si in the diamond structure.
ecut
13.005
ndtset
10
acell1
acell2
acell3
acell4
acell5
acell6
acell7
acell8
acell9
acell10
5.30
5.32
5.34
5.36
5.38
5.40
5.42
5.44
5.46
5.48
# Energy cutoff (Ha)
5.30
5.32
5.34
5.36
5.38
5.40
5.42
5.44
5.46
5.48
5.30
5.32
5.34
5.36
5.38
5.40
5.42
5.44
5.46
5.48
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Comparing the pseudopotential in SIESTA and
ABINIT: bulk Si (covalent semiconductor)
To be totally sure that we have run SIESTA and ABINIT with the same
peudopotential operator, i.e. with the same decomposition in local part and
Kleinman-Bylander projectors:
I.
2
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC FIELD.
1. Edit one of the
output files in SIESTA and search for the following lines:
1
Cite .
KBgen: Kleinman-Bylander
l= 0
rc= 1.936440
l= 1
rc= 1.936440
l= 2
rc= 1.936440
l= 3
rc= 1.936440
projectors:
el= -0.796617
4.661340
kbcos=
0.299549
I. ENERGYEkb=
FUNCTIONAL
FOR A
DIELECTRIC
INSIDE AN ELECTRIC FIELD.
el= -0.307040
Ekb= 1.494238
kbcos= 0.298299
Cite1 . 0.008863
el=
Ekb= -2.809034
kbcos= -0.000566
el= 0.013111
Ekb= -0.959387
kbcos= -0.000004
KBgen: Kleinman-Bylander projectors:
l= 0
2. Edit the
rc=
1.936440
el= -0.796617
Ekb=
4.661340
kbcos=
0.299549
2
EKB=
2.3306701075631429 l= 1 0.74711916272986667
-1.4045230577882524
rc= 1.936440
el= -0.307040
Ekb=
1.494238
kbcos= 0.298299-0.4796986487264058
l=
2
rc=
1.936440
el=
0.008863
Ekb=
-2.809034
kbcos= -0.000566 -3.9422142725424
COSKB=
0.30046489563475459
0.29921145875238120
-5.9527141624513349E-004
log file in ABINIT and l=search
for theel=following
lines:
3
rc=
1.936440
0.013111
Ekb=
-0.959387
kbcos= -0.000004
pspatm: epsatm=
9.82921496
--- l ekb(1:nproj)
EKB=--> 2.3306701075631429
The Kleinman-Bylander
0
2.330670
0.74711916272986667
-1.4045230577882524
1
0.747119
energies and cosines should
2
-1.404523COSKB= -0.47969864872640589
0.30046489563475459
be the same upto numerical
3
-0.479699
0.29921145875238120
-5.9527141624513349E-004
-3.9422142725424715E-006
pspatm: epsatm=
9.82921496
--- l ekb(1:nproj) -->
0
2.330670
1
0.747119
2
-1.404523
3
-0.479699
roundoff errors
Note: In SIESTA they are
written in Ry and in
ABINIT they are in Ha.
Comparing the energy versus lattice constant in
SIESTA and ABINIT: bulk Si (covalent semiconductor)
3
$ gnuplot
gnuplot> plot "Si.abinit.latcon.dat" u 1:2 w l, "../Runsiesta/Si.siesta.latcon.dat" u 1:2 w l
gnuplot> set terminal postscript colorTerminal type set to ’postscript’
Options are ’landscape noenhanced defaultplex \
leveldefault color colortext \
dashed dashlength 1.0 linewidth 1.0 butt \
palfuncparam 2000,0.003 \
"Helvetica" 14 ’
gnuplot> set output "Si.compar.latcon.ps"gnuplot> replot
#input for bulk Si in the diamond structure.
-215.78
13.005
ndtset
10
acell1
acell2
acell3
acell4
acell5
acell6
acell7
acell8
acell9
acell10
5.30
5.32
5.34
5.36
5.38
5.40
5.42
5.44
5.46
5.48
$ more Si.files
Si.input.latcon
Si.out
Sii
Sio
t1x
Si.xml
# Energy cutoff (Ha)
-215.79
5.30
5.30
5.32
5.32
5.34
5.34
5.36
5.36
5.38 -215.80
5.38
5.40
5.40
5.42
5.42
5.44
5.44
5.46
5.46
5.48 -215.81
5.48
Free Energy (eV)
ecut
-215.82
Siesta (DZP)
Abinit (13.005 Ha)
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
Angstrom
5.30
5.35
5.40
5.45
Lattice constant (Ang)
5.50
Comparing the pseudopotential in
SIESTA and ABINIT: bulk Al (sp metal)
To be totally sure that we have run SIESTA and ABINIT with the same
peudopotential operator, i.e. with the same decomposition in local part and
Kleinman-Bylander projectors:
I.
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC
Cite1 .
# Ekb (in Ha):
#
SIESTA
#
l= 0
1.556933
#
l= 1
0.439408
#
l= 2 -0.891462
#
l= 3 -0.308287
#
# Coskb
#
SIESTA
#
l= 0
0.285004
#
l= 1
0.270194
#
l= 2 -0.001042
#
l= 3 -0.000009
ABINIT
1.556933
0.439408
-0.891467
-0.308292
ABINIT
0.285875
0.271020
-0.001096
-0.000009
ElectronicTemperature 0.02 Ry
occopt
3
tsmear
0.01
# Occupation Option
# Fermi-Dirac smearing (finite-temperature metal)
# Temperature of smearing (in Ha)
Comparing
pseudopotential
in
I. the
ENERGY
FUNCTIONAL FOR A DIELECTRIC
INSIDE AN ELECTRIC FIELD.
SIESTA
and ABINIT: bulk Al (sp metal)
Cite .
I.
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRI
Cite1 .
1
# Ekb (in Ha):
#
SIESTA
ABINIT
#
l= 0
1.556933
1.556933
Ekb (in Ha):
#
l= 1
0.439408
0.439408
SIESTA
ABINIT # l= 2 -0.891462
-0.891467
-0.308292
l= 0
1.556933
1.556933 # l= 3 -0.308287
#
l= 1
0.439408
0.439408 # Coskb
I. ENERGY FUNCTIONAL
SIESTA
ABINIT
l= 2 -0.891462
-0.891467 #
#
l= 0
0.285004
0.285875
l= 3 -0.308287
-0.308292 # l= 1 0.270194
IESTA
BINIT
0.271020
Cite1 .
#
l= 2 -0.001042
-0.001096
#
l=
3
-0.000009
Coskb
I. ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE-0.000009
AN ELECTRIC FIELD.
For #the case of metallic system, besides the k-point sampling we have to pay
particular attention to the occupation option
#
#
#
#
#
S
#
#
#
SIESTA
ABINIT
1
# . l= 0
0.285004
0.285875
Cite
Default:
Fermi-Dirac
#
l= 1
0.270194
0.271020
# Ekb
(in
Ha):
#
l= 2 -0.001042
-0.001096
# # SIESTA
ABINIT
l= 3 -0.000009
-0.000009
FOR A DIELE
A
2
# Ekb (in Ha):
#
SIESTA
ABINIT
occopt #3
Option
l= 0# Occupation
1.432475
1.432472
# Fermi-Dirac smearing (finite-temperature metal)
#
l= 1# Temperature
0.827346
0.827336
tsmear 0.01
of smearing (in Ha)
#
l= 2 -3.174442
-3.174442
#
l= 3 -1.038935
-1.039329
#
l= 0
1.432475
1.432472
# ElectronicTemperature
l= 1
0.827346
0.827336
#
0.02 Ry
#
l= 2 -3.174442
-3.174442
# Coskb
# occopt
l= 3 -1.038935
#
SIESTA
ABINIT
3
# -1.039329
Occupation Option
#
l= 0
0.304116 metal) 0.305046
# Fermi-Dirac smearing #
(finite-temperature
# Coskb
Also, as explained
in
# the
l= Tutorial
1
0.208262
tsmear
0.01
# Temperature
of smearing
(in
Ha) 0.207625
#
SIESTA
ABINIT
l= 2 -0.703574
-0.705726
“Convergence
of electronic
and structural#properties
of a metal with
respect to
#
l= 0
0.304116
0.305046
#
l=
3
-0.000012
-0.000012
#
l= 1
0.207625
0.208262
#
#
l= 2 -0.703574
we should look
l= 3 -0.000012
ElectronicTemperature 0.02 Ry
the k-point sampling: bulk Al”
Energy and
notFreeEng
to the *.out
Kohn-Sham
energy
grep
> Al.siesta.latcon.dat
at-0.705726
the Free
-0.000012
grep FreeEng *.out > Al.siesta.latcon.dat
grep Total Al.out > Al.abinit.latcon.dat
grep Total Al.out > Al.abinit.latcon.dat
occopt
3
# Occupation Option
# Fermi-Dirac smearing (finite-te
Running the energy versus cutoff energy in ABINIT:
bulk Al (a sp metal)
Lattice constant 3.97 Å
-56,9
Free Energy (eV)
-57,0
-57,1
DZP
-57,2
0
5
10
20
15
Cutoff energy (Ha)
25
30
Comparing the energy versus lattice constant in
SIESTA and ABINIT: bulk Al (sp metal)
Basis set of SIESTA: DZP optimized with a pressure of 0.001 GPa at the theoretical
lattice constant of 3.97 Å)
Plane wave cutoff in ABINIT: 8.97 Ha
-57.14
Siesta (DZP)
Abinit (8.97 Ha)
Free Energy (eV)
-57.15
-57.16
-57.17
-57.18
3.85
3.90
4.00
3.95
4.05
Lattice constant (Ang)
4.10
Comparing the pseudopotential in SIESTA
and ABINIT: bulk Au (a noble metal)
To be totally sure that we have run SIESTA and ABINIT with the same
peudopotential operator, i.e. with the same decomposition in local part and
Kleinman-Bylander projectors:
I.
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC F
Cite1 .
# Ekb (in Ha):
#
SIESTA
#
l= 0
1.432475
#
l= 1
0.827346
#
l= 2 -3.174442
#
l= 3 -1.038935
#
# Coskb
#
SIESTA
#
l= 0
0.304116
#
l= 1
0.207625
#
l= 2 -0.703574
#
l= 3 -0.000012
ABINIT
1.432472
0.827336
-3.174442
-1.039329
ABINIT
0.305046
0.208262
-0.705726
-0.000012
grep FreeEng *.out > Al.siesta.latcon.dat
grep Total Al.out > Al.abinit.latcon.dat
occopt
3
tsmear
0.01
# Occupation Option
# Fermi-Dirac smearing (finite-temperature metal)
# Temperature of smearing (in Ha)
Running the energy versus cutoff energy in ABINIT:
bulk Au (a noble metal)
Lattice constant 4.08 Å
-903.6
-903.8
Free Energy (eV)
-904.0
-904.2
-904.4
DZP
-904.6
-904.8
-905.0
20
30
40
Cutoff energy (Ha)
50
60
Comparing the energy versus lattice constant in
SIESTA and ABINIT: bulk Au (a noble metal)
Basis set of SIESTA: DZP optimized with a pressure of 0.02 GPa at the theoretical
lattice constant of 4.08 Å
Plane wave cutoff in ABINIT: 17.432 Ha
-904.66
Siesta (DZP)
Abinit (17.432 Ha)
Free Energy (eV)
-904.67
-904.68
-904.69
-904.70
4.00
4.10
4.05
Lattice constant (Ang)
4.15
Comparing the pseudopotential in SIESTA and
ABINIT: bulk Fe (a magnetic transition metal)
To be totally sure that we have run SIESTA and ABINIT with the same
peudopotential operator, i.e. with the same decomposition in local part and
Kleinman-Bylander projectors:
I.
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC FIE
Cite1 .
# Ekb (in Ha):
#
SIESTA
#
l= 0
1.693954
#
l= 1
0.682661
#
l= 2 -7.130128
#
l= 3 -0.827854
#
# Coskb
#
SIESTA
#
l= 0
0.237303
#
l= 1
0.207724
#
l= 2 -0.711799
#
l= 3 -0.000008
XC.functional
XC.authors
SpinPolarized
ABINIT
1.693938
0.682647
-7.130129
-0.828135
ABINIT
0.238031
0.208363
-0.713975
-0.000008
GGA
PBE
.true.
nsppol
spinat
ixc
2
0.0 0.0 1.0
11
occopt
3
tsmear
0.01
# Number of spin polarizations
# Spin for atoms
# Integer for exchange-correlation choice
# Occupation Option
# Fermi-Dirac smearing (finite-temperature metal)
# Temperature of smearing (in Ha)
Comparing the pseudopotential in SIESTA and
ABINIT: bulk Fe (a magnetic transition metal)
I.
2
ENERGY FUNCTIONAL FOR A DIELECTRIC INSIDE AN ELECTRIC FIELD.
ENERGY FUNCTIONAL FOR A I.DIELECTRIC
INSIDE AN ELECTRIC FIELD.
For the case of metallic system,1 besides the k-point sampling we have to pay
Cite .
particular attention
to the occupation option.
Cite1 .
# Ekb (in Ha):
Now,
besides: ABINIT
#
SIESTA
# Ekb (in Ha):
2.129660
- # TheSIESTA
system is spin polarized
ABINIT # l= 0 2.129660
#
l= 1
1.425393
1.425393
l= 0use2.129660
2.129660 # l= 2 -6.283666
- # We
a GGA functional
-6.283666
l= 1
1.425393
1.425393 # l= 3 -0.825350
-0.825350
- ## We
include
non-linear
partial core corrections
in the pseudo
l= 2 -6.283666
-6.283666 #
#
l= 3 -0.825350
#
# Coskb
#
SIESTA
#
l= 0 S0.260207
IESTA
#
l= 1
0.174652
#
l= 2 -0.683298
#
l= 3 -0.000005
XC.functional
XC.authors
SpinPolarized
-0.825350 # Coskb
#
#
#
ABINIT #
0.261003 #
0.175186
-0.685387 XC.functional
XC.authors
-0.000005 SpinPolarized
GGA
PBE
.true.
nsppol
spinat
ixc
2
0.0 0.0 1.0
11
occopt
3
tsmear
0.01
l=
l=
l=
l=
SIESTA
0
0.260207
1
0.174652
2 -0.683298
3 -0.000005
ABINIT
0.261003
0.175186
-0.685387
-0.000005ABINIT
GGA
PBE
.true.
nsppol
spinat
ixc
2
0.0 0.0 1.0
11
occopt
3
# Number of spin polarizations
# Spin for atoms
# Integer for exchange-correlation choice
# Occupation Option
# Number of spin# polarizations
Fermi-Dirac smearing (finite-temperature metal)
#tsmear
Spin for
0.01 atoms# Temperature of smearing (in Ha)
# Integer for exchange-correlation choice
# Occupation Option
# Fermi-Dirac smearing (finite-temperature metal)
# Temperature of smearing (in Ha)
Running the energy versus cutoff energy in ABINIT:
bulk Fe (a magnetic transition metal)
Lattice constant 2.87 Å
-765.0
Free Energy (eV)
-770.0
-775.0
-780.0
20
DZP
30
40
50
Cutoff energy (Ha)
60
Comparing the energy versus lattice constant in
SIESTA and ABINIT: bulk Fe(magnetic transition metal)
Basis set of SIESTA: DZP optimized without pressure at the experimental lattice
constant of 2.87 Å
Plane wave cutoff in ABINIT: 34.82 Ha
-782.48
Siesta (DZP)
Abinit (34.82 Ha)
Free Energy (eV)
-782.50
-782.52
-782.54
-782.56
-782.58
-782.60
2.75
2.80
2.90
2.85
2.95
Lattice constant (Ang)
3,00