Ab initio computations of periodical (1D, 2D, and 3D) systems can be made much more efficient with application of pseudopotentials. The pseudopotentials replace nuclei and inner electrons. This reduces the number of electrons, which have to be described in simulations, and reduces size of basis sets needed to expand oneelectron wave functions leading to significant savings of computational time.
One of codes used in Materials and Process Simulation Center (MSC) is SeqQuest [1]. Pseudopotentials cannot be inputted to SeqQuest in usually applied parameterized form. SeqQuest code requires input of pseudopotentials in a specific format as numerical arrays. The last requirement makes more difficult application of pseudopotentials, described in literature. However, it is still possible to generate input data for pseudopotentials from published pseudopotentials, or from newly created pseudopotentials. The second approach has an advantage that one can generate pseudopotentials with the same computational approach, which then will be used for actual simulations. The most of published pseudopotentials were generated within HartreeFock approach. However, most of ab initio calculations are made with different flavors of Density Functional Theory. In particular, SeqQuest allows to produce computations in Local Density Approximation (LDA) [2] and with PerdewBerkeErzenkorf (PBE) exchangecorrelation functional [3]. Pseudopotentials based on these functionals can be generated with the code written in FritzHaber Institute [4].
SeqQuest code employs local basis sets from Gaussiantype atomic orbitals. Oneelectron wave functions are built as linear combinations of the atomic orbitals (LCAO). Before ab initio calculations with local basis sets can be performed, one has to optimize the basis sets for all atoms in a system under study.
Before
this project has started Hamman type [5] pseudopotentials were generated and
basis sets were optimized for majority of s, p and d elements using LDA
functional. For some elements the same kind of pseudopotentials and basis sets
were prepared with PBE functional as well.
The goal of this project
was to generate pseudopotentials and to optimize basis sets for the rest
of s, p and d atoms. Simultaneously, some of obtained pseudopotentials and
basis sets were tested by calculations of Equations of State for different
phases of several metals ( Sr, Zr, Pb, Y, Ru, Co ), metal oxides ( SrO, BaO,
PbO ) and lead titanate perovskite.
To generate Hamman type [5] normconserving pseudopotentials I applied the package of computer programs written in FritzHaber Institute [4]. I generated pseudopotentials for all s, p, and d elements (excluding noble gases), which have not been generated before this job started. Details, how to use this package at MSC can be found at [6]. Parameters used to generate the pseudopotentials are collected in Table 1 and also can be found at [7].
As a reference state I used neutral atoms in the ground state for halogen atoms and positively charged ions for all other cases. Charges of the ions sometimes were fractional. Such a choice was usually made for alkali and alkali earth elements. The choice to use positive ions as reference states for metals was made because in the most of chemical compounds atoms of these elements have positive charge. It is obvious, that the made choice is closer to situations, in which we need a good description of metal atoms and will not harm computations of metallic solids, where atoms of these elements are neutral.
Leaving outer core electron shells outside pseudopotentials allows significantly better description of bond lengths and binding energies. Oneelectron wavefunctions of these electrons are still sensitive to chemical environment of atoms because valence dorbitals have large magnitude in the same space, where wavefunction of the outer core electrons are. Valence electrons of alkali and alkali earth atoms usually are nearly completely ionized in chemical compounds. In such situations porbitals of the outer core electron shells are in direct contact with neighboring atoms. Therefore I tried, where it was possible, to leave porbitals from outer core electronic shells of atoms unfrozen. This approach was applied for transition elements and alkali and alkali earth elements, except of the first row atoms, where there are not such orbitals. It is also impossible to apply this approach for pelements, because the employed pseudopotential generation code does not allow for more than one free orbital in each orbital moment component of the pseudopotential. For example, generated pseudopotentials for atoms from Sr till Cd have free 4s, 3p, and 3d orbitals, and pseudopotential for atoms from Tl till Bi have free 5s, 4d, and 5p orbitals.
Matching radii determine distances from nuclei where pseudoorbitals and their logarifmic derivatives have to coincide with atomic orbitals. I applied usual choice for the matching radii. Their values have to be between the last nodes of the valence orbitals and its last maximums. Also, the radii cannot be chosen too close to the last node of valence oneelectron orbitals to avoid irregular behavior of pseudopotentials. Radii, used for generation of pseudopotentials are collected in Table 1 and at website [7].
Exchangecorrelation interactions
are not linear in respect to an electron density. This nonlinearity is more
pronounced in openshell systems and for alkali and alkali earth metal atoms
[8]. The regular semilocal
pseudopotentials have only linear interaction with pseudowavefunctions. The
nonlinear part of interactions between core and valence electrons is not
accounted for. To correct the situation, a model core electron density ρ_{m.core
}can be introduced [8] with the condition, that it coincides with core
electron density at distances from nuclei beyond of cutoff radius r_{nlc}. Then density ρ used for calculation of
exchangecorrelation energy is
ρ = ρ_{m.core } + ρ_{pseudo } ,
where
ρ_{pseudo} is electron density calculated with
pseudowavefunctions. Applied cutoff radii are shown in the Table 1 and at
[7].
Application
of new pseudopotentials in ab initio codes employing local Gaussiantype basis
sets requires a suitable basis sets. Such basis sets can be obtained by fitting
atomic pseudoorbitals, built during generation of pseudopotentials. Another
approach is directly to optimize basis sets in a representative system. The
technique utilized in this project combined both these approaches. At the
first, atomic pseudoorbitals were fitted with several Gaussian atomic orbitals.
Then the obtained expansions into Gaussian orbitals were put to input for
SeqQuest as contracted basis functions. One or two most diffuse Gaussian
exponents and their expansion coefficients in atomic orbitals were optimized in
elementary solids. At the final step the most diffuse Gaussian was set as a
separate basis function as well. In some cases, especially when the most
diffuse Gaussian function has very small exponent parameter, two most diffuse
Gaussians were splitted as a separate
basis function and exponent and expansion coefficients were optimized.
Basis
sets were optimized in elementary solids. In the majority of cases I used the
most stable phase of the elementary for basis set optimization. However, there
were several exceptions, where more symmetric phases with smaller number of
atoms in unit cell were used. Basis sets for magnetic metals ( Co .. Fe) were
optimized in ferroelectric states.
Also, basis sets for halide atoms were optimized in potassium halides KX
( X= Cl, Br, I ). Lattice constants in all employed crystals had experimental
values taken from [9].
Basis
set optimization was performed with help a separate optimization code, written
previously for different purposes. I modified interface of the code to make the
code to work with SeqQuest. This code performs Conjugated Gradients
optimizations. It basically contains two drivers. The one drives computations
of gradients with numerical differentiation over 3 or 5 points for each
variable. The second driver controls Conjugated Gradient optimization. The code
creates necessary inputs for another code ( for example, for SeqQuest ), using
prepared templates. Then it runs a single point calculation of another code and
extracts final energy from output of this code. The extracted energy is used
for calculation of gradients and in energy optimization procedure.
SeqQuest uses difference between actual electron density and electron density calculated as a sum of reference atomic densities to reach significant acceleration of computations. The reference atomic densities are given in SeqQuest input through occupation numbers of atomic orbitals. The criterion for reference state quality is the nonspherical part of electrostatic energy of interaction of electron density ( twoelectron interaction). When this contribution is small, approximations made in energy computation in SeqQuest work better. Therefore, to improve performance of computations with SeqQuest it was necessary to minimize the nonspherical electrostatic energy with respect to the occupation numbers. This optimization was done with help of the same optimization code, which was used for optimization of basis sets. A few modifications needed for this job were made to the code. The first, modification was introduced into interface of the code to extract nonspherical electrostatic energy from the code output instead of the total energy. The second modification was necessary to maintain sum of the occupation numbers constant. This was needed to keep atom charge constant during the optimization.
Results of
generation of pseudopotentials can be found at website [7] and in the
directory /source/seqquest/atoms_pbe in MSC’s computer system . These websites
contain links to .atm files with pseudopotentials, optimized basis sets, and
occupation numbers. These files can be also found in *_pbe_h subdirectories
of /source/seqquest/atoms_pbe directory.
The website [7] contain also detail information, how pseudopotentials
were generated.
Several
pseudopotentials ( Sr, Zr, Pb, Y, Ru, Co )
were tested in more details by
computations of Equations of States ( EoS) for Sr, Zr, Pb, Y, Ru, and Co
metals, SrO and ZrO_{2} metal oxides, and PbTiO_{3} perovskite.
EoS were calculated for bcc, fcc, simple cubic, diamond, a15, and hcp phases of
Sr, Zr, Pb, Y, and Ru metals, for bcc, fcc, simple cubic, diamond, and hcp
phases of ferromagnetic ( fully
spinpolarized ) Co, for rocksalt phase of SrO, for tetragonal phase of PbO,
for fluorite cubic, tetragonal, monoclinic and two orthorhombic phases of ZrO_{2},
cubic, tetragonal and rhombohedral perovskite phases of PbTiO_{3}. EoS
computations were performed with SeqQuest code with newly generated
pseudopotentials for Sr, Zr, Pb, Y, Ru, Co and previously generated
pseudopotentials for Ti [10] and O [11] . Integration grids and Monkhorst net applied for in the EoS calculations for
each material and phase is given at the web site [7]. Sizes of integration
grids were made to reach good convergence of energy, stresses, and forces, and
keep volume per a node of the grid approximately the same in the different
phases of the same materials. To
calculate EoS I first calculated equilibrium lattice constants and atomic
positions in unit cell at zero pressure. Then lattice constants and atomic
positions were optimized at several different volumes of crystal unit cell. At
the end EoS were written in form of energy
as function of WignerSeitz radius, where WignerSeitz radius is a
radius of a sphere with volume equal to volume per formula unit in the
respective crystals ( volume per atom in metals, volume per SrO, ZrO_{2}
in the oxides, volume per PbTiO_{3} in the studied perovskite).
On the next stage obtained EoS were
fitted with analytical expression in form of Rose’s EoS [12] :
_{} (1)
where E_{0}
is the atomization energy, a_{ws}
is WignerSeitz radius, a_{ws}^{(0)} is its equilibrium value at zero pressure, s is scale factor, and
k is a parameter. All parameters of the Rose’s EoS were optimized with the Least
Square Deviation technique to reproduce crystal energies at sets of volumes
where crystal energies were calculated with SeqQuest code. Representation of EoS in the analytical form
allows for further analysis of material properties. Bulk modulus B and its
derivative B’ over pressure can be calculated as
_{} (2)
_{}. (3)
The function of
pressure versus WignerSeitz radius can be calculated as
_{}, (4)
and then for
calculation of enthalpy one can use its definition:
_{}.
(5)
Together equations 4 and 5 define enthalpy function in respect to pressure in parametric form. The minimal possible pressure from equation 4 gives the pressure at which the crystal ruptures (rupture pressure). For the material, for which several phases were calculated, it was possible to consider relative stability of phases and to estimate phase transition pressures at zero temperature. Points of phase transition can be found by solving a system of equations
_{} (6)
in
respect to WignerSeitz radii a_{ws,1
}and a_{ws,2 } for both
phases. The described fitting and analysis of Equations of States was done with
Mathematica code [7].
Results of the
described studies can be found through links from the website [7], where all
described values tabulated and images of calculations with Mathematica code are
available.
Lattice
constants for calculated metals and oxides are in good agreement with
experiment. They are slightly overestimated usually, except Sr. For noncubic
phases of all crystals c/a ratio was usually overestimated as well. The best
coincidence with experimental values was found for yttrium metal. This is
typical results for PBE functional. Calculated bulk moduli of the crystals are
in very well consistent with experimental data.
In this project I generated pseudopotentials for 46 s, p and d elements, optimized basis sets and occupation numbers of reference states for these elements . The pseudopotentials were tested by computations of basic crystal properties.
Table 1. Parameters, used to generate pseudopotentials.
Atomic No 
Element 
Core configuration 
Reference Valence configuration 
Matching Radii, Bohrs 
Model density Cutoff radius, Bohrs 

