To Construct Database for Norm-Conserving Pseudopotentials
K. Kobayashi
National Institute for Research in Inorganic Materials
Namiki 1-1, Tsukuba-shi, Ibaraki 305, Japan()

Pseudopotential database has been prepared as for H, Li, B, C, N, O, Na, Mg, Al, Si,
P, Sc, Ti, V,Cr, Mn, Fe, Co, Ni, Cu, Zn, Y, Zr, Nb, Mo, Tc, Ru, Rh, Pd, Ag, Cd, In, Sn
and Sb. A partial core correction is considered to improve an overestimation of magnetic
moments of transition metals (Fe, Co and Ni). As for prepared pseudopotentials, bulk
properties (lattice constants, bulk moduli, electronic band structures etc.) are checked
and they agree well with other theoretical and experimental results. There are no ghost
bands in present electronic structure calculations using separable forms for non-local partsof pseudopotential.

I. INTRODUCTION
The

rst-principles molecular dynamics (FPMD) developed originally by Car and Parrinello[1] is a very

powerful tool for performing an optimization of internal coordinates of atoms in a unit cell and electronic structures. Appearing this method, a dramatic reduction in storage and computing time is realized. Already, a variety
of systems has been studied by many groups, and many remarkable achievements have been obtained.
A technique of total energy pseudopotential method is generally used in the FPMD. The basic concept ofpseudopotential is a frozen core approximation which assumes that the electronic states of core electrons areinsensitive to the neighboring atomic con guration. This assumption is valid in many cases of solid state physics.
Therefore, avoiding core electrons and smoothing potentials, it is possible to use a plane wave basis set in theelectronic structure calculation.

In addition, it is relatively easy to formulate interatomic forces or stresses

originated from Nielsen and Martin[2] in this approach.
Traditionally, empirical pseudopotentials [3] made in reproducing other (all-electron) theoretical or experimental
results were used. However, such empirical pseudopotentials have a weak point that calculated pseudo charge
densities do not agree with those of the real atoms even outside the core region.
In 1979, Hamann, Schluter and Chiang[4] proposed a norm-conserving pseudopotential. This is constructed
to reproduce eigen-values of valence states and norms of wave functions in all-electron calculations of atoms. A
reliability of norm-conserving pseudopotentials is almost same level as other kinds of the

rst-principles all-electron

calculation (LMTO, APW) in a variety of cases. Furthermore, Bachelet, Hamann and Schluter[5] (BHS) reported adetailed description of making the norm-conserving pseudopotential and tabulated parameters of pseudopotentials
tted by analytic functions. Recently, Rappe, et al. and Troullier and Martin (TM)[6,7] proposed a new normconserving pseudopotential scheme to overcome a problem, in which a huge number of plane waves is required in
the case of calculating

rst elements of atoms and transition metals. This is called 'optimized pseudopotential'.

Using pseudopotentials in the electronic structure calculation, one of the most serious problems is unphysicalghost bands which appear by using a separable form by Kleinman and Bylander[8] in order to reduce a costof the numerical calculation.

There is no way to overcome easily such a problem without any try and error.

Although some technical treatments avoiding ghost bands are given in reports[9,10], it requires a knowledge oftherst-principles calculation and a way of making pseudopotentials. The author thinks that it is ecient anduseful to prepare a database of norm-conserving pseudopotentials as the other way of approach to solve the above
problem. Already, pseudopotentials for over 30 kinds of atoms have been prepared.
In this proceeding, it is described the outline of pseudopotential database and calculational method of theelectronic structure calculation in section II. In section III, results of accuracy for the database (bulk properties,
electronic band structures for transition and noble metals) are presented. Section IV is devoted to summary ofthis study.

