- Email: [email protected]
An ab initio Hartree-Fock investigation of galena (PbS)
/
2 August 1996
ELSEVIER
CHEMICAL PHYSICS LETTERS Chemical Physics Letters 257 (1996) 627-632
An ab initio Hartree-Fock investigation of galena (PbS) M. Mian a, N.M. Harrison b, V.R. Saunders b, W.R. Flavell a a Department of Chemistry, UMIST, P.O. Box 88, Manchester M60 IQD, UK b CCLRC Daresbury Laboratory, Daresbury, Warrington WA4 4AD, UK
Received 23 January 1996; in final form 29 May 1996
Abstract Ab initio Hartree-Fock (HF) theory has been applied to galena (PbS). Pseudopotentials are used to describe the Pb core states and Gaussian basis sets have been developed for S (all electron) and Pb (valence only). The effect of the approximations made (HF theory, the basis set and the pseudopotential) on the structural energetics have been examined in detail. The computed lattice constant and elastic constants are in fair agreement with experiment. The effect of basis set on these properties is minor whereas different choices of pseudopotential and correlation corrections significantly affect the results. The density of states is insensitive to these approximations and in good agreement with spectroscopic data.
1. Introduction Galena (PbS) is a small band-gap semiconductor with the rocksalt structure. There is considerable uncertainty in the experimental elastic constants of PbS [1-6] and their temperature dependence [5,6] (see Table 1), motivating a theoretical study of these properties. We note that unexpected anomalies in its phonon spectrum [7] and a low temperature distortion of the structure at high pressure [8] occur, these issues are not addressed in the present work. Periodic Hartree-Fock (HF) theory [9] including correlation corrections based on density functional theory (DFT) [10] has been very successful in describing the structural and elastic properties of ionic and covalent materials containing light elements. The consequences of the approximations made in ab initio calculations on these materials have been studied in detail [9,33] establishing the accuracy of subsequent studies. Ab initio theory has been little used in
studies of materials containing heavy elements, the only periodic studies involving Pb being of [B-PbF2 I1 l] and PbTe [12] which both used pseudopotentials to describe the core electrons. Pseudopotentials are convenient for reducing the computational cost of calculations and, more critically for heavy elements, they remove the need to perform a relativistic treatment of the core electrons. As far as the author is aware no such treatment has been implemented within periodic Hartree Fock theory. It is timely to establish the reliability of ab initio calculations in these systems. The purpose of the present work is to determine the importance of the approximations made in I-IF calculations on these systems. Namely, the HF approximation, the basis set and the choice of pseudopotential. This has been achieved through the calculation of the lattice parameter and elastic constants of PbS using a variety of basis sets, pseudopotentials and a posteriori DFT estimates of electron correlation corrections.
0009-2614/96/$12.00 Copyright © 1996 Elsevier Science B.V. All fights reserved. Pll S 0009-2614(96)005 9 i -X
628
M. Mian et ai./ Chemical Physics Letters 257 (1996) 627-632
Table 1 A comparison of published experimental elastic constants (in GPa) of PbS at room temperature (RT) and (where temperature dependent data is available) the athermal-limit (AL) Reference no.
Year
Temp.
B
CI I
C44
CI 2
[1] [2] [3]
1951 1951 1970
RT RT RT
59 62 63
102 127 124
25 25 23
38 30 33
[4] [5] [6]
1976 1976 1976 1981
RT RT AL RT
53 57 48 62
126 127 i 55 126
17 25 29 25
16 23 - 5 29
1981
AL
73
152
32
34
[7] [8]
(31") basis sets respectively, was also investigated. DFT correlation corrections to these calculations [18] were examined, using the density functionals of Ceperley-Alder (CA) [19], ColleSalvetti (CS) [20], Perdew (P91) [21] and Becke [22] Each elastic constant was determined from the curvature, at equilibrium, of the energy versus elastic distortion, using a polynomial of 4 th or 5 th order. Each curve was sampled at a minimum of 32 points. We have referred our calculations to the observed elastic constants [6] and lattice parameter [23] extrapolated to the athermal limit.
A posteriori
2. Details o f the calculation
3. Results
Periodic HF theory [13] as implemented in the CRYSTAL92 program [14] was used. The sensitivity of the calculation to various sources of error has been explored in detail. The problem of accurately calculating the Coulomb and exchange series contributions to the Fock operator has been described by others [13], and was addressed by setting tight tolerances for the evaluation of these series (s c • t m = 10 - 7 a n d Sex----Pegx~- 10 - 7 a n d pelx = 1 0 - 1 4 ) , which converged the total energy to approximately 0.1 mhartree. For the Pb atom we used the Durand and Barthelat (DB) [15], and the Hay and Wadt smallcore(HW) [16] pseudopotentials. As there are no f orbitals in our basis set, the f projection of the HW pseudopotential may be expected to have negligible effect and was eliminated. Gaussian basis sets for both the S (86311) (all electron) and Pb (31) (valence only) atoms (see Table 2) were developed (the integers in parentheses refer to the number of Gaussian functions contracted to form a given radial function). These were first optimised in isolated atoms (the HW pseudopotential being used for Pb), the outer exponents of both atoms being subsequently re-optimised in PbS at the experimental lattice constant. For comparison, calculations were performed with Pople's 6-21G [17] S basis set with the exponent of the outer sp shell re-optimised in PbS to 0.134 bohr -2. The re-optimization of the Pb basis set for the DB pseudopotential yielded only a 0.1 mharuee gain in total energy and was therefore not applied. The addition of d orbitals to both the S and Pb atoms, to yield (86311 *) and
The variation of lattice constant ( a 0) and bulk modulus (B) at the HF level with pseudopotential and basis set is shown in Table 3. The calculated a 0 are somewhat larger, and B lower, than found experimentally. There is large variation in a 0 and B with pseudopotential. The DB pseudopotential is in better agreement with the experimental data than the HW pseudopotential. The low value of B(a o) can be seen to be a consequence of the overestimate of a o. B evaluated at the experimental a 0 is consistently 1523% higher than experiment. This is in accord with previous studies [9,10]. There is only a small variation in a 0 and B with basis set. The introduction of d-functions reduces a 0 and increases B slightly. The effect of the correlation corrections is largely independent of basis set, as has been observed in previous studies [24,25], with B increased and a 0 reduced (see Table 4). The strength of this binding term is greater in the P91 and Becke than in the CA and CS functionals. The wide variation in B(a o) is again attributable to the variation of a 0, with B at the experimental a 0 being consistently 12-18% too high. Table 5 shows how the calculated elastic constants compare with experiment. As stated above the choice of pseudopotential has a significant effect on a 0 and B. It can seen that this is also true for C H which increases when the DB pseudopotential is used. Ci2 and C , which do not vary strongly with volume are almost unaffected by changing pseudopotential or the addition of a correlation correction. The inclusion of d orbitals in the basis set has only a small effect. The addition of the correlation correc-
M. Mian et a l . / Chemical Physics Letters 257 (1996) 627-632
629
Table 2 The exponents ( a k in bohr -2 ) and coefficients (Cr, d.p) of the S (86311)/(86311 * ) and Pb (31)/(31 * ) basis sets Shell no.
type
1
s
2
sp
3
sp
4 5 6
sp sp d
1
sp
2 3
sp d
ak
Cs ord
109211.00 16235.206 3573.0286 943.23811 287.26179 99.914226 38.602137 15.531224 281.22171 67.106575 21.794135 8.2097646 3.4178289 1.5452225 4.3752432 1.809620 1 0.6833985 0.237 0.0814 0.3
! 35799 0.94753 0.23092 0.10 0.3
Cp
S (86311)/(86311" ) 0.0002520 0.0019934 0.0111177 0.0498945 0.1661455 0.3627018 0.4108787 0.1457875 - 0.0057780 - 0.0665855 - 0.1 203552 0.2741310 0.6463829 0.2925792 - 0.1750 - 0.5938952 0.8298996 1.0 1.0 1.0 Pb ( 3 1 ) / ( 3 1 " ) 0.32057 -0.81301 0.95198 !.0 1.0
0.0081427 0.0565570 0.2039582 0.3973328 0.3946313 0.1544345 - 0.0613439 0.127225 I 1.2215893 1.0 1.0 -
0.04285 -0.19167 0.64017 1.0 -
tion proposed by Perdew brings the calculations of B and a 0 closer to the observed values. We note that in a plane wave DFT calculation on the isovalent material PbTe the effects of the pseudopotential were
found to be significant. A full elastic tensor was not reported but the calculated value for C44 was 50 GPa in poor agreement with the observed value of 15 GPa [12].
Table 3 The effect of the different pseudopotentials and basis sets (the numbers in braces indicating the basis set used for Pb and S respectively) on the determination o f the lattice constant (a o in ~), and the bulk modulus (B in GPa) determined at the optimized and experimental geometries, of PbS at the I-IF level of approximation. The percentage deviation from the reference values is shown in parentheses
Table 4 The effect a posteriori DFT correlation corrections on the lattice constant (a o in/~), and bulk modulus ( B in GPa) determined at the optimized and experimental geometries, of PbS using the I-IW pseudopotential with the Pb (31 * ) and the S (8631 ! * ) basis sets. The percentage deviation from the reference values is shown in parentheses
Data type
ao
B(a o)
B(5.90 ,~,)
Experiment HW(31,86311}
5.90 6.191(+4.9)
73 51.1(-30)
87.0(+19)
HW{31,6-21G)
6.198(+5.0) 50.4(-30)
90.0(+23)
HW{31",86311"} DB (31, 8631 !}
6.160(+4.4) 6.073(+ 2.9)
82.7(+13) 83.7(+ 15)
51.8(-29) 6 1 . 1 ( - 16)
Damtype
ao
B(a o)
B(5.90~)
Experiment I-IF I-IF+CA HF+CS H F + Becke HF+P91
5.90 6.160(+4.4) 6.056(+2.6) 6.07.4(+2.1) 5.951(+0.9) 5.937(+0.6)
73 51.8(-29) 6 0 . 9 ( - 17) 6 4 . 9 ( - 11) 73.9(+ 1.2) 75.8(+3.8)
82.7(+ 85.8(+ 84.3(+ 82.3(+ 81.9(+
13) 18) 16) 13) 12)
630
M. Mian et a l . / Chemical Physics Letters 257 (1996) 627-632
I ~Pb/'x I P* I!it bt i / /-/ k L" ~ ~. ] ] |,] [ ' ,[ [ '
s UIIJl "/ H , l"[l lJllt~lt [ "xl " ",%'--~" 't
I <..%./; \K,K,~"~\~~~';~
/ff
: " / I
{ {
;',, ¢..'. .'-
.... ' " Pb
_
'"I/ffff
~ -,'.,/,. , ',',lllll :
illllll
Pb 6p
3
Illll
2 I 0
~S 3s
•
-'-.',,,
ru t ', 0
~ T ~ ' ~
x Itlllf/
I I*llll/ ,,11111
lUll
~
~N\~[, I I X~I I' ~
|IUl IIIIIII
,"
,
/" I
g -'"
Pu
Lattice Vector
I
Fig. 1. The deformation electron density of PbS (referred m the neuU'al atoms) in the (100) plane, determined using the (8631 ! * ) for S and the (31 * ) basis sets and HW pesudopotential for Pb at the optimised geometry. Dashed and solid lines represent negative and positive contours respectively. The contour values, whose meduli are indicated by the thickness of the line, are separated b y 10-4[el bohr -3, the narrowest solid line being the zero contour.
g
2 3 I
tl
O~
2 I I•
0L -20
A Mulliken population analysis shows the system to be rather ionic with charge transfer of 1.81el, the overlap population being - 0.0241 e[, indicating weak, if any, covalent bonding. The polarized nature of the solid is confirmed by examination of the electron density deformation (see Fig. 1). The approximately square shaped contours are typical of the addition of equally and oppositely charged spherical distributions of similar size on a rocksalt lattice (see Fig. lb in Ref. [26]). When interpreting I-IF eigenvalues one must be aware that in comparing the orbital energies of occu-
s3p
i
-15
-10 -5 0 5 10 Energy relative to the Fermi-level (eV)
Fig. 2. Projections of the density of states of PbS determined using the HW psedopotential and (86311 * ) and (31 * ) basis sets at the optimised geometry.
pied and unoccupied orbitals the former are stabilised by a self-interaction exchange term leading to a considerable over-estimate of the band gap. Differences in energy between orbitals of equal occupancy are more reliable and thus band widths are often in better agreement with spectroscopy [27]. The calculated density of states (DOS) (see Fig. 2) shows
Table 5 A comparison of the experimental athermal-limit lattice constant (a o in/~) and elastic constants (in GPa) of PbS with the computed values (with percentage deviations in parentheses) showing the effect of changing the basis set and pseudopotential and adding a !)91 a posteriori DFT correlation correction Type of data
ao
B
Ci t
C~
CI2
Experiment HF (HW[31, 86311]) HF (DB[31,86311 ]) HF (HW[31 *, 86311* ]) HF + !)91 (HW[3 I, 86311 ])
5.90 6.191(+ 4.9) 6.073( + 2.9) 6.160(+4.4) 5.975( + 1.3)
73 51( - 30) 61 ( - 16) 52(-29) 75( + 3)
152 106( - 30) 137( - I 0) 110(-28) 179( + 18)
32 32(0) 33( + 3) 33(+3) 32(0)
34 2 4 ( - 30) 23( - 33) 23(-33) 23( - 33)
M. Mian et a l . / Chemical Physics Letters 257 (1996) 627-632
some hybridization of the Pb 6s and 6p orbitals with the S 3p orbitals, indicating that there is weak covalent bonding in this system. The ultraviolet photoelectron spectroscopy (UPS) and photoabsorption of PbS has recently been measured by Santoni et al. [28]. Their spectra show features that can be recognised in the calculated DOS. As expected, the gap in the eigenvalues between the occupied and unoccupied states is considerably larger than the experimental band gap (6.1 eV compared with 0.4 eV). The widths of the fully occupied or unoccupied bands are in much closer agreement. The unoccupied Pb 6p band having a calculated width of 5.9 eV compared with the measured value of 6.0 eV. The occupied S 3p band width being computed as 6.4 eV and observed to be 6.0 eV. There is, however, a significant difference between the computed binding energy of the S 3s band ( = 17 eV) and that reported by Santoni et al. [28] ( = 12 eV). The observed peak may be due to surface contamination. Spectroscopic data for S containing molecules indicates that typical binding energies, relative to the highest occupied molecular orbital (HOMO), for the 3s orbital of S are 16-17 eV. The 3s orbital in H2S, however, has a binding energy, relative to the HOMO, of 11.7 [29]. It is possible that the peak observed by Santoni et al. is the result of the formation of HS-type surface species. UPS is a very surface sensitive technique due to the short inelastic mean free paths of electrons at these energies (typically 10-100 in the 20-100 eV range) [30]. The methods used to check for surface cleanliness (Auger and low energy electron diffraction) would not be able to measure the presence of adsorbed H on galena. Residual gas analysis for ultra-high vacuum systems typically shows H 2 as being the predominant gas present [31 ] and at the quoted pressure in the experimental chamber (2 × 10 -9 mbar) a monolayer of HS species could be formed within 5 minutes of the crystal being cleaved [32]. Unfortunately the region of the UPS spectrum examined by Santoni et al. did not include the 17 eV binding energy region required to validate this hypothesis.
4. Summary and conclusions We have performed a detailed study of the effects of basis set, differing pseudopotentials, and a poste-
631
riori DFF correlation corrections, on ab initio HF calculations of the elastic constants of PbS. The system is largely ionic but with some evidence of covalency in the projected density of states. Basis set details have only a small effect on the lattice parameter and elastic constants. The effects of the pseudopotential however are large. Two apparently similar potentials were found to produce large variations in the lattice parameter (2%) and related variations in B and C1~ (15-20%). Both of these potentials have been shown to perform well in molecular studies and are widely used. It is clear that the reliability of these pseudopotentials in solid state calculations needs to be examined further before these systems can be studied with confidence. Various DFT estimates of the effects of electron correlation also produce significant variations; however the more recently developed functionals (Becke, IX)l) are in encouragingly good agreement. Of the calculations performed the use of the HW pseudopotential and the P91 correlation correction produces the closest agreement with the experimental lattice parameter and elastic constants. The uncertainty in the measured C~2 is large indicating that for this property further experimental work is desirable. The density of states is insensitive to the approximations used and shows features that are consistent with recent spectroscopic measurements; however the S 3s peak may have been assigned incorrectly.
Acknowledgements The work was supported by the Human Capital and Mobility Programme of the European Community under contract CHRX-CT93-0155. MM thanks the EPSRC for financial support during this work.
References [1] G.N. Ramachandran and W.A. Wooster, Acta Cryst. 4(1951) 335. [2] S. Bhagavantam and R.T. Seshagiri, Nature 168 (1951) 42. [3] B.R.K. Gupta and V. Kumar, Sol. Stat. Comm. 45 (1983) 745. [4] G.I. Persada, E.G. Ponyatovskii and Zh. D. Sokoloskaya, Phys. Stat. Sol. 35 (1976) K177. [5] A. Sivaraman, V.C. Padke, E.S. Rajagopal and R.V. Raghavendra, Indian J. Cryo. 1 (1976) 277.
632
M. Mian et aL / Chemical Physics Letters 257 (1996) 627-632
[6] V.C. Padaki, S.T. Lakhshmikumar, S.V. Subrahmanyam and E.S.R. Gopal, Pramana 17 (1981) 25. [7] O.B. Maksimenko and A.S. Mishchenko, Solid State Commun. 92 0994) 797. [8] P.B. Littlewood, J. phys. C. 13 (1980) 4855. [9] R. Dovesi, C. Roetti, C. Freyria-Fava, E. Ap~, V.R. Sannders and N.M. Harrison, Phil. Trans. Roy. Soc. A 34 (1992) 203. [10] R. Orlando, R. Dovesi, C. Roetti and V.R. Saunders, J. Phys. Cond. Matt. 2 (! 990) 7769. [11] M. Nizam, Y. Bouteiiler, B. Silvi, C. Pisani, M. C.ausa and R. Dovesi, J. Phys. C 21 (1988) 5351. [12] K.M. Rabe and J.D. Joannopoulos, Phys. Rev. B 32 (1985) 2302. [13] C. Pisani, R. Dovesi and C. Roetti, ~ - F o c k ab initio treatment crystalline systems (Lecture notes in chemistry, Vol. 48), (Springer, Heidelberg, 1988). [14] R. Dovesi, V.R. Saunders and C. Roetti, CRYSTAL-92 Users Manual, (Universi~ di Torino, 1992). [15] J.C. Barthelat, P.H. Durand and A. Serafani, Mol. Phys. 33 (1977) 159. [16] P.J. Hay and W.R. Wadt, J. Chem. Phys. 82 (1985) 299. [17] W.J. Hehre, R. Ditchfield, R.F. Stewart and J.A. Pople, J. Chem. Phys. 4 (1 970) 2769. [18] M. Cans~. R. Dovesi, C. Pisani, C. Roetti and A. Fortunclli, Phys. Rev. B 36 (1987) 891. [19] D.A. Ceperley and BJ. Alder, Phys. Rev. Letters 45 (1980) 566.
[20] R. Colle and O. Salvetti, Theorct. Chim. Act& 37 (1975) 329. [21] J.P. Perdew, J.A. Chcvary, S.H. Vosko, K.A. Jackson, M.R. Pcderson, D. Singh and C. Fiolhais, Phys. Rev. B 46 0992) 6671. [22] A.D. Becke, J. Chem. Phys. 88 (1988) 1053. [23] Y. Noda, K. Matsomoto, S. Ohaba and Y. Saito, Acta CrysL 43 (1987) 1433. [24] N.M. Harrison and M.L. Leslie, Mo. Sim. 9 (1992) 171. [25] N.M. Harrison and V.R. Saunders, E. Apt& M. Cans~iand R. Dovesi, J. Phys. Cond. Matt. 4 (1992) L261. [26] V.R. Sannders, C. Freyria-Fava, R. Dovesi, L. Salasco and C. Roetti, Mol. Phys. 4 (1992) 629. [27] W.C. Mackrodt, N.M. Harrison, V.R. Saunders, N.L. Allan and M.D. Towler, Chem. Phys. Letters 250 (1996) 66. [28] A. Santoni, G. Paolucci, G. Santoro, K.C. Prince and N.E. Christenscn, J. Phys. Cond. Matt. 4 (1992) 6759. [29] Handbook of Spectroscopy, Vol. l, ed. J.W. Robinson, (CRC Press, 1974) p. 512. [30] M.P. Seah and W. Dench, Surface Intern. Anal. 1 (1979) 2. [31] DJ. Hucknall, Vacuum technology and applications, 1st Ed. (Butterworths-Heinemann Ltd, 1991). [32] D.P. Woodruff and T.A. Delchar, Modem techniques of surface science (Solid state science series), 1st Ed. (Cambridge Univ. Press, 1986), Section 1.2. [33] M.I. McCarthy and N.M. Harrison, Phys. Rev. B 49 (1994) 8574.
2 August 1996
ELSEVIER
CHEMICAL PHYSICS LETTERS Chemical Physics Letters 257 (1996) 627-632
An ab initio Hartree-Fock investigation of galena (PbS) M. Mian a, N.M. Harrison b, V.R. Saunders b, W.R. Flavell a a Department of Chemistry, UMIST, P.O. Box 88, Manchester M60 IQD, UK b CCLRC Daresbury Laboratory, Daresbury, Warrington WA4 4AD, UK
Received 23 January 1996; in final form 29 May 1996
Abstract Ab initio Hartree-Fock (HF) theory has been applied to galena (PbS). Pseudopotentials are used to describe the Pb core states and Gaussian basis sets have been developed for S (all electron) and Pb (valence only). The effect of the approximations made (HF theory, the basis set and the pseudopotential) on the structural energetics have been examined in detail. The computed lattice constant and elastic constants are in fair agreement with experiment. The effect of basis set on these properties is minor whereas different choices of pseudopotential and correlation corrections significantly affect the results. The density of states is insensitive to these approximations and in good agreement with spectroscopic data.
1. Introduction Galena (PbS) is a small band-gap semiconductor with the rocksalt structure. There is considerable uncertainty in the experimental elastic constants of PbS [1-6] and their temperature dependence [5,6] (see Table 1), motivating a theoretical study of these properties. We note that unexpected anomalies in its phonon spectrum [7] and a low temperature distortion of the structure at high pressure [8] occur, these issues are not addressed in the present work. Periodic Hartree-Fock (HF) theory [9] including correlation corrections based on density functional theory (DFT) [10] has been very successful in describing the structural and elastic properties of ionic and covalent materials containing light elements. The consequences of the approximations made in ab initio calculations on these materials have been studied in detail [9,33] establishing the accuracy of subsequent studies. Ab initio theory has been little used in
studies of materials containing heavy elements, the only periodic studies involving Pb being of [B-PbF2 I1 l] and PbTe [12] which both used pseudopotentials to describe the core electrons. Pseudopotentials are convenient for reducing the computational cost of calculations and, more critically for heavy elements, they remove the need to perform a relativistic treatment of the core electrons. As far as the author is aware no such treatment has been implemented within periodic Hartree Fock theory. It is timely to establish the reliability of ab initio calculations in these systems. The purpose of the present work is to determine the importance of the approximations made in I-IF calculations on these systems. Namely, the HF approximation, the basis set and the choice of pseudopotential. This has been achieved through the calculation of the lattice parameter and elastic constants of PbS using a variety of basis sets, pseudopotentials and a posteriori DFT estimates of electron correlation corrections.
0009-2614/96/$12.00 Copyright © 1996 Elsevier Science B.V. All fights reserved. Pll S 0009-2614(96)005 9 i -X
628
M. Mian et ai./ Chemical Physics Letters 257 (1996) 627-632
Table 1 A comparison of published experimental elastic constants (in GPa) of PbS at room temperature (RT) and (where temperature dependent data is available) the athermal-limit (AL) Reference no.
Year
Temp.
B
CI I
C44
CI 2
[1] [2] [3]
1951 1951 1970
RT RT RT
59 62 63
102 127 124
25 25 23
38 30 33
[4] [5] [6]
1976 1976 1976 1981
RT RT AL RT
53 57 48 62
126 127 i 55 126
17 25 29 25
16 23 - 5 29
1981
AL
73
152
32
34
[7] [8]
(31") basis sets respectively, was also investigated. DFT correlation corrections to these calculations [18] were examined, using the density functionals of Ceperley-Alder (CA) [19], ColleSalvetti (CS) [20], Perdew (P91) [21] and Becke [22] Each elastic constant was determined from the curvature, at equilibrium, of the energy versus elastic distortion, using a polynomial of 4 th or 5 th order. Each curve was sampled at a minimum of 32 points. We have referred our calculations to the observed elastic constants [6] and lattice parameter [23] extrapolated to the athermal limit.
A posteriori
2. Details o f the calculation
3. Results
Periodic HF theory [13] as implemented in the CRYSTAL92 program [14] was used. The sensitivity of the calculation to various sources of error has been explored in detail. The problem of accurately calculating the Coulomb and exchange series contributions to the Fock operator has been described by others [13], and was addressed by setting tight tolerances for the evaluation of these series (s c • t m = 10 - 7 a n d Sex----Pegx~- 10 - 7 a n d pelx = 1 0 - 1 4 ) , which converged the total energy to approximately 0.1 mhartree. For the Pb atom we used the Durand and Barthelat (DB) [15], and the Hay and Wadt smallcore(HW) [16] pseudopotentials. As there are no f orbitals in our basis set, the f projection of the HW pseudopotential may be expected to have negligible effect and was eliminated. Gaussian basis sets for both the S (86311) (all electron) and Pb (31) (valence only) atoms (see Table 2) were developed (the integers in parentheses refer to the number of Gaussian functions contracted to form a given radial function). These were first optimised in isolated atoms (the HW pseudopotential being used for Pb), the outer exponents of both atoms being subsequently re-optimised in PbS at the experimental lattice constant. For comparison, calculations were performed with Pople's 6-21G [17] S basis set with the exponent of the outer sp shell re-optimised in PbS to 0.134 bohr -2. The re-optimization of the Pb basis set for the DB pseudopotential yielded only a 0.1 mharuee gain in total energy and was therefore not applied. The addition of d orbitals to both the S and Pb atoms, to yield (86311 *) and
The variation of lattice constant ( a 0) and bulk modulus (B) at the HF level with pseudopotential and basis set is shown in Table 3. The calculated a 0 are somewhat larger, and B lower, than found experimentally. There is large variation in a 0 and B with pseudopotential. The DB pseudopotential is in better agreement with the experimental data than the HW pseudopotential. The low value of B(a o) can be seen to be a consequence of the overestimate of a o. B evaluated at the experimental a 0 is consistently 1523% higher than experiment. This is in accord with previous studies [9,10]. There is only a small variation in a 0 and B with basis set. The introduction of d-functions reduces a 0 and increases B slightly. The effect of the correlation corrections is largely independent of basis set, as has been observed in previous studies [24,25], with B increased and a 0 reduced (see Table 4). The strength of this binding term is greater in the P91 and Becke than in the CA and CS functionals. The wide variation in B(a o) is again attributable to the variation of a 0, with B at the experimental a 0 being consistently 12-18% too high. Table 5 shows how the calculated elastic constants compare with experiment. As stated above the choice of pseudopotential has a significant effect on a 0 and B. It can seen that this is also true for C H which increases when the DB pseudopotential is used. Ci2 and C , which do not vary strongly with volume are almost unaffected by changing pseudopotential or the addition of a correlation correction. The inclusion of d orbitals in the basis set has only a small effect. The addition of the correlation correc-
M. Mian et a l . / Chemical Physics Letters 257 (1996) 627-632
629
Table 2 The exponents ( a k in bohr -2 ) and coefficients (Cr, d.p) of the S (86311)/(86311 * ) and Pb (31)/(31 * ) basis sets Shell no.
type
1
s
2
sp
3
sp
4 5 6
sp sp d
1
sp
2 3
sp d
ak
Cs ord
109211.00 16235.206 3573.0286 943.23811 287.26179 99.914226 38.602137 15.531224 281.22171 67.106575 21.794135 8.2097646 3.4178289 1.5452225 4.3752432 1.809620 1 0.6833985 0.237 0.0814 0.3
! 35799 0.94753 0.23092 0.10 0.3
Cp
S (86311)/(86311" ) 0.0002520 0.0019934 0.0111177 0.0498945 0.1661455 0.3627018 0.4108787 0.1457875 - 0.0057780 - 0.0665855 - 0.1 203552 0.2741310 0.6463829 0.2925792 - 0.1750 - 0.5938952 0.8298996 1.0 1.0 1.0 Pb ( 3 1 ) / ( 3 1 " ) 0.32057 -0.81301 0.95198 !.0 1.0
0.0081427 0.0565570 0.2039582 0.3973328 0.3946313 0.1544345 - 0.0613439 0.127225 I 1.2215893 1.0 1.0 -
0.04285 -0.19167 0.64017 1.0 -
tion proposed by Perdew brings the calculations of B and a 0 closer to the observed values. We note that in a plane wave DFT calculation on the isovalent material PbTe the effects of the pseudopotential were
found to be significant. A full elastic tensor was not reported but the calculated value for C44 was 50 GPa in poor agreement with the observed value of 15 GPa [12].
Table 3 The effect of the different pseudopotentials and basis sets (the numbers in braces indicating the basis set used for Pb and S respectively) on the determination o f the lattice constant (a o in ~), and the bulk modulus (B in GPa) determined at the optimized and experimental geometries, of PbS at the I-IF level of approximation. The percentage deviation from the reference values is shown in parentheses
Table 4 The effect a posteriori DFT correlation corrections on the lattice constant (a o in/~), and bulk modulus ( B in GPa) determined at the optimized and experimental geometries, of PbS using the I-IW pseudopotential with the Pb (31 * ) and the S (8631 ! * ) basis sets. The percentage deviation from the reference values is shown in parentheses
Data type
ao
B(a o)
B(5.90 ,~,)
Experiment HW(31,86311}
5.90 6.191(+4.9)
73 51.1(-30)
87.0(+19)
HW{31,6-21G)
6.198(+5.0) 50.4(-30)
90.0(+23)
HW{31",86311"} DB (31, 8631 !}
6.160(+4.4) 6.073(+ 2.9)
82.7(+13) 83.7(+ 15)
51.8(-29) 6 1 . 1 ( - 16)
Damtype
ao
B(a o)
B(5.90~)
Experiment I-IF I-IF+CA HF+CS H F + Becke HF+P91
5.90 6.160(+4.4) 6.056(+2.6) 6.07.4(+2.1) 5.951(+0.9) 5.937(+0.6)
73 51.8(-29) 6 0 . 9 ( - 17) 6 4 . 9 ( - 11) 73.9(+ 1.2) 75.8(+3.8)
82.7(+ 85.8(+ 84.3(+ 82.3(+ 81.9(+
13) 18) 16) 13) 12)
630
M. Mian et a l . / Chemical Physics Letters 257 (1996) 627-632
I ~Pb/'x I P* I!it bt i / /-/ k L" ~ ~. ] ] |,] [ ' ,[ [ '
s UIIJl "/ H , l"[l lJllt~lt [ "xl " ",%'--~" 't
I <..%./; \K,K,~"~\~~~';~
/ff
: " / I
{ {
;',, ¢..'. .'-
.... ' " Pb
_
'"I/ffff
~ -,'.,/,. , ',',lllll :
illllll
Pb 6p
3
Illll
2 I 0
~S 3s
•
-'-.',,,
ru t ', 0
~ T ~ ' ~
x Itlllf/
I I*llll/ ,,11111
lUll
~
~N\~[, I I X~I I' ~
|IUl IIIIIII
,"
,
/" I
g -'"
Pu
Lattice Vector
I
Fig. 1. The deformation electron density of PbS (referred m the neuU'al atoms) in the (100) plane, determined using the (8631 ! * ) for S and the (31 * ) basis sets and HW pesudopotential for Pb at the optimised geometry. Dashed and solid lines represent negative and positive contours respectively. The contour values, whose meduli are indicated by the thickness of the line, are separated b y 10-4[el bohr -3, the narrowest solid line being the zero contour.
g
2 3 I
tl
O~
2 I I•
0L -20
A Mulliken population analysis shows the system to be rather ionic with charge transfer of 1.81el, the overlap population being - 0.0241 e[, indicating weak, if any, covalent bonding. The polarized nature of the solid is confirmed by examination of the electron density deformation (see Fig. 1). The approximately square shaped contours are typical of the addition of equally and oppositely charged spherical distributions of similar size on a rocksalt lattice (see Fig. lb in Ref. [26]). When interpreting I-IF eigenvalues one must be aware that in comparing the orbital energies of occu-
s3p
i
-15
-10 -5 0 5 10 Energy relative to the Fermi-level (eV)
Fig. 2. Projections of the density of states of PbS determined using the HW psedopotential and (86311 * ) and (31 * ) basis sets at the optimised geometry.
pied and unoccupied orbitals the former are stabilised by a self-interaction exchange term leading to a considerable over-estimate of the band gap. Differences in energy between orbitals of equal occupancy are more reliable and thus band widths are often in better agreement with spectroscopy [27]. The calculated density of states (DOS) (see Fig. 2) shows
Table 5 A comparison of the experimental athermal-limit lattice constant (a o in/~) and elastic constants (in GPa) of PbS with the computed values (with percentage deviations in parentheses) showing the effect of changing the basis set and pseudopotential and adding a !)91 a posteriori DFT correlation correction Type of data
ao
B
Ci t
C~
CI2
Experiment HF (HW[31, 86311]) HF (DB[31,86311 ]) HF (HW[31 *, 86311* ]) HF + !)91 (HW[3 I, 86311 ])
5.90 6.191(+ 4.9) 6.073( + 2.9) 6.160(+4.4) 5.975( + 1.3)
73 51( - 30) 61 ( - 16) 52(-29) 75( + 3)
152 106( - 30) 137( - I 0) 110(-28) 179( + 18)
32 32(0) 33( + 3) 33(+3) 32(0)
34 2 4 ( - 30) 23( - 33) 23(-33) 23( - 33)
M. Mian et a l . / Chemical Physics Letters 257 (1996) 627-632
some hybridization of the Pb 6s and 6p orbitals with the S 3p orbitals, indicating that there is weak covalent bonding in this system. The ultraviolet photoelectron spectroscopy (UPS) and photoabsorption of PbS has recently been measured by Santoni et al. [28]. Their spectra show features that can be recognised in the calculated DOS. As expected, the gap in the eigenvalues between the occupied and unoccupied states is considerably larger than the experimental band gap (6.1 eV compared with 0.4 eV). The widths of the fully occupied or unoccupied bands are in much closer agreement. The unoccupied Pb 6p band having a calculated width of 5.9 eV compared with the measured value of 6.0 eV. The occupied S 3p band width being computed as 6.4 eV and observed to be 6.0 eV. There is, however, a significant difference between the computed binding energy of the S 3s band ( = 17 eV) and that reported by Santoni et al. [28] ( = 12 eV). The observed peak may be due to surface contamination. Spectroscopic data for S containing molecules indicates that typical binding energies, relative to the highest occupied molecular orbital (HOMO), for the 3s orbital of S are 16-17 eV. The 3s orbital in H2S, however, has a binding energy, relative to the HOMO, of 11.7 [29]. It is possible that the peak observed by Santoni et al. is the result of the formation of HS-type surface species. UPS is a very surface sensitive technique due to the short inelastic mean free paths of electrons at these energies (typically 10-100 in the 20-100 eV range) [30]. The methods used to check for surface cleanliness (Auger and low energy electron diffraction) would not be able to measure the presence of adsorbed H on galena. Residual gas analysis for ultra-high vacuum systems typically shows H 2 as being the predominant gas present [31 ] and at the quoted pressure in the experimental chamber (2 × 10 -9 mbar) a monolayer of HS species could be formed within 5 minutes of the crystal being cleaved [32]. Unfortunately the region of the UPS spectrum examined by Santoni et al. did not include the 17 eV binding energy region required to validate this hypothesis.
4. Summary and conclusions We have performed a detailed study of the effects of basis set, differing pseudopotentials, and a poste-
631
riori DFF correlation corrections, on ab initio HF calculations of the elastic constants of PbS. The system is largely ionic but with some evidence of covalency in the projected density of states. Basis set details have only a small effect on the lattice parameter and elastic constants. The effects of the pseudopotential however are large. Two apparently similar potentials were found to produce large variations in the lattice parameter (2%) and related variations in B and C1~ (15-20%). Both of these potentials have been shown to perform well in molecular studies and are widely used. It is clear that the reliability of these pseudopotentials in solid state calculations needs to be examined further before these systems can be studied with confidence. Various DFT estimates of the effects of electron correlation also produce significant variations; however the more recently developed functionals (Becke, IX)l) are in encouragingly good agreement. Of the calculations performed the use of the HW pseudopotential and the P91 correlation correction produces the closest agreement with the experimental lattice parameter and elastic constants. The uncertainty in the measured C~2 is large indicating that for this property further experimental work is desirable. The density of states is insensitive to the approximations used and shows features that are consistent with recent spectroscopic measurements; however the S 3s peak may have been assigned incorrectly.
Acknowledgements The work was supported by the Human Capital and Mobility Programme of the European Community under contract CHRX-CT93-0155. MM thanks the EPSRC for financial support during this work.
References [1] G.N. Ramachandran and W.A. Wooster, Acta Cryst. 4(1951) 335. [2] S. Bhagavantam and R.T. Seshagiri, Nature 168 (1951) 42. [3] B.R.K. Gupta and V. Kumar, Sol. Stat. Comm. 45 (1983) 745. [4] G.I. Persada, E.G. Ponyatovskii and Zh. D. Sokoloskaya, Phys. Stat. Sol. 35 (1976) K177. [5] A. Sivaraman, V.C. Padke, E.S. Rajagopal and R.V. Raghavendra, Indian J. Cryo. 1 (1976) 277.
632
M. Mian et aL / Chemical Physics Letters 257 (1996) 627-632
[6] V.C. Padaki, S.T. Lakhshmikumar, S.V. Subrahmanyam and E.S.R. Gopal, Pramana 17 (1981) 25. [7] O.B. Maksimenko and A.S. Mishchenko, Solid State Commun. 92 0994) 797. [8] P.B. Littlewood, J. phys. C. 13 (1980) 4855. [9] R. Dovesi, C. Roetti, C. Freyria-Fava, E. Ap~, V.R. Sannders and N.M. Harrison, Phil. Trans. Roy. Soc. A 34 (1992) 203. [10] R. Orlando, R. Dovesi, C. Roetti and V.R. Saunders, J. Phys. Cond. Matt. 2 (! 990) 7769. [11] M. Nizam, Y. Bouteiiler, B. Silvi, C. Pisani, M. C.ausa and R. Dovesi, J. Phys. C 21 (1988) 5351. [12] K.M. Rabe and J.D. Joannopoulos, Phys. Rev. B 32 (1985) 2302. [13] C. Pisani, R. Dovesi and C. Roetti, ~ - F o c k ab initio treatment crystalline systems (Lecture notes in chemistry, Vol. 48), (Springer, Heidelberg, 1988). [14] R. Dovesi, V.R. Saunders and C. Roetti, CRYSTAL-92 Users Manual, (Universi~ di Torino, 1992). [15] J.C. Barthelat, P.H. Durand and A. Serafani, Mol. Phys. 33 (1977) 159. [16] P.J. Hay and W.R. Wadt, J. Chem. Phys. 82 (1985) 299. [17] W.J. Hehre, R. Ditchfield, R.F. Stewart and J.A. Pople, J. Chem. Phys. 4 (1 970) 2769. [18] M. Cans~. R. Dovesi, C. Pisani, C. Roetti and A. Fortunclli, Phys. Rev. B 36 (1987) 891. [19] D.A. Ceperley and BJ. Alder, Phys. Rev. Letters 45 (1980) 566.
[20] R. Colle and O. Salvetti, Theorct. Chim. Act& 37 (1975) 329. [21] J.P. Perdew, J.A. Chcvary, S.H. Vosko, K.A. Jackson, M.R. Pcderson, D. Singh and C. Fiolhais, Phys. Rev. B 46 0992) 6671. [22] A.D. Becke, J. Chem. Phys. 88 (1988) 1053. [23] Y. Noda, K. Matsomoto, S. Ohaba and Y. Saito, Acta CrysL 43 (1987) 1433. [24] N.M. Harrison and M.L. Leslie, Mo. Sim. 9 (1992) 171. [25] N.M. Harrison and V.R. Saunders, E. Apt& M. Cans~iand R. Dovesi, J. Phys. Cond. Matt. 4 (1992) L261. [26] V.R. Sannders, C. Freyria-Fava, R. Dovesi, L. Salasco and C. Roetti, Mol. Phys. 4 (1992) 629. [27] W.C. Mackrodt, N.M. Harrison, V.R. Saunders, N.L. Allan and M.D. Towler, Chem. Phys. Letters 250 (1996) 66. [28] A. Santoni, G. Paolucci, G. Santoro, K.C. Prince and N.E. Christenscn, J. Phys. Cond. Matt. 4 (1992) 6759. [29] Handbook of Spectroscopy, Vol. l, ed. J.W. Robinson, (CRC Press, 1974) p. 512. [30] M.P. Seah and W. Dench, Surface Intern. Anal. 1 (1979) 2. [31] DJ. Hucknall, Vacuum technology and applications, 1st Ed. (Butterworths-Heinemann Ltd, 1991). [32] D.P. Woodruff and T.A. Delchar, Modem techniques of surface science (Solid state science series), 1st Ed. (Cambridge Univ. Press, 1986), Section 1.2. [33] M.I. McCarthy and N.M. Harrison, Phys. Rev. B 49 (1994) 8574.