s 
p 
d 
f 

3 
Li 
0s^{2} 
1s^{0.2}1p^{0} 
1.63 
0.35 
 
 
0.20 
17 
Cl 
0s^{2}1s^{2}1p^{6} 
2s^{2}2p^{5} 
0.744 
0.78 
0.84 
 
0.235 
19 
K 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{0.2}2d^{0} 
2.1 
0.65 
0.7 
 
0.40 
20 
Ca 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{0.2}2d^{0} 
1.77 
0.59 
0.3 
 
0.37 
21 
Sc 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{1}2d^{1} 
1.72 
0.58 
0.43 
 
0.57 
23 
V 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{1}2d^{3} 
1.40 
0.45 
0.30 
 
0.40 
24 
Cr 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{1}2d^{4} 
1.47 
0.47 
0.33 
 
0.46 
25 
Mn 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{1}2d^{5} 
1.41 
0.45 
0.30 
 
0.40 
26 
Fe 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{2}2d^{5} 
1.20 
0.40 
0.25 
 
0.40 
27 
Co 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}2d^{6}3s^{2} 
1.20 
0.40 
0.25 
 
0.40 
29 
Cu 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{1}2d^{9} 
1.20 
0.40 
0.25 
 
0.40 
30 
Zn 
0s^{2}1s^{2}1p^{6}2s^{2} 
2p^{6}3s^{1}2d^{10} 
1.15 
0.35 
0.23 
 