1II. METHOD OF CALCULATIONA. Outline of pseudopotential database
At present, norm-conserving pseudopotentials for H, Li, B, C, N, O, Na, Mg, Al, Si, P, Sc, Ti, V, Cr, Mn,
Fe, Co, Ni, Cu, Zn, Y, Zr, Nb, Mo, Tc, Ru, Rh, Pd, Ag, Cd, In, Sn and Sb are available. Pseudopotentials for
Ca, Ga, Ge, Rb and Sr are preparing and testing. Li, Na, Mg, Al, Si and P are constructed by using the wayof BHS[5], the other is for the optimized pseudopotentials by TM[7]. In this database, it is checked that there
are no ghost bands as a result of the electronic structure calculation for solid even by using separable forms[8].
As for transition and noble metals, non-local parts of s, p and d pseudopotentials are used. In this case, a corepotential Vcore (r) is used as a local pseudopotential based on BHS[5]. Parameters of Vcore (r) are appropriately

adjusted to avoid ghost bands in this work. As for the other pseudopotentials(H, Li, Na, B, C, N, O, Mg, Al, Si,
P, In, Sn and Sb), non-local parts are s and p in which the d part of pseudopotential is treated as a local one.
Pseudopotentials for empty states (for examples 4p states of 3d transition metals) are made in a condition of ioniccon guration for valence electrons in the same way of BHS[5]. A scallar relativistic e ect is considered to thesepseudopotentials[11]. This e ect is important in the case of 4d transition metals.

The name of this database is 'NCPS95'(Norm-Conserving PSeudopotential 1995). Calculated results of bulk
properties for transition and noble metals and alloys will be presented in the next section in detail. With regards
to the results of other atoms, another papers of the author and colleagues are available for references[12,13,14,15].

B. Electronic band structure calculation
The present electronic structure calculation is based on the local (spin) density approximation(L[S]DA)[16,17]
in the density functional theory[18,19]. In this study, two types of local (spin) density approximation for exchangecorrelation are used. The electronic part is optimized by using the steepest descent type of algorithm which ismodi ed to accelerate a convergence of the total energy[20].
Although the optimization of only the electronic part is performed in this work, if necessary, it is possible tooptimize interatomic coordinates in a unit cell and unit cell shapes using forces acting on atoms and stresses [2]on unit cell surfaces, respectively, in this computer code.

III. RESULTS AND DISCUSSION
A. Bulk properties and electronic band structures of transition and noble metals
Calculated bulk properties and electronic band structures will be presented as for Sc, Ti, V, Cr, Mn, Fe,
Co, Ni, Cu, Zn, Y, Zr, Nb, Mo, Tc, Ru, Rh, Pd, Ag and Cd. The bulk properties (equilibrium lattice constants
and bulk moduli) are tabulated in Table I. In this table, results of second column are present work, third is the
all-electron band structure calculation[21] and forth is the experimental data of Kittel[22]. Left values of every
 and right ones are bulk moduli(Mbar). Almost results are calculated
column are equilibrium lattice constants(A),in a bcc or fcc structure(see Table I). As for only Zn and Cd, bulk properties and electronic band structures arecalculated in a hcp structure, in which c/a ratios of Kittel[22] are employed. Although the hcp structure is morestable than the fcc structure in many cases, the fcc structure is adopted in order to compare the reference dataof calculation[21]. In addition, it is dicult to optimize c/a ratios of every hcp structure because of computerresources, although it is possible to optimize the unit cell shape (c/a ratio) by using stresses.

20 % for bulk moduli without a few exceptions. In the cases of Ti, V, Fe

In this table, the bulk properties agree with other theoretical and experimental data within 2 3 % for equilibrium lattice constants and 10

and Pd, disagreement is more than 2% for equilibrium lattice constants. Particularly, the lattice constant ofFe(paramagnetic) in this work is 3.9 % shorter than that of experiment. Since the values of bulk moduli depend
tightly on a way of
di erence of about 15

tting and calculational condition(energy cuto , a number of sampling k points, etc), the 20 % is not so serious.

Conditions of calculation for every case are also shown in table I. The energy cutois usually 81Ry. In a few

