19254
J. Phys. Chem. B 2006, 110, 19254-19263
First-Principles Electronic Structure Study of the Monoclinic Crystal Bismuth Triborate BiB3O6 Jun Yang* and Michael Dolg* Institute of Theoretical Chemistry, UniVersity of Cologne, Greinstrasse 4, Cologne, D-50939, Germany ReceiVed: June 2, 2006; In Final Form: August 1, 2006
Monoclinic BiB3O6 is an excellent nonlinear optical material with many advantages compared to other borate crystals. The origins of the optical effects and the chemical stability of BiB3O6 are studied with gradientcorrected hybrid B3PW density functional theory within the Gaussian-orbital-based CO-LCAO scheme. Including spin-orbit coupling, the B3PW hybrid functional provides an estimate of the indirect band gap of 4.29-4.99 eV closer to the experimental value of 4.3 eV than HF, LDA, or GGA. The crystal orbital overlap population to give a detailed first-principles analysis of chemical bonding and the density of optical absorptions by convoluting the occupied density of states and the virtual density of states have been calculated. Obvious Bi-O covalent bonds have been found with different energy ranges for 6s-2p and 6p-2p interactions. The reason that [BiO4]5- units are mainly responsible for the optics of BiB3O6 in the long-wavelength region is due to the electronic transfer from occupied O 2p to empty Bi 6p orbitals favored by the Bi-O covalent bonds. The relativistic and correlation effects lead to fundamental differences of the band structure, chemical bonds, and optical effects for BiB3O6 compared with nonrelativistic and uncorrelated calculations.
1. Introduction During the past two decades, increasing scientific interest focused on the crystalline growth1 and excellent physical properties of polar monoclinic bismuth triborate BiB3O6.2 Numerous reports have confirmed that BiB3O6 is an outstanding nonlinear optical (NLO) material3 due to its large effective SHG (Second Harmonic Generation) coefficient, which is higher than that of other borate-based NLO materials such as LiB3O5 and BaB2O4. BiB3O6 is more advantageous than other borates in terms of several aspects such as a wider transparency range, a short ultraviolet absorption wedge, a lower laser threshold, a higher damage threshold, a small beam convergence, and a larger angular acceptance.4,5 Among the series of metal-containing triborates MB3O6, only BiB3O6 crystallizes in space group C2 without inversion centers in the asymmetric crystallographic cells. Its 3-dimensional structure consists of c direction-alternating layers of [B3O6]3rings forming sheets of corner-sharing distorted [BO3]3- triangles and [BO4]5- tetrahedra (Figure 1). These are linked by 6-fold coordinated Bi(III) cations forming [BiO6]9- units in which four oxygen atoms at the same side of Bi(III) are the nearest neighbors (2.390 Å off O4 and O5, 2.086 Å off O6 and O7) and the other two oxygens at the other side of Bi(III) remain relatively far away (2.632 Å off O6′ and O7′) (Figure 2). Nevertheless, the two Bi-O6′ (Bi-O7′) bonds are too long to be accepted as a chemical bond in the usual sense. Obviously O6′ and O7′ are coordinated to the next translated Bi along the b lattice vector as the shortest Bi-O6-like bonds. Therefore, in our paper Bi is considered to be 4-fold coordinated in [BiO4]5units. While the optical properties of BiB3O6 have been widely investigated experimentally, the electronic structure has not been * Address correspondence to the authors at the following: phone 0049221-4706893 and 0049-221-4706886; fax 0049-221-4706896; e-mail
[email protected] and
[email protected].
Figure 1. The asymmetric crystallographic cell of monoclinic BiB3O6 (C2 space group) with the lattice parameters of a ) 7.116 Å, b ) 4.993 Å, c ) 6.508 Å, and β ) 105.62°. The polyhedral [B3O6]3- sheet (dark gray areas) along the c direction contains triangular [BO3]3- (∆) and tetrahedral [BO4]5- (T) units with a ratio of 2∆:1T.
examined in great detail partly because of the tremendously demanding computational complexity of the first-principles
10.1021/jp0634151 CCC: $33.50 © 2006 American Chemical Society Published on Web 09/14/2006
First-Principles Electronic Structure Study of BiB3O6
Figure 2. The experimental 6-fold coordination of Bi(III) in monoclinic BiB3O6. Bismuth is indicated by the sequence number 1. Six-folded oxygens are numbered by 4, 5, 6, 7, 6′, and 7′. Borons in the centers of tetrahedral units [BO4]5- are labeled by 8 (8′) and those in triangular units [BO3]3- are labeled by 9 and 10. Translationally identical atoms are indicated by the same number with a prime as the atoms in the reference primitive cell. The bond distances of Bi-O6′ (Bi-O7′), BiO4 (Bi-O5), and Bi-O6 (Bi-O7) are 2.632, 2.390, and 2.087 Å.
calculations for a crystal with low symmetry such as BiB3O6. Only a few theoretical reports focusing on the evaluation of the optical tensors of BiB3O6 based on different semiempirical models appeared so far.6-8 In general it is simply assumed that the central Bi cation holds a +3 charge by transferring the 6p3 electrons to the [B3O6]3- rings. The remaining 6s2 electrons form on Bi the so-called “inert lone pairs” which are said to stay nonbonding. The 6s2 lone pair electrons at Bi(III) have been suggested to point along the b lattice vector in the opposite direction against the four nearest oxygens around Bi(III) according to the distorted square-pyramidal structure of [BiO4]5units. Without knowing the reasons for the optical properties resulting from electronic structure in detail, the large NLO effects of BiB3O6 are usually explained by the main contributions from the Bi 6s lone pair electrons which are suggested to have more pronounced effects than the [BO3]3- and [BO4]5units. Beyond the optical studies for BiB3O6, several other contributions have attempted to explore the topological rules for the structure stability of polyborate anion type compounds via a statistical analysis of existing structural data.9-12 The [BO4]5tetrahedrons have been considered as favored units in polyborates under high pressures. Unfortunately BiB3O6 turns out to fall into the type for which the structure occurrence rules are least reliable. To our knowledge no first-principles electronic structure calculations on BiB3O6 appeared so far in the literature, except for a “PW-LDA” evaluation of optical tensor components by using a wave function cutting approximation.8 Therefore a detailed depiction of the electronic structure of BiB3O6 is timely and would bring us important insights, as one example of polyborates, to understand the electronic origins of both the structure stability and optical properties by employing quantum chemical calculations. We note that the plane-wave methods in conjunction with pseudopotential and DFT techniques are prevailing in the literature against the application of localized basis functions (Gaussian orbitals) in connection with the CO-LCAO (Crystal Orbital Linear Combination of Atomic Orbitals) scheme in bulk periodic solids. However, several advantages of Gaussian orbitals motivated us to apply the CO-LCAO scheme besides the plane-wave approaches. First of all, a relatively small number of Gaussian orbitals is sufficient for an accurate description of the periodic system, and makes the approach attractive in terms of the computational costs. Due to the closely packed nature of atoms in cells, the most diffuse functions of standard basis sets can be eliminated or replaced by tighter orbitals which significantly shorten the computational time of integral evalu-
J. Phys. Chem. B, Vol. 110, No. 39, 2006 19255 ation without significant loss of accuracy. Second, some popular hybrid functionals such as B3PW and B3LYP, well established for molecules13 to be a crucial improvement in terms of accuracy over pure LDA and GGA functionals, can be used since Gaussian orbitals enable the CO-LCAO scheme to evaluate the HF (Hartree-Fock) exchange interactions analytically. Nevertheless, we performed calculations with different methods from LDA through GGA toward hybrid functionals to establish their accuracies by comparing with the experimental values. Although we applied DFT to calculate BiB3O6 in this paper, it is worthwhile pointing out that by following the HF calculation within CO-LCAO one is in principle able to include the dynamic correlation energy beyond DFT by some highly correlating postHF method such as the recent reports on linear scaling local MP214 or even more advanced FCI (Full Configuration Interaction) and CCSD (Coupled Cluster Single and Double Excitations) in an incremental scheme at the controlled accuracy of several kcal/mol per unit cell.15-17 In practice the BiB3O6 crystal due to its large primitive unit cell containing 10 atoms and its low symmetry can currently not be treated with the wave function-based correlation schemes. The goal of the present work was not to evaluate optical tensors of BiB3O6 but to concentrate on its first-principles electronic structure including correlation effects and relativistic effects originating from Bi cations related to chemical bonding and optical effects. The paper is basically trying to answer the following questions: (1) Where does the optical effect of BiB3O6 originate from the electronic structure point of view? (2) How do the [BiO4],5- [BO3]3- triangles, and [BO4]5- tetrahedrons contribute to the structural stability of BiB3O6? (3) How are the correlation and relativistic effects involved in the electronic structure of BiB3O6? 2. Theoretical Methods Different first-principles methods (i.e., HF, LDA, GGA, and hybrid functionals) were used to optimize the atomic coordinates of BiB3O6 and the accuracies were compared with the experimental bond lengths and angles by running the CO-LCAO-based package CRYSTAL03.18 The optimizations were performed by using a modified conjugate gradient algorithm.19 Analytical gradients were evaluated for both HF and DFT and the optimization convergence was achieved by fulfilling three criteria, i.e., between optimization steps the root-mean-square (RMS) of energy change, gradient change, and displacement change should stay less than 10-7 au, 0.0003 au, and 0.0012 au, respectively. The following tolerances were employed in the evaluation of infinite Coulomb and Hartree-Fock exchange series: 10-7 for the Coulomb overlap, HF exchange overlap, Coulomb penetration, and the first exchange pseudo-overlapand 10-14 for the second exchange pseudo-overlap. The Fock matrix has been diagonalized at 24 k-points within the irreducible Brillouin zone corresponding to a shrinking factor of 4 in the Monkhorst net.20 The larger number of 150 k-points reduces the energy by only 6.4 × 10-6 au. To improve the convergence, a negative energy shift of 1.0 au to the diagonal Fock/KS matrix elements of the occupied orbitals was added to reduce their coupling to the unoccupied set and maintained after the diagonalization. A very accurate extra-large grid consisting of 75 radial points and 974 angular points was employed in the DFT calculations, where Becke grid point weights21 were chosen. To investigate the relativistic effects resulting from the heavy Bi cations, we have used energy-consistent scalar-relativistic as well as nonrelativistic ab initio ECP (effective core pseudo-
19256 J. Phys. Chem. B, Vol. 110, No. 39, 2006
Yang and Dolg
TABLE 1: The Dependence of the Geometry and Band Gap of BiB3O6 on HF, LDA, GGA, and Hybrid Functionals for Scalar-Relativistic Calculations with Lattice Constants Fixed to Their Experimental Values [BiO4]5- (Å) [BO3]3- (Å) [BO4]5- (Å) EnDash∠O6′-Bi-O7′ (deg) ∠O4-Bi-O5 (deg) ∠O6-Bi-O7 (deg) |∆d|Bi-O (Å) |∆d|B-O (Å) |∆-| (deg) band gap (eV)
exptl1b
HF
2.390 2.086 1.411 1.365 1.339 1.487 1.436 102.15 152.20 90.58
2.380 2.061 1.394 1.354 1.336 1.489 1.446 103.18 154.66 89.92 0.018 0.009 1.38 12.54
4.3h
VWNa-LDAb
PWc-LDAb
PBEd
PW91e
B3LYPf
B3PWg
2.405 2.140 1.404 1.368 1.345 1.490 1.444 100.79 149.69 89.18 0.034 0.005 1.63 3.96
2.404 2.140 1.404 1.368 1.345 1.490 1.443 100.89 149.75 89.22 0.034 0.005 1.70 3.96
2.394 2.135 1.412 1.376 1.352 1.505 1.451 99.85 148.80 90.04 0.026 0.011 2.08 4.05
2.394 2.133 1.411 1.375 1.350 1.504 1.450 100.02 149.03 90.01 0.026 0.010 1.96 4.07
2.389 2.105 1.404 1.367 1.344 1.500 1.446 101.17 151.20 90.15 0.010 0.007 0.80 5.79
2.388 2.109 1.405 1.367 1.344 1.498 1.446 101.40 151.32 89.88 0.012 0.007 0.78 5.69
a
LSD. Vosko-Wilk-Nusair parametrization of the Ceperley-Alder free electron gas correlation results.41 b LSD. Dirac-Slater exchange.42 LSD. Perdew-Wang parametrization of the Ceperley-Alder free electron gas correlation results.43-45 d GGA. Perdew-Becke-Ernzerhof.46 e GGA. Perdew-Wang.43-45 f Becke’s three-parameter functional combined with the nonlocal correlation LYP.47,48 g Becke’s three-parameter functional combined with the nonlocal correlation PW.43-45,47 h Reference 3b.
c
potentials) of the Stuttgart-Koeln variety for Bi.22 The 1s-5d shells were included in the ECP core and others (i.e., 6s6p... shells) were left in the valence. The explicit quantum chemical treatment is restricted to the valence electrons and scalarrelativistic effects are usually implicitly accounted for by a proper adjustment of free parameters in the valence model Hamiltonian. The reference data used to determine the spinorbit averaged relativistic potential have been taken from relativistic all-electron (AE) calculations by using the so-called Wood-Boring (WB) scalar-relativistic Hartree-Fock (HF) approach. Both AE WB and ECP calculations have been performed with an atomic finite-difference HF scheme to avoid basis set effects in the determination of the ECP parameters. A (4s4p1d)/[2s2p1d] double-ζ valence basis set22,23 was applied to the Bi atoms in all calculations with energy-optimized outmost orbital exponents of 0.0900 (p orbital) for the scalar-relativistic ECP as well as 0.1063 (s orbital) and 0.0944 (p orbital) for the nonrelativistic ECP. The basis sets of boron and oxygen atoms were originally taken from the public EMSL database.24 Modifications of the exponents of uncontracted outmost orbitals are necessary if they are too diffuse. The triple-ζ Dunning contraction (11s, 6p)/[5s, 3p]25 and 6-311G* basis sets were employed for oxygens and borons, respectively. One additional d polarization function with the optimized exponent of 0.8700 and 0.8268 was added to (11s, 6p)/[5s, 3p] contractions of oxygens for scalar-relativistic and nonrelativistic calculations, respectively. The following process was executed for the basis set optimization: Each exponent was varied with a step size of 0.05 and a single point calculation was completed with CRYSTAL03. At every single optimization iterative stage, each optimal exponent was obtained by a three-point parabolic fitting. The above step was recycled until the variation of SCF energy stayed below 10-6 au. All the basis set optimizations were performed at the hybrid B3PW scalar-relativistic level and the optimized orbitals were directly transferred to HF and other DFT calculations. We assume that these optimized sets for B3PW are also quite close to the optimum for other methods. Fox example, the exponent variation of the oxygen d polarization function from 0.8700 to 0.8200 ends up with the scalarrelativistic HF energy difference of only 3.5 × 10-5 au. The spin-orbit effects were calculated with the plane-wave code WIEN2K26 only to correct the scalar-relativistic LDA and GGA calculations for band gaps. The muffin-tin radii of Bi, O, and B are 2.36, 1.26, and 1.26 au, respectively. The cutoff
energy is -6.0 Ry separating the core and valence configurations of 5d106s26p3, 2s22p4, and 2s22p1 for Bi, O, and B, respectively. The calculations were performed with 1418 LAPWs (Linearized Augmented Plane-wave) and converged with 20 k-points due to the 2-fold convergence criteria of the energy less than 10-6 Ry and the charge center distance less than 10-4 au. The occupied core states were treated with fully relativistic effects (mass-velocity, Darwin, and Spin-Orbit (SO) corrections) by solving the relativistic Dirac-Fock equation,27 and the valence states were individually calculated with scalar relativistic corrections and SO perturbations in a second variational procedure28 by using the scalar-relativistic eigenfunctions as a basis. The SO corrections were applied only to Bi. 3. Results and Discussions 3.1. Accuracies of HF, LDA, GGA, and Hybrid Functionals. Before the discussion of the electronic structure, we present results of scalar-relativistic calculations with different Hamiltonians (i.e., HF, LDA, GGA, and hybrid functionals, see Table 1) to establish their accuracies with respect to experimental data.1a To make parallel comparisons between different methods, the lattice constants should stay invariant with respect to different methods. Thus in all our subsequent calculations, only the inner-cell atomic coordinates have been fully optimized while fixing the lattice constants to their experimental values. To be sure of its accuracy by the experimental lattice constants for the electronic structure calculations, we also performed the lattice optimization by calculating the fully atom-relaxed energies for B3PW with respect to forty cell volumes corresponding to the lattice constants with the variations from +3% to -3% of the experimental ones. The optimized crystallographic lattice constants deviate from the experimental values by less than 0.6% with an energy reduction of only 0.0076 eV. Therefore we assume that the calculations with the experimental lattice constants are sufficiently accurate to describe the band structure of BiB3O6 if other approximations (ECPs, limited basis sets, etc.) are taken into account. The calculated average absolute deviations |∆d|Bi-O for Bi-O bond lengths, |∆d|B-O for B-O bond lengths, and |∆-| for bond angles between the experimental and theoretical values are listed in Table 1. We found that |∆d|Bi-O is considerably greater than |∆d|B-O for all methods especially for LDA, which may be due to the incomplete error cancellation between the Bi ECP and the functionals. For |∆d|B-O there seems to be no evident method
First-Principles Electronic Structure Study of BiB3O6
J. Phys. Chem. B, Vol. 110, No. 39, 2006 19257
Figure 3. The band structures of BiB3O6 over the valence band region obtained from scalar-relativistic B3PW (left), nonrelativistic B3PW (middle), and scalar-relativistic HF (right) calculations.
bias. The hybrid functionals B3LYP and B3PW indeed reduce the deviations |∆d|Bi-O and |∆-| by a factor of 2 with respect to LDA and GGA. Although HF gives similar |∆d|B-O and |∆d|Bi-O values to B3PW and B3LYP with differences of only a few thousandths of an angstrom, the |∆-| value is almost two times larger than that for the hybrid functionals. From a theoretical point of view, HF is not suitable for the electronic structure investigation of BiB3O6 due to its exclusion of correlation effects. Moreover, our previous report29 suggests that, for some cases such as hexagonal ionic lanthanide sesquioxides Ln2O3 (Ln ) La-Pm), B3PW agrees with the experimental geometry (i.e., both lattice constants and internal atomic coordinates) even notably better than B3LYP. Therefore, the B3PW functional is recommended and utilized for all the subsequent calculations of BiB3O6. Nevertheless, the small overestimations of some Bi-O bond lengths even at the B3PW level are at least partially caused by the omission of a corepolarization potential29,30 in our applied ECPs for Bi, i.e., the neglect of static and dynamic core-valence correlation. The scalar-relativistic HF, LDA, GGA, and hybrid functionals give the indirect band gaps of 12.54 eV, 3.96 and 3.96 eV, 4.05 and 4.07 eV, and 5.79 and 5.69 eV (cf. Table 1). Whereas HF overestimates the band gap by three times due to the absence of correlation effects, the LDA and GGA routines underestimate the band gap by about 0.34 and 0.24 eV, respectively, which is in line with the general observation that band gaps are usually underestimated within the DFT scheme. The B3PW and B3LYP functionals overshoot the band gap by 1.39 eV, which is, however, not surprising since B3PW and B3LYP are not pure DFT approaches and mix exact HF and pure DFT exchange with empirically determined weights. Accordingly these hybrid approaches should produce the band gaps located between the results of HF and pure DFT methods. The good performance of GGA and LDA compared to the hybrid approaches has to be reviewed with great care, however. If we take into account the spin-orbit effects of the Bi cations, as we will discuss in the following, the hybrid functionals would bring the calculated band gaps much closer to the experimental result than the LDA and GGA approaches. 3.2. Electronic Structure of BiB3O6. The band structures for BiB3O6 are illustrated in Figure 3 for scalar-relativistic B3PW, nonrelativistic B3PW, and scalar-relativistic HF calcula-
tions. The bands below ca. -2.0 eV refer to the occupied valence states (valence band, VB) and the bands above refer to the virtual states (conduction band, CB). The scalar-relativistic B3PW gives the indirect band gap of 5.69 eV determined between the Z point at the bottom of conduction bands and the E point at the top of valence bands. This suggests the photonphonon coupling for a nonvertical E-Z transition in order to meet the energy conservation mainly fulfilled by photons and the momentum conservation mainly fulfilled by phonons. The band dispersions are significant along the specified k-point paths. As pointed out by Hoffmann,31 band dispersions are determined by the inter-unit-cell overlapping of atomic orbitals, i.e., the overlapping degree sets the bandwidth and the topology of that overlapping provides the way that the bands run up or down. Therefore, the band dispersions along B-D and Z-Γ (i.e., the intercell orbital overlapping toward the direction of the c lattice vector) imply additional notable interactions besides the electrostatic attraction between Bi cations and coordinating oxygens anions, since the layers of [B3O6]3- polyanions are alternating along the c lattice direction connecting the Bi cations between (Figure 1a). Thus we do not support the finding that BiB3O6 has the electronic structure typical for a molecular crystal as was concluded in the literature.8 Other k-point paths in Figure 3 correspond to the strong B-O covalent bonds in both [BO3]3and [BO4]5- units and thus their bands are much more dispersive than the bands for the Bi-O interactions. We explore the density of states projected out for the atomic valence orbitals of bismuths and oxygens in the [BiO4]5- units (Figure 4) as well as for the triangular [BO3]3- and tetrahedral [BO4]5- units (Figure 5). For the scalar-relativistic B3PW calculation, the lower part of the valence states below -24 eV mainly corresponds to the occupied O 2s states, which are well separated by ∼7.7 eV from the occupied O 2p states dominating the upper part of the valence states until the top of the valence states at ca. -5.6 eV. The valence states have contributions also from the admixture of Bi 6s and 6p orbitals contrary to the separation of O 2s and 2p orbitals. Although Bi 6s6p states appear relatively weak, they can be neglected by no means because of their important overlapping with the occupied O 2p states widely spreading over the same region. This indicates that the covalent interactions between the central Bi cations and the coordinating oxygens may take place from the energy
19258 J. Phys. Chem. B, Vol. 110, No. 39, 2006
Yang and Dolg
Figure 4. The l-momentum-projected density of states for the valence orbitals of bismuths and oxygens of [BiO4]5- units of BiB3O6 for scalarrelativistic B3PW (top), nonrelativistic B3PW (middle), and scalar-relativistic HF (bottom) calculations. The projected density of states for Bi 6s, Bi 6p, O 2s, and O 2p orbitals are plotted in solid, dash, dot, and dash-dot lines, respectively.
Figure 5. The atom-projected density of states for the triangular [BO3]3- (dash line) and tetrahedral [BO4]5- (solid line) units of BiB3O6 for scalar-relativistic B3PW (top), nonrelativistic B3PW (middle), and scalar-relativistic HF (bottom) calculations.
matching point of view, which is consistent with the band structure feature. The theoretical measurement of this Bi-O covalent degree will be quantitated in the next section. The lower parts of the conduction states between 0 and 3 eV result from the strong origin of the unoccupied Bi 6p orbitals with only small contributions from other orbitals. The region above 3 eV can be divided into two sectors (Figure 5), i.e., from ∼3 to ∼7 eV consisting of the states of [BO3]3- units
and above ∼7 eV consisting of the states of both [BO3]3- and [BO4]5- units. These unoccupied states result mainly from the boron empty orbitals since the oxygen 2s2p states show little contribution to the conduction bands (Figure 4). 3.3. Chemical Bonds in BiB3O6. We have applied the COOP (crystal orbital overlap population) scheme following the classic introduction by Hughbanks and Hoffmann32 to analyze the chemical bonds in BiB3O6. We note that several other methods
First-Principles Electronic Structure Study of BiB3O6
J. Phys. Chem. B, Vol. 110, No. 39, 2006 19259
Figure 6. Crystal orbital overlap population (COOP) plots for the interactions of Bi-O in [BiO4]5- units (top), B-O in [BO3]3- (middle), and [BO4]5- (bottom) for the scalar-relativistic B3PW calculation. The COOP plots for the nonrelativistic B3PW and scalar-relativistic HF calculations can be found in Figures 3 and 4 in the Supporting Information.
such as COHP33 (crystal orbital Hamilton population) and BCOOP34 (balanced crystal orbital overlap population) are alternative options. COOP describes the density of bonding and antibonding interactions between specific orbitals at a given energy in solids. Regions with positive COOP contributions are bonding, regions with negative COOP contributions are antibonding, and zero COOP contributions are nonbonding. We g projected out the overall density matrix element PuV (E) up to the given energy E for the atomic orbital u belonging to the reference cell (0, 0, 0) and V belonging to the cell g. Those g (E) were evaluated by eq 1. matrix elements puV g puV (E) )
g dPuV (E) dE
(1)
g The COOPRβ between two atom-centered orbital sets R in the reference cell (0, 0, 0), and β in the cell g for BiB3O6 can be calculated by
g (E) ) COOPRβ
g g puV (E)suV dRβ(E) ∑ ∑ u∈R V∈ β
(2)
g where suV is the element of the overlap matrix between the atomic orbital u in the reference cell (0, 0, 0) and V in the cell g, and dRβ is the projected density of states for all orbitals belonging to R and β. Equation 2 is partially different from the original COOP definition since we replace the band-resolved total density of states by the projected density of states dRβ so as to avoid some arbitrary enhancement of COOP, for instance, between the orbital u and V with only small density of states. COOP between the Bi 6s (6p) and O 2p was computed by overlapping corresponding orbitals at Bi1 with those at O4, O5, O6, and O7 in the cell (0, 0, 0), COOP for the triangular [BO3]3units at B9 with O5 in the cell (0, 0, 0) and O2′, O6′ in the cell (1, 0, 0), and COOP for the tetrahedral [BO4]5- units at B8 with O4, O3 in the cell (0, 0, 0) and O2′, O5′ in the cell (0, 0, 1).
The strength and dispersion of Bi-O covalent bonds are quantitatively manifested by scalar-relativistic B3PW COOP in Figure 6. Bi 6s and O 2p orbitals form both separate bonding and antibonding 6s-2p interactions, while only Bi 6p-O 2p bonding interactions spread continuously below the top of the valence band. Therefore it turns out that the Bi 6s orbitals do host active lone pair electrons instead of “inert” ones. The Bi 6s-O 2p bonding should contain dominant Bi 6s components and the Bi 6s-O 2p antibonding should contain dominant O 2p components. The Bi-O bonding and antibonding are located over different energy intervals. For example, for O6 and O7, 6s-2p bonding stays within a narrow energy range from -16.6 to -14.5 eV; 6s-2p antibondings spread from -13.5 eV to the top of valence bands at ca. -5.6 eV. For the nearest O6 and O7 atoms around the central Bi, the 6s-2p and 6p-2p bond intensities are remarkably larger than those for the more distant O4 and O5 atoms because of the fast decay of orbital overlapping with increasing Bi-O distance. This is consistent with the classic chemical concept of short-range covalent bonds. From Table 2, the 6p-2p bonding is found to be significantly important to the stabilization of the [BiO4]5- units in the way that it compensates to a great extent for the 6s-2p antibondings. However, since the 6p-2p compensation is not complete by leaving the net antibonding of -5.29 in [BiO4],5- it is reasonable to conclude that the electrostatic attraction between Bi cations and O anions is still indispensably necessary to assemble together the two layers of [B3O6]3- rings at both sides of Bi (Figure 1a). We will point out in the next section that this unstable net antibonding in [BiO4]5- is caused by the scalarrelativistic effect of the Bi cations. The electronic densities in the O6-Bi-O7 and O4-Bi-O5 planes projected between -16.6 and -14.5 eV are visualized in Figure 7 for Bi 6s-O 2p bondings. We do not visualize the electronic densities at higher energy levels to avoid ambiguity because this region is contributed by both Bi 6p-O 2p bonding and Bi 6s-O 2p antibonding interactions. It is clearly observed that the electrons are concentrated in the middle between Bi
19260 J. Phys. Chem. B, Vol. 110, No. 39, 2006
Yang and Dolg
TABLE 2: Scalar-Relativistic and Correlation Effects on the Bond Densities over the Energy for the Valence States of BiB3O6 bonds
bondinga
antibondinga
net bondsb
scalar-relativistic B3PW Bi-O 6s-2p 6p-2p total in [BiO4]5B-O in [BO3]3B-O in [BO4]5total in BiB3O6c
5.64 18.38 24.02 168.12 23.94 216.08
-29.31 conduction band -29.31 conduction band -7.56 -36.87
-23.67 18.38 -5.29 168.12 16.38 179.21
nonrelativistic B3PW Bi-O 6s-2p 6p-2p total in [BiO4]5B-O in [BO3]3B-O in [BO4]5total in BiB3O6 c
7.25 17.24 24.49 170.78 24.25 219.52
-17.69 conduction band -17.69 conduction band -10.52 -28.21
-10.44 17.24 6.80 170.89 13.73 191.31
scalar-relativistic HF Bi-O 6s-2p 6p-2p total in [BiO4]5B-O in [BO3]3B-O in [BO4]5total in BiB3O6c
4.84 11.61 16.45 100.33 13.70 130.48
-25.78 conduction band -25.78 conduction band -19.02 -44.80
-20.94 11.61 -9.33 100.33 -5.32 85.68
a Bonding and antibonding were given by corresponding integral intensities of COOP over the valence region. b Net bonds were evaluated by the summation of bonding and antibonding. c Totals in BiB3O6 were evaluated by the summation of the bonding, antibonding, and net bonds for Bi6s-O2p, Bi6p-O2p, and B-O in [BiO4],5- [BO3]3-, and [BO4]5units, respectively.
and O atoms, which typically represents covalent bonds between Bi 6s and O 2p orbitals. The covalent bonds in both triangular [BO3]3- and tetrahedral [BO4]5- units are much stronger than those in [BiO4]5- units. From Figure 6 and Table 2, the B-O bonds in [BO3]3- show a bonding character across the valence band and the antibonding stays in conduction bands, which totally leaves the strong net bonding intensity of 168.12. The B-O bonds in [BO4]5- stay bonding from -15.2 to -11.0 eV and antibonding from -8.5 to -6.3 eV ending up with the net bonding intensity of 16.38. Therefore the total net bonding in BiB3O6 adds up to 179.21 so that it is sufficient to achieve the stability of BiB3O6 at ambient conditions significantly contributed by the B-O bondings in polyanions. Becker10 pointed out according to structure data analysis for numerous anhydrous polyborates that the tetrahedral [BO4]5- units are not favored units in polyborates at ambient conditions. Most crystallographic structures of anhydrous polyborates contain more triangular [BO3]3- units than tetrahedral [BO4]5- (cf. Figure 2 in ref 10). One reason can be suggested here due to our calculations for BiB3O6 that a large ratio of [BO4]5- to [BO3]3- may greatly intensify the B-O antibonding interactions in [BO4]5- close to the Fermi level to some extent that only the application of external pressure can stabilize the structure. 3.4. Electronic Origins for Optical Effects of BiB3O6. Equation 3 is formulated to quantitatively describe the density of optical absorptions gauging the contributions from different units in BiB3O6, which is basically the mathematical convolution at a given energy shift pω between the density of states dVB u (E) CB for atomic orbital u in VB and the density of states dV (E) for atomic orbital V in CB.
GuV(ω) )
CB ∫-∞+∞ dVB u (E) dV (E - pω) dE
(3)
Therefore the contributions for optical absorptions from [BiO4]5-, [BO3]3-, and [BO4]5- units can be calculated by eq 4 running over all VB orbitals and those CB orbitals building up [BiO4],5[BO3]3-, and [BO4]5- units, respectively. Θ
GI(ω) )
V
∑u ∑V GuV(ω)
(4)
where Θ ∈ VB, V ∈ {[BiO4]5-CB, [BO3]3-CB, [BO4]5-CB}. The plots with respect to the wavelength are presented in Figure 8 for the scalar-relativistic B3PW calculation. It is found that the absorption of BiB3O6 begins at ∼200 nm. In the long wavelength region from 200 to 130 nm, [BiO4]5- units dominate the optical absorption mainly due to the electronic transfer from occupied O 2p to virtual Bi 6p orbitals at the bottom of the conduction bands between 0 and 3 eV. The contribution of single Bi cations due to the inner-shell electronic transfer from 6s to 6p is much less than that of the [BiO4]5- units since the majority of Bi 6s states stay deep in energy at ca. -16 eV far away from the top of the valence band. The electronic transitions from occupied O orbitals to virtual B orbitals in borate anions, however, have to demand higher photon energies since the virtual B states stay above ∼3 eV. For example, the absorption of triangular [BO3]3- units fast becomes increasingly important if the wavelength is below 120 nm and reaches the maximum at ∼50 nm. In the short wavelength region below 80 nm, the tetrahedral [BO4]5- units even become one of the significant absorbers to ∼50 nm where the maximum absorption is achieved. To study the effect of optical polarization of BiB3O6, we define eq 5 following the Kramers-Kro¨nig relation. We call this quantity density of optical polarizations since it depicts the number of polarized states per volume at a given energy interval ω + dω.
PuV(ω) ) lim
δf0
{∫
ω-δ
0
ω′GuV(ω′) ω′2 - ω2
dω′ +
ω′G (ω′)
+∞ uV dω′ ∫w+δ 2 ω′ - ω2
}
(5)
As revealed in eq 5, the density of optical polarization should have a similar trend as the density of optical absorptions against the energy although Puv(ω) is explicitly calculated. Therefore [BiO4]5- units also dominate the optical polarization of BiB3O6 in the long wavelength region. The polarization of [BO3]3- and [BO4]5- units requires higher photon energies than that of [BiO4]5- units. These findings actually explain the calculated SHG coefficients (e.g., d22BiO4-only ) -2.829 pm/V, d22BO3-only ) -0.233 pm/V, and d22BO4-only ) -0.118 pm/V)8 in the static limit based on the wave function cutting approximation. The crystal BiB3O6 shows a quite different electronic origin for optical effects compared to other IA and IIA metalcontaining borates such as β-BaB2O4,35 LiB3O5, CsB3O5, and CsLiB6O10 36 where the optical effects are mainly accounted for by the contributions of [BO3]3- and [BO4]5- borate units. This special optical feature of BiB3O6 is, on one hand, due to the presence of the empty 6p orbital at Bi as an electron acceptor and, on the other hand, due to the high spatial-overlapping between the atomic orbital on the Bi and O atoms which results in the covalent Bi-O bonds and large optical transition matrix elements. However, the relatively difficult polarization of [BO4]5- units is the consequence of the 3-D network of
First-Principles Electronic Structure Study of BiB3O6
J. Phys. Chem. B, Vol. 110, No. 39, 2006 19261
Figure 7. The projected electronic densities in the O6-Bi-O7 and O4-Bi-O5 planes for Bi 6s-O 2p bondings. The contour lines are plotted between the minimum of 0.0 and the maximum of 0.1 with the step size of 0.001.
Figure 8. The density of optical absorptions for the total BiB3O6 (solid line), [BiO4]5- (dash line), [BO3]3- (dot line), and [BO4]5- (dash-dot line) at the scalar-relativistic B3PW level. The plots for the nonrelativistic B3PW and scalar-relativistic HF calculations can be found in Figures 1 and 2 in the Supporting Information.
tetrahedral structures which prohibits the electrons from being conjugated in a plane and thus higher photon energies are required. 3.5. Relativistic and Correlation Effects. We first examine the relativistic effect on the band gap of BiB3O6. The nonrelativistic B3PW calculation gives an indirect band gap of 4.92 eV. The scalar-relativistic effect enlarges the gap to 5.69 eV. One may note that the scalar-relativistic band gap is considerably larger than the experimental value of 4.3 eV,3b while the nonrelativistic band gap appears to agree quite well with only half of the deviation from the experimental value. In Table 3, it also seems that LDA and GGA approaches give band gaps closer to the experimental value than B3PW and B3LYP. We explain this result as the consequence of neglecting the spinorbit splitting of the Bi 6p states due to the averaged SO potential in the applied ECP of Bi, whereas SO effects originating from Bi have been theoretically found to significantly reduce the band gaps of some Bi-containing crystals such as CsBi4Te637 and Bi2Te3.38 The SO effects were therefore explored by recalculating the band gaps of BiB3O6 without and with SO correction, using the plane-wave code WIEN2K. The calculated band gaps without SO corrections are 3.90 eV (PW-LDA), 4.08 eV (PBE), and 4.08 eV (PW91) (Table 3), which are quite close to the corresponding results obtained with the CO-LCAO code CRYSTAL03. It can further be seen from Table 3 that SO interactions considerably lower the band gaps by ∼0.7 eV to 3.21 eV (PW-LDA), 3.37 eV (PBE), and 3.36 eV (PW91), which are left ∼1.0 eV lower than the experimental result and
thus are not sufficiently accurate. Although the plane-wave method cannot derive a SO correction for hybrid functionals, it is reasonable to infer that the reduction of the gap for B3PW due to SO couplings would be similar to that for LDA and GGA. On the other hand, at the all-electron state-averaged multiconfiguration Dirac-Hartree-Fock level with the Dirac-Coulomb Hamiltonian (GRASP),39 the Bi 6p1/2-6p3/2 atomic spinor splitting is 2.10 eV, i.e., 6p1/2 is lowered by 1.40 eV compared to the scalar relativistic 6p. We conclude that around 50% of the atomic SO effect is transferred to the solid BiB3O6. Although without SO corrections the B3PW functional overrates the band gap significantly, a SO lowering of the conduction band by this magnitude (ca. between 0.70 and 1.40 eV) would bring the B3PW estimate (ca. between 4.29 and 4.99 eV) into much better agreement with the experimental value of 4.3 eV than the LDA and GGA results. In Figure 4, there is an obvious state redistribution in the valence band due to the scalar-relativistic and correlation effects. A large part of Bi 6s states stay at ca. -5 eV for the nonrelativistic calculation. However, the scalar-relativistic effect strongly stabilizes the majority of Bi 6s states at ca. -16 eV and leaves only small contributions at the top of the valence band. It is mainly this negative shift of Bi 6s states that increases the band gap from 4.92 to 5.69 eV by scalar-relativistic effects. The relativistic stabilization of Bi 6s orbitals on one hand forces more electrons in Bi 6s orbitals (i.e., 6s2.17) than in the nonrelativistic case (i.e., 6s1.89), but on the other hand impairs the energy level matching between Bi 6s and O 2p orbitals. The net bonds in [BiO4]5- units are 6.80 and -5.29 for the nonrelativistic and scalar-relativistic calculations, respectively, which obviously indicates that by the scalar-relativistic effect the Bi-O bonds are partially destroyed in the [BiO4]5- units and stay antibonding. This bond instability mainly originates from the significantly intensified Bi 6s-O 2p antibonding interaction from -17.69 to -29.31, while the Bi 6s-O 2p and Bi 6p-O 2p bonding states are only slightly affected. It is this intensified antibonding repulsion between Bi 6p and O 2p orbitals that results in the elongations of the bonds Bi-O4 (O5) and Bi-O6 (O7) by 0.010 and 0.032 Å, respectively (see Table 4), compared to the nonrelativistic case, which is in contradiction to the expected relativistic effect that seems to always be a bond contraction effect in previous reports.40 BiB3O6 could not have come into being at ambient conditions if no dynamic electron-electron correlations have been taken into account, since the tetrahedral [BO4]5- units connecting the [BiO4]5- and [BO3]3- units would collapse. This is expected from the HF result given in Table 4 that the net antibonding in [BO4]5- is notably repulsive (-5.32). The reason for the occurrence of BiB3O6, from the chemical bonds point of view, is that the correlation effect fully removes the inner orbital
19262 J. Phys. Chem. B, Vol. 110, No. 39, 2006
Yang and Dolg
Figure 9. The total electronic densities in the (100) (left) and (001) (right) planes. The contour lines are plotted between the minimum of 0.0 and the maximum of 0.1 with the step size of 0.003.
TABLE 3: The Spin-Orbit (SO) Effect on the Band Gap of BiB3O6
c
band gap (eV)
HF
VWN-LDA
PW-LDA
PBE
no SO CRYSTAL03 ∆gap (eV)a
12.54 +8.24
3.96 -0.34
3.96 -0.34
4.05 -0.25
no SO WIEN2K ∆gap (eV)a
noc no
no no
3.90 -0.40
SO WIEN2K ∆gap (eV)a
no no
no no
3.21 (-0.69)b -1.09
PW91
B3LYP
B3PW
4.07 -0.23
5.79 +1.49
5.69 +1.39
4.08 -0.23
4.08 -0.22
no no
no no
3.37 (-0.71)b -0.93
3.36 (-0.72)b -0.94
no no
no no
a The deviations ∆ b gap compared with the experimental value (4.30 eV). The numbers in parentheses are the SO corrections calculated by WIEN2K. No HF, B3PW, and B3LYP calculations are available for the plane-wave code WIEN2K.
TABLE 4: The Scalar-Relativistic Effect on Bond Lengths of BiB3O6 bond lengths
scalar-relativistic
nonrelativistic
variations
2.388 2.109
2.378 2.077
0.010 0.032
1.405 1.367 1.344
1.409 1.365 1.342
-0.004 0.002 0.002
1.498 1.446
1.506 1.443
-0.008 0.003
[BiO4]5-
(Å) Bi-O4 (O5) Bi-O6 (O7) [BO3]3- (Å) B-O B-O B-O [BO4]5- (Å) B-O B-O
repulsion within [BO4]5- in the way that the B-O antibonding is weakened from -19.02 to -7.56 and the B-O bonding is intensified simultaneously from 13.70 to 23.94. In addition, the B-O bond intensity within [BO3]3- is further reinforced by the correlation effect from 100.33 to 168.12, which is another important factor that primarily contributes to the stability of BiB3O6. For [BiO4]5- units, in contrast to the scalar-relativistic effect, electron correlation brings an evident intensification on 6p-2p interactions rather than the intensification of 6s-2p for Bi-O bonds so as to partially counteract the instability in [BiO4].5Scalar-relativistic and correlation effects are important not only in the aspect of BiB3O6 stability, but also to its optical performance. If one compares the densities of states in Figure 4, the absence of either the scalar-relativistic effect or the correlation effect leads to the appearance of Bi 6s states at the high end of the valence band. In that case single Bi cations would exhibit a contribution to the optical response of BiB3O6 in the long wavelength region due to the inner-shell electronic transfer from 6s to 6p orbitals. 3.6. Asymmetric 6s Lone Pair Lobes at Bi Cations. If we take a look at the total electronic densities in the (100) and (001)
planes in Figure 9, at one side of Bi cations against oxygen atoms there is an asymmetric lobe-like distribution of the electronic density along the b lattice vector. The lobe is more spherical and expanded in the (100) plane than in the (001) plane, which implies that this electronic density is highly anisotropically distributed in 3-dimensional space. If one attributes the origin of the lobe electronic densities only to Bi 6s orbitals one cannot explain its asymmetry since a 6s orbital is shaped spherical. Some nonspherical orbitals must be hybridized with 6s orbitals. Since there are evident covalent bonds between Bi and O atoms, i.e., 6s-2p and 6p-2p covalent interactions as concluded above, we assume that some of those nonspherical orbitals are 2p orbitals at nearby O. Nevertheless, the detailed mechanism of the lone-pair formation at Bi in BiB3O6 will be systematically investigated in the future. 4. Conclusions Systematic first-principles electronic structure calculations (HF, LDA, GGA, and hybrid density functionals) for the geometry and the band gap of the monoclinic crystal BiB3O6 within the CO-LCAO scheme using Gaussian basis sets revealed that the hybrid functionals, especially B3PW, yield the best results. Estimates for the band gap of BiB3O6 (4.29-4.99 eV) proved to be only reliable if spin-orbit effects are accounted for. Two important consequences are due to the covalency between Bi and O as studied with a first-principles COOP analysis. First Bi 6p-O 2p bonding interactions promote its spatial-overlapping and electronic densities in the middle of the Bi-O bond, which accordingly favors the O 2p-Bi 6p electronic transitions and leads to the dominant contribution from [BiO4]5- units to optical effects of BiB3O6 in the long wavelength region. Second, the Bi-O covalent interactions bring more asymmetric O 2p orbitals into originally symmetric Bi 6s orbitals and thus partly explain the formation of the Bi
First-Principles Electronic Structure Study of BiB3O6 lone-pair lobe. The single Bi cations according to our calculations do not play a critical role in determining the optical properties of BiB3O6 in contrast to explanations given in some experimental reports. BiB3O6 thus shows quite a different electronic origin for optical effects from other IA and IIA metalcontaining borates such as β-BaB2O4, LiB3O5, CsB3O5, and CsLiB6O10 where the [BO3]3- and [BO4]5- units contribute most. The scalar-relativistic effect results in a Bi-O net antibonding repulsion and the Bi-O bond lengths are therefore elongated opposite to the typical relativistic bond contraction effect. The occurrence of BiB3O6 at ambient conditions is due to the correlation effect that significantly stabilizes tetrahedral [BO4]5- units from unstable B-O net antibondings to stable B-O net bondings. The covalent B-O bonds in [BO3]3- and Bi-O bonds in [BiO4]5- are enhanced as well due to the correlation effect. Acknowledgment. This work was financially supported by the project Graduiertenkolleg 549, “Azentrische Kristalle” at University of Ko¨ln. The authors are grateful to Prof. Dr. L. Bohaty´ for valuable discussions on BiB3O6. We also thank the Computer Center at University of Cologne (RRZK) for providing us access to the high performance computing CLIO SunOpteron cluster. Supporting Information Available: Additional figures for reference of scalar-relativistic and correlation effects. This material is available free of charge via the Internet at http:// pubs.acs.org. References and Notes (1) (a) Bohaty´, L.; Liebertz, J. Z. Kristallogr. 1985, 170, 18. (b) Becker, P.; Liebertz, J.; Bohaty´, L. J. Cryst. Growth 1999, 203, 149. (c) Becker, P.; Wickleder, C. Cryst. Res. Technol. 2001, 36, 27. (d) Teng, B.; Wang, J. Y.; Wang, Z. P.; Hu, X. B.; Jiang, H. D.; Liu, H.; Cheng, X. F.; Dong, S. M.; Liu, Y. G.; Shao, Z. S. J. Cryst. Growth 2002, 240, 331. (2) (a) Becker, P.; Wickleder, C. Cryst. Res. Technol. 2001, 36, 27. (b) Teng, B.; Wang, Z. P.; Jiang, H. D.; Cheng, X. F.; Liu, H.; Hu, X. B.; Dong, S. M.; Wang, J. Y.; Shao, Z. S. J. Appl. Phys. 2002, 91, 3618. (c) Haussu¨hl, S.; Bohaty´, L.; Becker, P. Appl. Phys. A 2006, 82, 495. (3) (a) Hellwig, H.; Liebertz, J.; Bohaty´, L. Solid State Commun. 1998, 109, 249. (b) Hellwig, H.; Liebertz, J.; Bohaty´, L. J. Appl. Phys. 2000, 88, 240. (c) Du, C.; Wang, Z.; Liu, J.; Xu, X.; Teng, B.; Fu, K.; Wang, J.; Liu, Y.; Shao, Z. Appl. Phys. B 2001, 73, 215. (d) Wang, Z. P.; Teng, B.; Fu, K.; Xu, X. G.; Song, R. B.; Du, C. L.; Jiang, H. D.; Wang, J. Y.; Liu, Y. G.; Shao, Z. S. Opt. Commun. 2002, 202, 217. (e) Kityk, I. V.; Majchrowski, A. Opt. Mater. 2004, 25, 33. (f) Gossling, A.; Moller, T.; Stein, W. D.; Becker, P.; Bohaty´, L.; Gru¨ninger, M. Phys. Status Solid B 2005, 242, R85. (4) Xia, H. R.; Li, L. X.; Teng, B.; Zheng, W. Q.; Lu, G. W.; Jiang, H. D.; Wang, J. Y. J. Raman Spectrosc. 2002, 33, 278. (5) Ghotbi, M.; Ebrahim-Zadeh, M. Opt. Express 2004, 12, 6002. (6) Xue, D.; Betzler, K.; Hesse, H.; Lammers, D. Phys. Status Sol. A 1999, 176, R1. (7) Xue, D.; Betzler, K.; Hesse, H.; Lammers, D. Solid State Commun. 2000, 114, 21. (8) Lin, Z. S.; Wang Z. Z.; Chen C. T.; Lee M. H. J. Appl. Phys. 2001, 90, 5585.
J. Phys. Chem. B, Vol. 110, No. 39, 2006 19263 (9) Burns, P. C.; Grice, J. D.; Hawthorne, F. C. Can. Mineral. 1995, 33, 1131. (10) Becker, P. Z. Kristallogr. 2001, 216, 523. (11) Parthe, E. Z. Kristallogr. 2002, 217, 179. (12) Barber, M. J.; Becker, P. Z. Kristallogr. 2002, 217, 205. (13) Koch, W.; Holthausen, M. C. A Chemist’s Guide to Density Functional Theory, 2nd ed.; Wiley-VCH Verlag GmbH: Weinheim, Germany, 2001. (14) Pisani, C.; Busso, M.; Capecchi, G.; Casassa, S.; Dovesi R.; Maschio, L.; Zicovich W. C.; Schutz, M. J. Chem. Phys. 2005, 122, 094113. (15) Shukla, A.; Dolg, M.; Fulde, P.; Stoll, H. Phys. ReV. B 1999, 60, 5211. (16) Stoll, H.; Paulus, B.; Fulde, P. J. Chem. Phys. 2005, 123, 144108. (17) Mukhopadhyay, A. B.; Dolg, M.; Oligschleger, C. J. Chem. Phys. 2004, 120, 8734. (18) Saunders: V. R.; Dovesi, R.; Roetti, C.; Orlando, R.; Wilson, C. M. Z.; Harrison, N. M.; Doll, K.; Civalleri, B.; Bush, I. J.; Arco, P. D.; Llunell, M. CRYSTAL2003 1.0 User’s Manual, 2003, August 7. (19) Schlegel, H. B. J. Comput. Chem. 1982, 3, 214. (20) Monkhorst, H. J.; Pack, J. D. Phys. ReV. B 1976, 13, 5188. (21) Becke, A. D. J. Chem. Phys. 1988, 88, 1053. (22) Ku¨chle, W.; Dolg, D.; Stoll, H.; Preuss H. Mol. Phys. 1991, 74, 1245. (23) Energy-Consistent Pseudopotentials of the Stuttgart/Cologne Group: http://www.theochem.uni-stuttgart.de/pseudopotentiale/clickpse.en.html. (24) EMSL Gaussian Basis Set Order Form: http://www.emsl.pnl.gov/ forms/basisform.html. (25) Dunning, T. H. J. Chem. Phys. 1971, 55, 716. (26) Blaha, P.; Schwarz, K.; Madsen, G. K. H.; Kvasnicka, D.; Luitz, J. WIEN2k An Augmented Plane-waVe + Local Orbitals Program for Calculating Crystal Properties (Karlheinz Schwarz, Techn. Universita¨t Weinheim), 2001, ISBN 3-9501031-1-2. (27) Desclaux J. P. Comput. Phys. Comm. 1975, 9, 31. (28) Ko¨lling, D. D.; Harmon, B. J. Phys. C 1980, 13, 6147. (29) Yang, J.; Dolg, M. Theor. Chem. Acc. 2005, 113, 212. (30) Wang, Y. X.; Dolg, M. Theor. Chem. Acc. 1998, 100, 124. (31) Hoffmann, R. Angew. Chem., Int. Ed. Engl. 1987, 26, 846. (32) Hughbanks, T.; Hoffmann, R. J. Am. Chem. Soc. 1983, 105, 3528. (33) Dronskowski, R.; Blo¨chl, P. E. J. Phys. Chem. 1993, 97, 8617. (34) Grechnev, A.; Ahuja R.; Eriksson, O. J. Phys.: Condens. Matter 2003, 15, 7751. (35) Lin, J.; Lee, M. H.; Liu, Z. P.; Chen, C. T.; Pickard, C. J. Phys. ReV. B 1999, 60, 13380. (36) Lin, Z. S.; Wang, Z. Z.; Chen, C. T.; Lee, M. H. Phys. ReV. B 2000, 62, 1757. (37) Larson, P.; Mahanti, S. D.; Chung, D. Y.; Kanatzidis, M. G. Phys. ReV. B 2002, 65, 045205. (38) Larson, P.; Mahanti, S. D.; Kanatzidis, M. G. Phys. ReV. B 2000, 61, 8162. (39) Dyall, K. G.; Grant, I. P.; Johnson, C. T.; Parpia, F. A.; Plummer, E. P. Comput. Phys. Commun. 1989, 55, 425. GRASP is a relativistic atomic electronic structure code. (40) Pyykko¨, P. Chem. ReV. 1988, 88, 563. (41) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200. (42) Dirac, P. A. M. Proc. Cambridge Philos. Soc. 1930, 26, 376. (43) Perdew, J. P. Phys. ReV. B 1989, 40, 3399. (44) Perdew, J. P.; Wang, Y. Phys. ReV. B 1992, 45, 13244. (45) Perdew, J. P. Electronic structure of solids; Akademie Verlag: Berlin, Germany, 1991. (46) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865. (47) Perdew, J. P.; Yue, W. Phys. ReV. B 1986, 33, 8800. (48) Lee, C.; Yang, W.; Parr, R. G. Phys. ReV. B 1988, 37, 785.