0.35 
31 
Ga 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6} 
3s^{1}2d^{10}3p^{1} 
1.20 
1.10 
0.10 
 
0.265 
33 
As 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6} 
3s^{2}2d^{10}3p^{2} 
1.05 
1.05 
0.10 
 
0.25 
34 
Se 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6} 
3s^{1}2d^{10}3p^{3} 
0.98 
1.03 
0.10 
 
0.245 
35 
Br 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6} 
3s^{2}2d^{10}3p^{5} 
0.98 
0.99 
0.10 
 
0.235 
37 
Rb 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10} 
3p^{6}3d^{0}4s^{0.2} 
2.20 
0.80 
1.35 
 
0.8 
38 
Sr 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{0.1}4s^{0.5} 
2.00 
0.8 
1.0 
 
2.5 
39 
Y 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{1}4s^{1} 
1.88 
0.74 
1.04 
 
1.20 
40 
Zr 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{2}4s^{2} 
1.14 
1.53 
0.635 
 
2.02 
41 
Nb 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{3}4s^{2} 
1.85 
0.60 
0.85 
 
1.6 
43 
Tc 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{5}4s^{2} 
1.50 
0.55 
0.70 
 
1. 
44 
Ru 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{6}4s^{2} 
1.47 
0.54 
0.69 
 