cases, it is 144Ry(see Table I). A number of sampling k points are 55 (bcc), 89 (fcc) and 95(hcp) in an irreducible
rst Brillouin zone(BZ). The Wigner interpolation formula for exchange-correlation[16] in a paramagnetic phaseis used.
There are no ghost bands as a result of checking calculated electronic band structures. Logarithmic derivatives
for pseudo and real atomic wave functions in all most cases are also checked. Agreements of each other are well2

in almost cases. About a few cases, although there is no problem in electronic band structures, it is necessary toinvestigate a further analysis of calculated results and logarithmic derivatives in future study.
In the next stage, ferromagnetic phases of Fe, Co and Ni are calculated. Recently, some papers[23,24] suggestthat it is necessary to consider a partial core correction (PCC)[25] in order to improve overestimated values ofmagnetic moments. With the PCC (Fe, Co and Ni), calculated results of equilibrium lattice constants, bulk
moduli and magnetic moments are tabulated in Table II. The calculated magnetic moment of Ni without PCC

is 0.8B which is overestimated on the comparison of the other calculated value[21] (0.59B ), and improved with

PCC (0.603B ).

In the paramagnetic phase of Fe, the calculated value of the equilibrium lattice constant isunderestimated in the di erence order of 3.9% on the comparison of the experimental value. On the other hand,the di erence order of the equilibrium lattice constant is improved to 2.3% in the ferromagnetic case. Electronic
band structures for Fe, Co and Ni(ferromagnetic) are shown in Fig. 1 (a),(b) and (c), respectively. The overall
trend of the dispersion curves for above band structures agrees well with the other theoretical results[21,24].
Calculational conditions are as follows: values of energy cuto

and a number of sampling k points in an

irreducible BZ are 81Ry (for Fe, Co, Ni), 121Ry(for Ni), 55 (for Fe[bcc]) and 89(for Co,Ni[fcc]), respectively. The
LSDA formula by Moruzzi, Janak and Williams[17] is used for Fe, Co and Ni.
It is possible to calculate stresses[2] acting on unit cell surfaces in this code of program.

This lead to an

automatic optimization of a unit cell shape under a constant pressure condition. The di erence order of theequilibrium lattice constants from total energies and stresses is about 0.6% in the case of 81Ry. Calculated valuesof stresses and total energies as a function of volume are shown in Fig. 2 for Ni(121Ry, ferromagnetic). It isfound that the agreement with the minimum of total energy and the zero point of stress is improved within 0.3%
for the lattice constant (1.0% for volume). This di erence (0.3%) would be improved in a condition of increasingenergy cuto .
It is necessary to consider a correction term of PCC in the stress calculation. A detailed formulation of thiscorrection will be presented in Appendix.

B. Calculation for alloys (NbMo, CoAl and NiAl)
Alloys of NbMo, CoAl and NiAl are calculated as a CsCl type crystal structure in order to check accuracyof pseudopotentials (related to transferability).

Calculational conditions are as follows.

An energy cutois81Ry, a number of sampling k-points is 84 in an irreducible BZ. Calculated bulk properties (NbMo, CoAl andNiAl) and electronic band structure (CoAl) are shown in Table I and Fig. 3, respectively. In all cases, they aretreated as paramagnetic. Equilibrium lattice constants agree quite well with other theoretical and experimental
results[26,27]. The agreement of the overall trend of dispersion curves between this work and that by the electronic
band structure calculation for CoAl[26] is quite well.

IV. SUMMARY
Database of norm-conserving pseudopotentials for over 30 kinds of atoms has been constructed. The bulk
properties and electronic band structures of them are calculated and checked accuracy. Calculated results givenno ghost bands agree well with other theoretical and experimental results without a few exceptions.
The purpose of this database is as follows. It is possible to perform the optimization of electronic, interatomic
and unit cell shape structures in a system of an arbitrary combination of atoms in a variety of environments(complex compounds, surfaces, interfaces, impurities[defects] and clusters, etc. ), using this database with therst-principles electronic structure calculation. This procedure will lead to prediction (or creating) of new materials,
or support to synthesize new materials in experiment. Also, it will be possible to search a new chemical reaction
process and solid state characteristics. Therefore, thenal purpose of this study is to design new materials byusing this database in the FPMD. The author expects that using this database is the