0.9 
45 
Rh 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{7}4s^{2} 
1.45 
0.55 
0.65 
 
0.9 
46 
Pd 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{8}4s^{2} 
1.40 
0.50 
0.60 
 
0.85 
47 
Ag 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{9}4s^{2} 
1.37 
0.50 
0.57 
 
0.80 
48 
Cd 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}3p^{6}2d^{10} 
3p^{6}3d^{10}4s^{2} 
1.34 
0.47 
0.54 
 
0.75 
49 
In 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6} 
3d^{10}4s^{1}4p^{1} 
1.20 
1.50 
0.50 
 
0.75 
50 
Sn 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6} 
3d^{10}4s^{1.8}4p^{1.2} 
1.20 
1.50 
0.50 
 
0.40 
51 
Sb 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6} 
3d^{10}4s^{2}4p^{2} 
1.10 
1.34 
0.51 
 
0.43 
53 
I 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6} 
3d^{10}4s^{2}4p^{5} 
1.01 
1.20 
0.47 
 
0.41 
55 
Cs 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}3d^{10}4s^{2} 
4p^{6}5s^{0.5}4d^{0} 
2.50 
1.05 
1.35 

0.8 
56 
Ba 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}3d^{10}4s^{2} 
4p^{6}5s^{0.5}4d^{0} 
2.20 
0.95 
1.25 
 
0.7 
57 
La 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}3d^{10}4s^{2} 
4p^{6}5s^{1}4d^{1} 
2.15 
0.95 
1.30 
 
2.05 
72 
Hf 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{2} 
1.60 
1.90 
0.90 
1.70 
0.69 
73 
Ta 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{3} 
1.55 
1.85 
0.85 
1.25 
0.67 
74 
W 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{4} 
1.50 
1.80 
0.85 
1.25 
0.67 
75 
Re 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{5} 
1.40 
1.70 
0.80 
1.25 
0.66 
76 
Os 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{6} 
1.38 
1.60 
0.75 
1.25 
0.64 
77 
Ir 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{7} 
1.33 
1.65 
0.70 
1.20 
0.62 
79 
Au 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{9} 
1.30 
1.60 
0.70 
1.10 
0.598 
80 
Hg 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{10} 
1.25 
1.60 
0.60 
1.10 
0.591 
81 
Tl 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{10}5p^{1} 
1.20 
1.50 
0.60 
1.10 
0.58 
82 
Pb 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{10}5p^{3} 
2.00 
0.8 
1.0 
1.53 
3.5 
83 
Bi 
0s^{2}1s^{2}1p^{6}2s^{2}2p^{6}3s^{2}2d^{10}3p^{6}4s^{2}
3d^{10}4p^{6}3f^{14} 
5s^{1}4d^{10}5p^{3} 
1.10 
1.30 
0.60 
1.00 
0.6 
1.
http://www.wag.caltech.edu/group/quest
2. D. M. Ceperley, B.J. Alder, Phys.Rev.Lett. 45 (1980),
566.
3.J.P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77 (1996), 3865.
4. M. Fuchs, M.
Scheffler, Comput. Phys. Commun., Comput. Phys. Commun. 119, 6798
(1999); http://www.fhiberlin.mpg.de/th/fhi98md/fhi98PP/
5. D.
R. Hamann, Phys. Rev. B 40, 2980 (1989).
6. http://www.wag.caltech.edu/group/pseudopotentials.
7. http://www.wag.caltech.edu/home/heifets/MSCpseudo/Reports/PseudopotentialsProjectHome.htm
8. S.G.Louie, S.Froyen, M.L. Cohen, Phys.Rev.B 26 (1982), 1738.
9. CRC Handbook of Chemistry and Physics, 81th edition, Ed. in chief D. R. Lide, CRC Press, 2000 ( Boca Raton, London, New York, Washington, DC).
10. File /source/seqquest/atoms_pbe/ti_pbe/ti_pbe.atm
11. File /source/seqquest/atoms_pbe/o_tm_pbe/o_tm_pbe.atm
12. P. Vinet, J. Ferrante, J.R. Smith, J.H. Rose, J. Phys. C: Sol. St. Phys. 19 (1986), L467L473