rst step of this purpose.

ACKNOWLEDGEMENTS
The author wishes to thank Dr. Y. Morikawa (JRCAT) for helpful assistance in the calculation partiallyfor pseudopotential(Li, C, O) optimization. A Computer program for TM pseudopotentials was coded by him.
Helpful discussions of the correction of stresses in PCC with Dr. K. Kusakabe (ISSP) are gratefully acknowledged.
The numerical calculations were performed at National Institute for Research in Inorganic materials (AlphaServer2100 4/200, DEC), and the author also thanks the Supercomputer Center, Institute for Solid State
Physics(ISSP), University of Tokyo for the facilities and use of the FACOM VPP500.

Appendix:Derivation of the correction term of stress for PCC
3

Here, the derivation of the correction term of stress for PCC is described in detail. The following equationis the starting point in a non-spin polarized case.

@c (G)@ (G)@ ( xc (G)(G))= c xc (G) + (G) + (G) xc
;
(1)
@
@ xc
@
where xc is an exchange-correlation energy,  a total charge density, c a charge density of partial core,avolume,  a distortion and G reciprocal lattice vectors. In this equation, v is invariant in the lattice distortion
, in which v is a charge density of valence electrons  = v + c .
The third term at the right hand side of equation (1) is derived as follows.

@xc (r) 1
=
@
=

1

Z

v

v +  cZ(Z@xc (r) 0iGr
e
dr;
@

(2)

@v @c xc 0 xc 0iGr
+
)
e
dr;
@
@ v + c

(3)@xc (G) 1=@

0 xc )e0iGr dr + 1

(xc

Z

@c (r) xc 0 xc 0iGredr:@ v + c

(4)Z@c (r) X @ 1=f(c (r)e0iGr dr)eiGr g
@
@
G
Z
X
X 1
=0
c (G)eiGr +f (02G G )(0 21G ) c (r)(j1(Gr))r3drgeiGr ;

G

(5)

(6)

G

where xc is an exchange-correlation potential.
Therefore,@xc (G)
@
=

1

Z

(xc ((r ))0 xc ((r)))e0iGr dr + 1

where
Pc (G) =The

02G G ( 21G )

Z

Z X
(

G

Gr ) xc 0 xc e0iGr dr;i
Pc (G)ev +  cc (r)(j1 (Gr))r 3 dr:

(7)

(8)

nal formula is@(  (G)(G))@ xc=(xc (G)c (G)0 xc (G))(v (G) + c (G)) + Pc (G)(xc (G) 0 xc (G)) + xc (G) @@:(9)

It is possible to expand this formula in a spin polarized case. The charge density  is divided "" and # , and#also the exchange-correlation potential xc and xc .
In the"#rst term of right hand side at above equation, xc is divided as 2 (xc + xc ).
1

The second term is reformulated as the same way,
(The second term) =Pc (G)(

xc " (G)2+

xc # (G)2References55, 2741(1985).
32, 3792(1985)[3] B. J. Austin, V. Heine and L. J. Sham, Phys. Rev.,127,276(1962)[1] R. Car and M. Parrinello, Phys. Rev. Lett.

[2] O. H. Nielsen and R. M. Martin, Phys. Rev. B

V. Heine, Solid State Physics. Volume 24, New York, Academic Press, 1970.
[4] D. R. Hamann, M. Schluter and C. Chiang, Phys. Rev. Lett.,

43, 1494(1979)
26, 4199(1982)[5] G. B. Bachelet, D. R. Hamann and M. Schluter, Phys. Rev. B

40 xc (G)):(10)

[6] A. M. Rappe, K. M. Rabe, E. Kaxiras and J. D. Joannopoulos, Phys. Rev. B[7] N. Troullier and J. L. Martins, Solid State Commun.,; and Phys. Rev. B,

43, 1993(1991)

74, 13(1990)

41, 1227(1990)

48, 1425(1982)
44, 8503(1991)[8] L. Kleinman and D. M. Bylander, Phys. Rev. Lett.,[9] X. Gonze, R. Stumpf and M. Scheer, Phys. Rev.

51, 16608(1995)[10]A. Khein, Phys. Rev. B[11] D. D. koeling and B. N. Harmon, J. Phys. C

10,3107(1977)

45, 3469(1992)[12] K. Kobayashi, Y. Morikawa, K. Terakura and S. Blugel, Phys. Rev. B[13] K. Kobayashi, K. Kokko, K. Terakura and M. Matsui, in Proceedings of Computer Aided Innovation of New Materials
II 1993 Elsevier Science Publishers B. V.
[14] K. Kobayashi, in Proceedings of Advanced in New Diamond Science and Technology, MYU, Tokyo 1994
[15] K. Kobayashi, in Extended Abstract of JRCAT International Symposium on Nanoscale Self-Organization 1995
[16] W. Kohn and L. J. Sham, Phys. Rev.

140, A1133(1965).
136, B864(1964).

[17] P. Hohenberg and W. Kohn, Phys. Rev.
[18] E. Wigner, Phys. Rev.

46, 1002(1934).

12, 1257(1975)
32,562(1987)[19] J. F. Janak, L. Moruzzi and A. R. Williams, Phys. Rev. B[20] A. R. Williams and J. Soler, Bull. Am. Phys. Soc.,[21] V. L. Moruzzi, J. F. Janak and A. R. Williams, 'Calculated Electronic Properties of Metals'(Pergamon, New York1978)[22] C. Kittle, 'Introduction to Solid State Physics', John Wiley and Sons, Inc., New York, London, Sydney[23] J. Cho and M. Kang, Phys. Rev. B

52, 9159(1995)

52, 12760(1995)
26, 1738(1982)[24] T. Sasaki, A. M. Rappe and S. G. Louie, Phys. Rev. B[25] S. G. Louie, S. Froyen and M. L. Cohen, Phys. Rev. B gu[26] S. Q t and K. M. Rabe, Phys. Rev. B

50, 2075(1994)[27] P. Villars and L. D. Calvert, Pearson's Handbook of Crystallographic Data for Intermetallic Phases (American Societyof Metals, Metals Park, OH, 1985)[28] Calculated by Dr. S. Takizawa in Hokkaido univ. (private communication)

(eV)

(eV)

G

H

N
(a)

P

G

N

(eV)

a = 6.55a.u.

a = 5.305a.u.

W

L

G
(b)

a = 6.59 a.u.

X

W

K

W

L

G
(c)

X

W

K

Fig.1. Calculated electronic band structures for Fe, Co and Ni in ferromagnetic
are plotted at (a), (b) and (c), respectively. The Fermi level is indicated byhorizontal line. Values of lattice constants are indicated in every figure.
Dashed curves are for minority spin and solid curves for majority spin.
(eV)

a = 5.41a.u.

G(a.u.)

Fig. 2. The total energy andstress of Ni(121Ry,ferromagnetic)
as a function of lattice constants
are plotted. Solid curve is forthe total energy, and dashed curvefor the stress.

X

M

G

R

X

Fig. 3. The dispersion curve for CoAlalloy. The Fermi level is indicated byhorizontal line. The lattice constantis 5.41a.u.

5

TABLE I.

Calculated bulk properties (equilibrium lattice constants and bulk moduli) and calculational condi-

tions (energy cuto and crystal structure) for transition and noble metals are tabulated. Unit of lattice constant , bulk modulus Mbar and energy cuto Ry, respectively. Left values are for lattice constants and right values
is A
for bulk moduli in every column. Numerical values at parentheses denoted in a unit of % means the di erenceof order for present values from other theoretical or experimental results. The values of block for hcp Zn and
Cd are lattice constants of a and c, respectively, bulk moduli are not calculated. 'ASW' denotes the all-electron
calculation (ASW), see reference[28].

present work

MJW

Ti(bcc)3.254 1.183

Ti(bcc)3.228 1.152(pcc)

Ti(fcc)4.064

4.00 (1.6%) 1.16

Exp
3.327 (2.2%)
3.327 (3.0%)

V(bcc)

2.999 1.833

2.932 (2.3%) 1.64

3.03 (1.0%) 1.619

Cr(bcc)2.81 2.931

2.805 (0.2%) 2.70

2.88 (2.4%) 1.901

Cu(fcc 81Ry)

3.622 1.62

3.577 (1.3%) 1.55

3.61 (0.3%) 1.37

Cu(fcc 144Ry)

3.64 1.43

3.577 (1.8%) 1.55

3.61 (0.8%) 1.37

Zr(fcc)4.448 0.954

4.403 (1.0%) 0.94

Nb(bcc)3.266 1.81

3.281 (0.5%) 1.68

3.3 (0.1%) 1.702

Mo(bcc)3.144 2.78

3.117 (0.9%) 2.51

3.15 (0.2%) 2.725

3.837 (2.2%) 0.82

Zn(fcc)3.921 0.76

Zn(hcp)(2.66,4.95)

Cd(fcc 81Ry)

4.519 0.506

4.445 (1.7%) 0.461

Cd(fcc 144Ry)

4.507

4.445 (1.4%) 0.461

Cd(hcp 81Ry)(3.032,5.718)

Pd(fcc)3.99 1.7

3.93 (1.6%) 1.7

3.89 (2.6%) 1.808

Ag(fcc)4.157 1.04

4.122 (0.8%) 1.02

4.09 (1.6%) 1.007

Sc(fcc)4.499 0.52

4.493 (0.1%) 0.57

Y(fcc)

4.896 0.44

4.884 (0.2%) 0.33

Y(fcc)

4.880 0.46

4.884 (0.1%) 0.33

Tc(fcc)3.884 3.26

3.852 (0.8%) 2.93

Fe(bcc)2.759 3.14

2.725 (1.2%) 3.06

Co(fcc)3.453 2.61

3.412 (1.2%) 2.84

(within 1%)(within 1.75%)(Ef - 7 eV ghost)2.87 (3.9%) 1.68

Ni(fcc)3.511 2.28

3.466 (1.3%) 2.20

Ru(fcc)3.833 3.38

3.810 (0.6%) 2.89
3.831 (0.9%) 2.61

3.800 (1.8%) 2.70
2.861 (1.6%) 1.62

Rh(fcc)3.867 2.93

NbMo(CsCl)

3.195 2.52

CoAl(CsCl)

2.808 1.89

2.805 (0.1%) (1.99)

NiAl(CsCl)

2.854 1.64

2.868 (0.5%) (ASW)

3.519 (0.2%) 1.86

3.203 (0.2%)

6

TABLE II.

Calculated bulk properties (equilibrium lattice constants[left], bulk moduli[middle] and magnetic

moments[B ,right]) for ferromagnetic Fe, Co and Ni are tabulated.

In the case of Ni(81Ry), the equilibrium

lattice constant and bulk modulus are not calculated. The magnetic moment (Ni,81Ry) is obtained at 3.466 A.
'pcc' indicates considering the partial core correction. The form of table is the same of Table 1, exception for
appending values of magnetic moments. In the column of experimental results, lattice constants (left) and bulk
moduli (right) are only presented.

present work

MJW

Exp

Fe(bcc,pcc)

2.805 2.18 2.25

2.789 (0.6%) 2.17 2.15

2.87 (2.3%) 1.683

Co(fcc,pcc)

3.484 2.64 1.61

3.461 (0.7%) 2.40 1.56

Ni(fcc,81Ry)

0.8

Ni(fcc,81Ry,pcc)

3.487 2.35 0.602

3.466 (0.6%) 2.27 0.59

3.52 (0.9%) 1.86

Ni(fcc,121Ry,pcc)

3.485 2.37 0.603

3.466 (0.6%) 2.27 0.59

3.52 (1.0%) 1.86

7