The Relevance of Dispersion Interactions for the Stability of Oxide

Dec 8, 2010 - An implementation of artificial neural-network potentials for atomistic materials simulations: Performance for TiO 2. Nongnuch Artrith ,...
0 downloads 10 Views 143KB Size
22718

J. Phys. Chem. C 2010, 114, 22718–22726

The Relevance of Dispersion Interactions for the Stability of Oxide Phases Jose´ C. Conesa* Instituto de Cata´lisis y Petroleoquı´mica, CSIC, Marie Curie 2, Cantoblanco, 28049 Madrid, Spain ReceiVed: September 23, 2010; ReVised Manuscript ReceiVed: October 28, 2010

The total energies of TiO2 and Al2O3 allotropic forms computed using density functional theory (DFT) methods that include dispersion interaction effects are compared. For TiO2, adding energy terms of r-6 form with coefficients derived from atomic polarizabilities leads to the correct result that rutile is more stable than brookite, anatase, and other forms, while previous DFT studies without this correction wrongly predicted rutile to be less stable. The magnitude of the correction is significant because of the high polarizability of the oxide ion, and produces energy differences between the titania phases that are reasonably close to available experimental values. The van der Waals density functional does not yield the correct result, but the error is significantly decreased. For Al2O3 the experimental energy difference between R and γ forms is also approached better when including dispersion corrections in DFT calculations. The results show that dispersion interaction should not be ignored when computing with DFT energy differences between oxide structures and that for the latter the dispersion coefficients appropriate for neutral species should not be used. It is proposed that the correct reproduction of these differences be included in the benchmarks used for testing van der Waals-type functionals. Introduction Density functional theory (DFT) is used today extensively with great success to compute many aspects of the atomic and electronic structure of solids and molecules, including their electromagnetic and optical properties, their energies, and their vibrational characteristics and chemical reactivity. In the prediction and understanding of these properties the approximations inherent in DFT lead however to inaccuracies; sometimes the results are systematically biased and in other cases they vary more randomly. Typical deficiencies in DFT are, for example, the underestimation of the separation between filled and empty electron levels (the gap problem, which may even lead to the prediction of a metallic character for some semiconducting solids) and the imperfect cancellation of electron self-interaction (which leads to excessive electronic delocalization and frequently wrong magnetic properties). These two problems are related to the inadequate representation of the electron exchange term in the functional, and several “beyond DFT” treatments have been devised trying to correct these deficiencies. This work is concerned however with the possible effect of another wellknown deficiency of DFT, i.e., its lack of a suitable treatment of dispersion interactions (hereinafter called DI), also known as London interactions, which are related to correlated fluctuations of electron densities in different parts of the same system. These are responsible for the weak van der Waals forces between nonpolar or scarcely polar molecules, between these and solid surfaces or between the layers of certain lamellar compounds such as graphite or some sulfides. The energies involved in these interactions usually follow, at least to the first order, an r-6 dependence on the distance r between the molecules that interact as a result of this effect. Standard DFT, with only local or semilocal terms in the correlation functional, is unable to account correctly for these fully nonlocal effects. * To whom correspondence should be addressed. E-mail: jcconesa@ icp.csic.es.

Some first-principles treatments of these interactions within the DFT framework have been attempted.1-4 Particular attention is paid currently to the so-called van der Waals density functional (vdW-DF) method,1 which uses a fully nonlocal functional combined with the correlation LDA functional and a GGA exchange functional (revPBE in its initial formulation) and has already been applied to a good number of cases.5 This type of method is currently under discussion and development,6,7 and there is ongoing work to determine which exchange correlation functional should be used with it.8 On the other hand, it has been indicated that completely omitting correlation functional gradient corrections, as formulated in the original vdW-DF work,1 might not be the best option.9 Thus the details of the vdW-DF method are still under discussion. While these first-principles methods are not yet consolidated, several authors have proposed taking DI into account in DFT calculations by adding to the self-consistently determined electronic energy ESCF a sum of terms of r-6 type (plus other r-n terms, with n g 8, in some cases) that reflect the DI existing between the constituent atoms. This type of approach has been followed in recent works by Grimme10,11 who proposes to add such a term as follows

E ) ESCF - s0

∑ f(rij)Cijrij-6

(1)

i>j

Here i and j are indices running over all atoms in the system and rij is the distance between atoms i and j. s0 is a scaling parameter that depends on the functional used to obtain ESCF and reflects the fact that typical GGA functionals may already include some long-range attractive component; in most cases s0 < 1. f is a damping function that screens (decreases) the interaction between atoms that are rather close and therefore begin to experience Pauli repulsion. In some works the Cij coefficients of these terms are derived for each particular case through a more or less accurate first-

10.1021/jp109105g  2010 American Chemical Society Published on Web 12/08/2010

Relevance of Dispersion Interactions

J. Phys. Chem. C, Vol. 114, No. 51, 2010 22719

principles approach,12 but in Grimme’s original method (or its newer version13) the coefficient Cii corresponding to the interaction between atoms of the same element, as well as the radii Ri entering the damping function, are taken from properties of neutral atomic species and tabulated. The Cij coefficients for the interaction between dissimilar atoms are also computed as the geometric mean of the corresponding homoatomic coefficients:

Cij ) (CiiCjj)1/2

(2)

Grimme’s original method is being used at present by a number of authors and is being incorporated to widely available DFT codes such as Quantum Espresso, CRYSTAL09, VASP (in a modification of the current 5.2 version, not yet released publicly at the time of writing this work), or fhi-aims (although in this code Cii and Ri values are adjusted as a function of charge redistribution). The above-mentioned more recent version of Grimme’s method13 has not yet been incorporated into these codes. All of these approaches to evaluate the DI energy (DFT-D or vdW-DF type) have been used until now mainly when modeling systems where the interaction of interest is a weak one (the van der Waals interaction), which affects neutral nonpolar or lightly polar molecules and/or some solid surfaces. These approaches are normally not used when modeling inorganic bulk solids or chemically reactive surfaces with quantum methods. This is surprising however, because dispersion terms of the r-n type (damped or not) have been very generally used for decades in the methods that model solids and surfaces with interatomic potentials, especially in the case of ionic oxides in which the oxide ions are known to be highly polarizable and thus subjected to non-negligible correlated fluctuations in the electronic density.14 The work presented here tries to evaluate, at least in a semiquantitative way, how large the effect of DI on the energetics of oxides can be, paying particular attention to the relative stabilities and structural characteristics of different polymorphs and taking TiO2 and Al2O3, prototypical oxides of high technological interest, as case examples. These have several crystalline forms with special applications in catalysis and photocatalysis, and for these some data on their relative stabilities have been determined experimentally.15-17 It turns out that DI may have a significant effect on the energy differences between them. Theoretical Calculation Methods In the present calculations the DI will be mostly modeled, as in Grimme’s original work, by eq 1 (this will henceforth be called the DFT-D method). However, it is unlikely that the Cii coefficients given in that work, based on the neutral state of the elements, are suitable for oxides such as TiO2 or Al2O3. Indeed these coefficients are known to be directly related to the (high frequency) polarizabilities of the atomic species involved, which in turn are different for different charge states of the latter. For the present study, in particular, the Cii values tabulated in Grimme’s original paper10 are not to be relied upon. In fact the latter gives Cii coefficients of 10.8 and 0.7 J nm6/ mol for Ti and O, respectively, while it is clear that the polarizability of the oxide ion is much higher than that of the Ti cation. Cij coefficients can be computed approximately from the polarizabilities of the species involved using the well-known Slater-Kirkwood equation, which can be expressed as

Cij ) K

RiRj 1/2

(Ri/Ni)

+ (Rj/Nj)1/2

(3)

Ri is the high-frequency polarizability of atomic species i, and Ni represents the effective number of electrons contributing to that polarizability. If R values are given in Å3 units, eq 3 with K ) 15.7 gives Cij in eV Å6 units.18 For the interaction between identical atoms (or ions) this is reduced to

Cii ) (K/2)(Ri3Ni)1/2

(4)

On the other hand, the polarizabilities of the individual atoms (ions) forming a compound are related to the refractive index n of the latter through the Lorentz-Lorentz relationship19

RM ) (3VM/4π)

n2 - 1 ) n2 + 2

∑ miRi

(5)

i

VM being the volume per formula unit and RM the polarizability per formula unit. As reflected in the right side of eq 5, polarizabilities are assumed to be additive,20 i.e., RM would be equal to the sum of the atomic (or ionic) polarizabilities Ri the latter being weighted by the number mi of i-type species in the formula. As regards the Ni value, an approximation to it in terms of the number of external and internal orbital electrons Next and Nint has been proposed by Cambi et al.18

N/Next ) 1 + (1 - Next/Nint)(Nint/Ntot)2

(6)

where Ntot ) Next + Nint. Here the 2sp electrons will be considered as the external ones for O2- and Al3+ ions, while for Ti4+ ions the 3sp electrons will be considered external. Equations 3 and 4 lead to the combination formula

Cij )

2[(CiiCjj)2NiNj]1/3 (CjjNii2)1/3 + (CiiNjj2)1/3

(7)

It is easy to see that Grimme’s geometric mean (eq 2) predicts values for Cij that are greater than those in eq 7 by a factor τ ) (x + x-1)/2 where

x ) [(Cii/Cjj)(Nj2/Ni2)]1/6

(8)

In this work all heteroatomic Cij coefficients will be obtained through eq 7, while eq 4 will be used to estimate homoatomic Cii coefficients using atomic polarizabilities determined from experimental refractive indices (where these are available) through eqs 6, 3, and 5. Since the latter formula gives only the total molecular polarizability RM directly, it is necessary to choose some rule for decomposing it into the Ri values. For this, it can be assumed20,21 that the polarizability of cations of the type considered here is hardly affected by the environment. For Ti4+ and Al3+ the values included in the compilation of Dimitrov et al.,20 obtained using the ionic refraction data of Kordes,22 will be used here. The average anion polarizability can then be directly derived for each oxide phase from the index of refraction of the latter. This index can be anisotropic, but

22720

J. Phys. Chem. C, Vol. 114, No. 51, 2010

Conesa

only its average value is used in this formulation. This is an approximation; actually, the ion polarizability itself can be anisotropic,23 while the equations above assume that Ri and Cij are isotropic. These anisotropies will be ignored here. For the damping function f the expression first proposed by Wu and Yang,24 used as well by Grimme, will be assumed

{

[ (

f(rij) ) 1 + exp -d

)]}

rij -1 R _

-1

(9)

where R ) (Ri + Rj) is taken to be the sum of the van der Waals radii of the ions; parameter d (having typically a value d ) 20) governs the steepness of the damping function. It is worth noting that this damping function is much steeper than that used in the Tang-Toennies potential25 frequently used in interatomic potentials-based modeling methods

TABLE 1: Parameters Used for Computing the DI Correction parameter 〈n〉 (averaged from the anisotropic values) RM (Å3)a RO (Å3)b NOc NMc RO (Å)d RM (Å)d COO (eV Å6)e CMM (eV Å6)e a

6

f (rij) ) 1 - exp(-y) TT

∑ y /k! k

(10)

for TiO2 (AN, RU, and BR phases)

for Al2O3 (R phase)

2.677 (RU) 2.622 (BR) 2.537 (AN) 0.184 2.415 (RU) 2.44 (BR) 2.527 (AN) 6.7 8.5 1.36 0.605 86.1 (for RU) 87.7 (for BR) 92.2 (for AN) 2.035

1.765

b

0.054 1.365 6.7 6.7 1.38 0.675 36.6 0.29 c

Reference 20. Computed through eq 5. Obtained through eq 6. d Shannon’s ionic radii27 appropriate for the corresponding coordination numbers. e Computed through eq 4.

k)0

where y ) rij/R (with R defined in some specific way) and that other damping functions have been proposed as well.26 For the Ri parameter, used in the damping function of eq 9, the standard Shannon ionic radii27 have been adopted here. The calculations were carried out with the plane-wave program VASP,28,29 widely used for modeling periodic solids, modified to include r-6 terms having the form given by eq 1. The forces analytically derived from these terms were also taken into account when relaxing both atomic positions and unit cell dimensions. The PBE functional has been used, and correspondingly in eq 1 a scaling parameter s0 ) 0.75 has been assumed.11 The core region was represented using the PAW method,28c,30 with the standard options of the VASP program and keeping the 2s, 3s, and 3p electron shells (and higher) for O, Al, and Ti atoms, respectively, in the valence state. The planewave expansion cutoff used (500 eV) and the density of the reciprocal space sampling (with Γ-centered Monkhorst-Pack grids, having, for example, 10 × 10 × 14, 14 × 14 × 14, and 6 × 8 × 8 points for rutile, anatase, and brookite, respectively) were verified to be high enough for the energy differences relevant for this work to converge within less than 2 meV per MOx formula unit. Some test calculations were made using for the damping function formula the radii Ri resulting from the Hirshfeld-type charge analysis proposed by Tkatchenko et al.12c To obtain these Ri values, the all-electron DFT program fhi-aims,31,32 which expands the wave functions in a basis of atom-centered numeric functions, was used. The basis adopted was of the “tier 2” class for both anions and cations, which ensures a good convergence in the wave function representation. In some other tests on TiO2 an alternative evaluation of the DI energy has been carried out using the fully nonlocal vdWDF functional1a through the implementation of the latter, based in a Fourier transform integration method33 and reported by Thonhauser et al.,1b in the DFT planewave-based program Quantum-Espresso.34,35 The PBE functional (in the exchange and local correlation parts) was used, and the ion cores were represented by ultrasoft pseudopotentials taken from the QuantumEspresso web page. These pseudopotentials kept in the valence space the 3s and 2s electron shells (and higher) for Ti and O, respectively. In the currently available form the code allows the automatic relaxation of atomic coordinates but not cell

TABLE 2: Crystal Cell Parametersa and Energy Differences ∆Eb in Respect to Rutile, Computed for the Three Main TiO2 Phases and Compared to Experimental Data crystal form

DFT results

DFT-D results

exptl datac

RU

a ) 4.654 c ) 2.974 VF ) 32.20 c/a ) 0.639 ∠(O-Ti-O) ) 98.36°

a ) 4.654 c ) 2.955 VF ) 32.0 c/a ) 0.653 ∠(O-Ti-O) ) 97.45°

a ) 4.587 c ) 2.954 VF ) 31.07 c/a ) 0.644 ∠(O-Ti-O) ) 98.76°

BR

a ) 9.269 b ) 5.512 c ) 5.184 VF ) 33.11 c/a ) 0.559 c/b ) 0.940 ∆E ) -4.8

a ) 9.214 b ) 5.528 c ) 5.123 VF ) 32.62 c/a ) 0.556 c/b ) 0.927 ∆E ) +0.6

a ) 9.180 b ) 5.430 c ) 5.164 VF ) 32.18 c/a ) 0.562 c/b ) 0.951 ∆E ) +0.716

AN

a ) 3.821 a ) 3.874 a ) 3.782 c ) 9.639 c ) 9.084 c ) 9.502 VF ) 35.18 VF ) 34.09 VF ) 33.98 c/a ) 2.523 c/a ) 2.345 c/a ) 2.514 ∠(Ti-O-Ti) ) 156.1° ∠(Ti-O-Ti) ) 163.1° ∠(Ti-O-Ti) ) 156.3° ∆E ) -8.8 ∆E ) +4.25 ∆E ) +1.715

a Vector lengths are in Ångstroms; volumes VF are per TiO2 formula unit in Å3. b kJ/mol. c Structural data are at 15 K for RU and AN,36 at ambient temperature for BR.

dimensions; therefore the relaxation of these latter was carried out manually, the final cell shape being obtained through a polynomial fit to energy values obtained for a set of cell dimensions around the energy minimum. Results and Discussion Main TiO2 Phases. Calculations based on either standard DFT or the DFT-D method were made first for the three main phases of TiO2 (rutile, brookite, and anatase, here labeled as RU, BR, and AN), using the parameter values given in Table 1 for the second method. It can be noted that, with the latter, the aforementioned τ factor through which Grimme’s geometric mean (eq 2) would overestimate the Cij values of eq 7 would amount to ca. 1.26. The structural parameters obtained after geometry optimization are summarized in Table 2, together with experimental data available for these phases. For anatase and rutile experimental values taken at temperatures close to 0 K36 are used in the comparison, considering that the calculations

Relevance of Dispersion Interactions include no thermal effects. For brookite no such low temperature data were available; therefore ambient temperature data are used.37 As expected, when the attractive DI is included the crystal lattices contract somewhat, as can be seen in the values of the volume per formula unit VF. This effect is small for rutile, which is already a rather compact oxide, and is higher for brookite and anatase, particularly in the latter, for which a value rather closer to the experimental value is reached; this occurs however mostly through a contraction in the c axis (connected to an increase in the Ti-O-Ti angle), which distorts the cell shape somewhat. This volume contraction is of a magnitude similar to that which would be produced, in the case of anatase, by an added external isotropic pressure of 60 kBar (for rutile the equivalent pressure is rather smaller, ca. 10 kBar), as was proven by a separate standard DFT calculation including this pressure. The resulting cell shape distortion obtained in the latter calculation is however different (data not shown); probably, the directionality of the DI forces is not well reproduced with the present DFT-D approach. Indeed the latter is based on assuming scalar polarizability values for the ions, while as said above these are known to be anisotropic,23 which is not surprising taking into account that the coordination of oxygen in these lattices is triangular planar. More importantly, the DI effects modify the total energies by different amounts for the different phases studied. The difference is due not only to the different lattice topologies but also to the different interionic distances (these tend to be smaller in denser forms, leading to a higher DI stabilization, although there can be also an effect on the damping function f) and to the differences in oxide ion polarizabilities and, hence, in COO coefficients (which are seen to be higher in less dense forms). This is reflected in the excess energy, as compared to rutile, of the other phases; this excess energy is called ∆E here. The net consequence of these differences is that, as shown in Table 2, while standard DFT wrongly predicts that anatase and brookite are (in bulk phase) more stable than rutile, i.e., that their ∆E values are negative (a well-known deficiency of DFT, even when using hybrid functionals38-40), the DI-corrected values recover the correct order of stability, leading to positive ∆E values. This effect of increasing the (computed) stability of rutile vs the other phases can probably be traced to its higher density, which increases the magnitude of the r-6 terms. Note that its lower COO value (due to a lower RO) works in the opposite direction; it obviously only partially compensates this effect. The variations in ∆E (that will be called δ∆E here) produced by including the DI-related correction (eq 1) amount thus to +13.45 and +5.4 kJ/mol for anatase and rutile, respectively. These values can be compared with the variations that would equate the standard DFT ∆E values with the experimental ones: +10.45 and +5.2 kJ/mol. They are clearly of a very similar magnitude. Although eq 1 is just a crude way of including DI effects (so that the coincidence between the experimental and computed ∆E within 0.1 kJ/mol found for brookite is surely fortuitous) and is affected by a number of approximations and uncertainties (e.g., ignoring the r-n terms with n > 6, method of estimating Cij, value of the s0 factor, shape and parameters of the damping function), this similarity is noteworthy and strongly suggests that the inadequacy of standard DFT to reproduce correctly the relative stabilities of these phases could arise, at least to a significant extent, from its inability to account properly for DI.

J. Phys. Chem. C, Vol. 114, No. 51, 2010 22721 If so, the correction computed through eq 1 is, for anatase, somewhat larger than it should be, since the resulting ∆E (having now the correct sign) is significantly larger than the experimentally measured value, although for brookite it is much closer. One could obtain somewhat lower ∆E values, for example, by using a smaller s0 value: indeed test calculations using s0 ) 0.65 instead of 0.75 lead to a closer prediction for anatase (∆E ) +3.1 kJ/mol) while still obtaining a positive value (∆E ) +0.2 kJ/mol) for brookite, although the latter is lower than the experimental one. Different radii Ri could also be used in the damping formula. In this respect a calculation of Hirshfeld partition-corrected radii, as proposed by Tkatchenko et al.,12c was undertaken using the fhi-aims31 program, but rather large values, probably unrealistic, were then obtained (e.g., RTi ) 2.275 Å, RO ) 1.65 Å for anatase), most probably due to the fact that this approach is based on adjusting the atomic radii of neutral species. This led to a much higher damping, and correspondingly very small δ∆E values were obtained in the resulting calculation with the said code. These approaches, looking for a better agreement with experiment through ad hoc changes in the parameters of the method, seem unjustified and were not retained for the rest of this work, given the approximations involved. Another test was carried out using the same Cij and Ri parameters given in the original Grimme’s method. It was found that anatase was correctly predicted to be less stable than rutile by ∆E ) +2.0 kJ/mol, but brookite was predicted to be more stable than rutile by ∆E ) -1.5 kJ/mol. This recovery of the correct ∆E sign for anatase but not for brookite when using Grimme’s parameters seems to have been found by some other author.41 The volumes per formula unit also decrease (VF ) 31.55, 32.61, and 34.87 Å3 for RU, BR, and AN, respectively) without reaching the experimental ones. Thus even if these coefficients are probably not to be relied upon, as already said, because they correspond to the neutral atoms and not to the ionic species, they provide some adjustment to the ∆E and VF values in the correct direction; still, they do not lead to full qualitative agreement with experiment. As an additional test, an alternative estimation of the DI effect was carried out also for anatase and rutile using the vdW-DF functional. In this case the standard DFT calculation, using ultrasoft pseudopotentials, led to ∆E ) -10.1 kJ/mol, somewhat higher than in the PAW-based calculation done with VASP. When using vdW-DF the calculation predicted, with both atomic positions and lattice dimensions relaxed, ∆E ) -6.95 kJ/mol; thus it falls short of recovering the experimentally found higher energy of anatase. It may be noted here that the cell volumes relaxed with vdW-DF, amounting to 35.81 and 32.57 Å3 for anatase and rutile, respectively, are higher than those obtained with the same code and pseudopotentials using standard DFT (34.49 and 31.54 Å3), while one would expect that the inclusion of DI should contract this volume as said before. It is therefore unclear whether including the DI with the vdW-DF functional is suitable to model the interaction in this type of highly polar solid compounds properly. Note, in this respect, that as mentioned above there is at present some discussion concerning which exchange functional should be used with the vdW-DF method and whether some semilocal GGA contribution to the correlation should be kept. Thus it might be that, by modifying the rest of the functional, better results would be obtained, but this is not possible with the code available here. Still, the δ∆E ) +3.25 kJ/mol value obtained in the present vdW-DF calculations is again larger than the experimental ∆E value. Thus this approach, although not proving that including DI is the

22722

J. Phys. Chem. C, Vol. 114, No. 51, 2010

Conesa TABLE 3: Crystal Cell Parametersa and Energy Differences ∆Eb in Respect to Rutile, Computed for Other TiO2 Phases, and Compared to Experimental Data crystal form

DFT-D results

B42

a ) 12.300 b ) 3.763 c ) 6.626 β ) 106.95° VF ) 36.67 c/a ) 0.539 c/b ) 1.761 ∆E ) -10.6

a ) 12.154 b ) 3.815 c ) 6.556 β ) 106.86° VF ) 36.36 c/a ) 0.539 c/b ) 1.719 ∆E ) +3.05

a ) 12.163 b ) 3.735 c ) 6.513 β ) 107.29° VF ) 35.31 c/a ) 0.535 c/b ) 1.744

II47

a ) 4.583 b ) 5.577 c ) 4.938 VF ) 31.55 c/a ) 1.077 c/b ) 0.885 ∆E ) -3.0

a ) 4.552 b ) 5.574 c ) 4.939 VF ) 31.33 c/a ) 1.085 c/b ) 0.886 ∆E ) +1.4

a ) 4.532 b ) 5.502 c ) 4.906 VF ) 30.58 c/a ) 1.083 c/b ) 0.892

R44

a ) 9.566 b ) 2.983 c ) 4.951 VF ) 35.32 c/a ) 0.518 c/b ) 1.660 ∆E ) +6.7

a ) 9.221 b ) 2.949 c ) 4.828 VF ) 32.83 c/a ) 0.524 c/b ) 1.637 ∆E ) +8.4

a ) 9.459 b ) 2.963 c ) 4.902 VF ) 34.30 c/a ) 0.518 c/b ) 1.657

H45

a ) 10.296 c ) 2.982 VF ) 39.52 c/a ) 0.290 ∆E ) +6.9

a ) 10.180 c ) 2.950 VF ) 38.21 c/a ) 0.290 ∆E ) +11.0

a ) 10.164 c ) 2.963 VF ) 38.26 c/a ) 0.292

BD52

a ) 4.845 b ) 4.925 c ) 5.065 β ) 100.02° VF ) 29.76 c/a ) 1.045 c/b ) 1.028 ∆E ) +8.6

a ) 4.806 b ) 4.931 c ) 5.061 β ) 98.58° VF ) 29.65 c/a ) 1.053 c/b ) 1.026 ∆E ) +46.7

a ) 4.606 b ) 4.896 c ) 4.933 β ) 99.17° VF ) 27.46 c/a ) 1.071 c/b ) 1.008

OI50

a ) 9.334 b ) 4.947 c ) 4.785 VF ) 27.62 c/a ) 0.513 c/b ) 0.967 ∆E ) +22.8

a ) 9.303 b ) 4.981 c ) 4.771 VF ) 27.64 c/a ) 0.513 c/b ) 0.958 ∆E ) +52.8

a ) 9.052 b ) 4.836 c ) 4.617 VF ) 25.26 (VF0 ) 27.23) c/a ) 0.510 c/b ) 0.955

CT51

a ) 5.213 b ) 3.197 c ) 6.343 VF ) 26.43 c/a ) 1.217 c/b ) 1.984 ∆E ) +73.5

a ) 5.284 b ) 3.113 c ) 6.274 VF ) 25.8 c/a ) 1.187 c/b ) 2.016 ∆E ) +130

a ) 5.163 b ) 2.989 c ) 5.966 VF ) 23.02 (VF0 )26.24) c/a ) 1.156 c/b ) 1.996

Figure 1. Correlation between the oxide ion polarizability computed (from the mean refractive index) for the three main TiO2 phases and their volume per TiO2 formula unit. The straight line is a linear fit.

solution to obtaining the correct stability order of titania phases, does support the assumption that this interaction is relevant for determining the aforementioned order. Other TiO2 Phases. To see what the effect of DI could be on other known titania phases DFT and DFT-D calculations were carried out for them following the same method that leads to the results in Table 2. The set of phases studied includes several which have 6-coordinated Ti: the TiO2(B) phase, first known through synthesis42 and later also detected in mineral specimens;43 TiO2(R) and TiO2(H) with ramsdellite and hollandite/priderite structures,44,45 respectively (both synthesizable at ambient temperature), and TiO2(II) of R-PbO2 (columbite) structure46 (which is frequently obtained under high pressure but can also be obtained at 1 bar47 and has been found in mineral samples as well48). These structures will be labeled (B), (R), (H), and (II), respectively. It also includes two phases having 7-coordinated Ti, the baddeleyite (monoclinic zirconia) structure49 (here labeled BD), and the OI phase,50 as well as one with Ti in (7 + 2) coordination, the cotunnite-type phase51 (here labeled CT). These three latter phases with Ti in greater than 6-fold coordination have been made in the laboratory only under pressure; although it may be noted that the BD phase has been detected recently in a natural mineral sample52 and that the CT phase could be quenched at ambient pressure at 77 K.51 There are reports about a cubic high pressure phase,53 but its structure, conjectured to be of a pyrite or fluorite type, has not been determined; it will not be considered here. In modeling these other titania phases a difficulty arises since to this authors’ knowledge their refractive indices have not been reported, and thus polarizabilities and Cij coefficients cannot be estimated directly. However, examination of the oxygen polarizability in rutile, brookite, and anatase revealed a remarkable linear correlation with the experimental molar volume (Figure 1), and so these data were fitted to a straight line, leading to the equation

RO ) 1.177 + 0.03964VF

(11)

(with RO and VF in Å3). This was used to estimate RO values for the other phases. Other details of the method, as well as the damping function shape, were kept the same as in the previous calculations. The DFT-D results obtained for these other phases, compared to those given for them by standard DFT and to the available experimental values (which unfortunately do not include en-

exptl datac

DFT results

a Vector lengths are in Ångstroms; volumes VF are per TiO2 formula unit in Å3. b kJ/mol. c Structural data are from the literature references quoted after the crystal form label. For BD the data, measured at ambient pressure, are for the natural form (akaogiite). For OI and CT phases, the data are those measured at P ) 28 and 60 GPa, respectively, and VF0 refers to the VF value extrapolated to zero pressure from a fitted third-degree Birch-Murnaghan equation of state.

thalpy data), are presented in Table 3. It can be seen that, as was the case in anatase and brookite, some other phases (in particular, B and II) are also erroneously predicted by standard

Relevance of Dispersion Interactions

Figure 2. Value of the Grimme-type correction term (as from eq 1) computed for the different TiO2 crystal forms, plotted against the square of the density (both data after full relaxation of the structure at zero pressure and temperature). The straight line, drawn as guide to the eye, highlights the correlation apparent for most of the phases containing 6-coordinated Ti.

DFT to have ∆E < 0, i.e., to be more stable than rutile; including the DI corrections leads however to the correct ∆E > 0 result for all phases. Note that these two phases (B and II) remain the less unstable of all of those listed in Table 3; this can be linked to the fact that, as previously stated, they have been found in natural minerals. Note also that the BD phase, which according to the standard DFT calculations would have an energy only slightly above those of the R and H phases (both of which can be prepared at ambient pressure), is predicted by the DFT-D method to have an ∆E value much higher than the latter, in accordance with the fact that it is much less stable and can be made in the laboratory only under high pressure. Still, the ∆E of the BD phase remains rather lower than that of phases OI and CT, which may be related to the fact that among the three high-pressure phases only BD has been found (until now) in a natural mineral. Thus, the order of stabilities predicted with the DFT-D method roughly agrees with the experimental observations, and the set of results obtained again supports the idea that DI are relevant in determining the real relative stabilities of the different titania phases. The total values obtained for the DI correction term in all the TiO2 phases studied are plotted in Figure 2 vs the square of the density obtained in the calculation. This plot is expected to be of possible interest, since according to eq 1 this term scales as r-6, i.e., roughly the inverse of the squared volume. A correlation can be seen mainly within the phases having 6-coordinated Ti, especially if one excludes the R and (especially) H phases, which appear with a relatively higher correction term probably because their structures contain intracrystalline voids that lead to higher VF values and hence to the prediction of high (possibly too high) RO values through eq 11. Phases with higher Ti coordination numbers clearly fall outside the correlation; the effect of their lower VF on the estimated RO values, and maybe also a more important damping through function f, may be the reason. The reliability of the results for these phases might therefore be lower. The relevant fact that remains, however, is the high value of the correction term (larger than 1.5 eV per TiO2 unit in absolute value) predicted for all these phases; it is clear that a term of this magnitude could be relevant in determining the overall cohesion energy of these oxides.

J. Phys. Chem. C, Vol. 114, No. 51, 2010 22723 Al2O3 Phases. A similar study to above was also undertaken for the experimentally known alumina polymorphs. In addition to the most stable R phase (corundum, with all Al octahedral), transition aluminas were studied. Two of these have a wellestablished structure: θ, with octahedral and tetrahedral Al in ratio O:T ) 1:1, and κ, where O:T ) 3:1. There is one additional transition alumina of great practical interest: the γ polymorph. This is a disordered structure, formed by decomposition of the boehmite phase AlO(OH), in which the O anions form a facecentered cubic lattice, and the cations have been frequently assumed in the literature to occupy spinel-type positions leaving some of the latter more or less randomly unoccupied. For the calculations one must opt for an ordered model, and it is not clear which one to choose (besides, the energy of this model might well be different from that of the experimentally handled disordered material). Furthermore, there have been in recent years strong indications that some nonspinel sites are occupied,54 in agreement with earlier suggestions.55,56 Here the ordered model obtained by Toulhoat et al.57,58 (labeled γT) that specifies the presence of cations in nonspinel positions and has a ratio O:T ) 6:2, will be used, as well as a monoclinic model proposed by Mene´ndez-Proupin and Gutie´rrez59 (labeled γG) derived from a tripled primitive cell of the spinel structure in which octahedral cation vacancies are located leading to a ratio O:T ) 5:3. A third model (labeled γO) that has been built here has, like γG, a spinel-type structure with minimum primitive cell volume and only octahedral vacancies but is established by starting from a tetragonally distorted spinel cell, tripling it along the c axis (having a centered tetragonal Bravais lattice, with space group I41/amd, at this moment), and then locating two octahedral cation vacancies in the resulting primitive cell in the most stable way (which then contains 40 atoms and has its symmetry lowered to monoclinic space group C2/c).60 Another known relevant transition alumina is the tetragonal δ phase that usually forms during the transformation from the γ to the R phase. Although similar to the γ phase, it is normally rather better crystallized, and as in the γO or γG models, its structure seems to be derived from that of a spinel by locating cation vacancies in the octahedral sites. The structure has a clearly tetragonal character with a tripled spinel cell and space group P-4m2 (Nr. 115);61 recent X-ray diffraction data confirm the same symmetry.62 It should be noted that it probably coincides with the so-called γ′-Al2O3 phase,63 indexed within the same space group and cell volume and which according to diffraction data has a rather complex distribution of partial occupancies. As a model for this tetragonal δ (or γ′) phase, the structure of γ-Fe2O364 (labeled δF) with a tripled spinel cell having cation vacancies ordered in octahedral sites according to a tetragonal symmetry (space group P41212, Nr. 92) will be used here. A number of other intermediate phases of orthorhombic or monoclinic symmetry (with diverse names such as δ′, orthorhombic δ, monoclinic θ′, θ′′, λ, or U), identified with transmission electron microscopy and/or electron diffraction, have been reported65,66 to appear in a dense form or in thin film specimens grown on aluminum metal and other surfaces but will not be considered here since their detailed atomic structures have not been clarified. The other two known Al2O3 phases, which are prepared only at high pressure, will be examined as well: that with Rh2S3type structure (where all Al is 6-fold coordinated),67 which will be labeled RH, and that with a CaIrO3-type structure (also with 6-fold coordinated Al, although 50% of it has trigonal prismatic geometry),68 which will be labeled CI. For this latter phase, the

22724

J. Phys. Chem. C, Vol. 114, No. 51, 2010

molecular volume extrapolated to zero pressure has been determined through its equation of state.69 In this study on alumina phases an additional difficulty arises from the fact that, among all of them, reliable refractive index values are known only for corundum (〈n〉 ) 1.765). Besides this, there is only one report of an experimentally determined refractive index for thin films of (supposedly pure) γ-Al2O3, for which 〈n〉 ) 1.61 is claimed;70 but since neither a dense character nor the absence of amorphous components have been proven for this specimen, such a value cannot be taken into consideration. In this situation, and to keep things simple, the same COO values and damping function shape used for corundum will be assumed for all phases. By use of the RAl value given by Dimitrov,20 the parameters obtained for the DI correction term are those given in Table 1. It can be noted that, with these values, the aforementioned τ factor of Cij overestimation by Grimme’s geometric mean would amount to 1.34. The calculations carried out using both standard DFT and DFT-D methods yield the structural and energetic values given in Table 4, which also includes the experimental data available. Here the DFT-D method does not lead uniformly to contracted volumes; lower values of the damping function f might be responsible for this. With regard to energies, both methods correctly predict that corundum is the most stable phase; thus standard DFT has no qualitative error needing correction. For the Al2O3 system only one experimental energy datum seems to be available: the γ phase is less stable than corundum by ∆E ) +29 kJ/mol.17 Which of the structures studied here would be the best model for the γ specimen used in this experimental work may be questionable, but one observes that although the computed ∆E is, for the models of γ-Al2O3 that appear most stable (γT and γO), lower than the experimental value with both DFT and DFT-D, the latter method gives higher ∆E values, being thus closer to that of the experiment. Also, although the difference δ∆E between ∆E values that is given by both methods (δ∆E ) 4.4 and 4.8 kJ/mol for γT and γO, respectively) is smaller than that found in the case of TiO2, it is still not negligible in comparison with the experimental ∆E value. Thus DI appears again to be significant for determining the relative energies of the different alumina polymorphs. Apart of this, it may be remarked that the δF model of δ-Al2O3, with all cation vacancies ordered in octahedral spinel sites, has lower energy than any of the models used for γ-Al2O3, in spite of the fact that among the latter the nonspinel γT, having higher O:T ratios, appears to be the most stable one. Finally, the θ and κ phases are predicted to be less stable than R, and overall θ is seen to be the most stable of the non-R forms, while significantly higher energies are predicted for the CH and CI models in agreement with the higher difficulty in their experimental synthesis. The relative stability order of all these phases is in general predicted to be the same by the DFT and DFT-D methods and consistent with the experimental observations. No attempt to model these systems with the vdW-DF functional was undertaken, as it is not clear which model would be the most suitable one for representing the γ-Al2O3 specimens from which the experimental ∆E value of ref 17 was obtained. Conclusions Differences in the energies involved in DI, modeled with correction terms of r-6 form (Grimme’s method DFT-D), are seen to affect significantly the relative stabilities of different TiO2 polymorphs. It must be noted that the parameters used to evaluate these terms that are given in the original Grimme paper, which are appropriate for neutral species, are likely to be

Conesa TABLE 4: Crystal Cell Parametersa and Energy Differences ∆Eb in Respect to Corundum, Computed for Al2O3 Phases, and Compared to Experimental Datac crystal form/model DFT results DFT-D results

exptl data

R

a ) 4.8042 c ) 13.1109 VF ) 43.68 c/a ) 2.729

a ) 4.806 c ) 13.162 VF ) 43.87 c/a ) 2.739

a ) 4.762 c ) 12.999 VF ) 42.54 c/a ) 2.730

θ

a ) 11.910 b ) 2.938 c ) 5.667 β ) 75.98° VF ) 48.10 c/a ) 0.476 c/b ) 1.929 ∆E ) +4.63

a ) 11.919 b ) 2.932 c ) 5.667 β ) 76.15° VF ) 48.07 c/a ) 0.475 c/b ) 1.933 ∆E ) +10.7

a ) 11.795 b ) 2.910 c ) 5.621 β ) 76.21° VF ) 46.84 c/a ) 0.477 c/b ) 1.932

κ

a ) 4.882 b ) 8.390 c ) 9.018 VF ) 46.17 c/a ) 1.847 c/b ) 1.075 ∆E ) +8.2

a ) 4.885 b ) 8.392 c ) 9.015 VF ) 46.20 c/a ) 1.845 c/b ) 1.074 ∆E ) +11.5

a ) 4.834 b ) 8.310 c ) 8.935 VF ) 44.86 c/a ) 1.848 c/b ) 1.075

γD

a ) 5.576 b ) 8.398 c ) 8.068 β ) 90.53° VF ) 47.23 c/a ) 1.447 c/b ) 0.961 ∆E ) +13.3

a ) 5.568 b ) 8.393 c ) 8.055 β ) 90.63° VF ) 47.055 c/a ) 1.447 c/b ) 0.960 ∆E ) +17.7

data for (cubic disordered) γ phase: a ) 7.952; VF ) 47.14; ∆E ) +2917

γG

a ) 9.751 b ) 5.677 c ) 13.667 β ) 89.38° VF ) 47.28 c/a ) 1.402 c/b ) 2.408 ∆E ) +20.3

a ) 9.7468 b ) 5.6762 c ) 13.664 β ) 89.37° VF ) 47.24 c/a ) 1.402 c/b ) 2.407 ∆E ) +25.0

γO

a ) 8.003 b ) 7.988 c ) 12.558 β ) 109.75° VF ) 47.23 c/a ) 1.569 c/b ) 1.572 ∆E ) 15.9

a ) 7.999 b ) 7.989 c ) 12.555 β ) 109.75° VF ) 47.19 c/a ) 1.570 c/b ) 1.572 ∆E ) 20.7

δF

a ) 7.981 c ) 23.887 VF ) 46.67 c/a ) 2.993 ∆E ) +7.8

a ) 7.974 c ) 23.877 VF ) 47.45 c/a ) 2.994 ∆E ) +13.8

a ) 5.599 c ) 23.657 VF ) 46.35 c/(2a) ) 2.988

RH

a ) 7.088 b ) 4.841 c ) 4.981 VF ) 42.72 c/a ) 0.703 c/b ) 1.029 ∆E ) +45.0

a ) 7.0892 b ) 4.8425 c ) 4.9836 VF ) 42.771 c/a ) 0.703 c/b ) 1.029 ∆E ) +47.9

a ) 6.393 b ) 4.362 c ) 4.543 VF ) 31.67 c/a ) 0.711 c/b ) 1.041

CI

a ) 2.675 b ) 8.938 c ) 7.036 VF ) 42.06 c/a ) 2.630 c/b ) 0.787 ∆E ) +135

a ) 2.677 b ) 8.952 c ) 7.050 VF ) 42.24 c/a ) 2.633 c/b ) 0.787 ∆E ) +137

a ) 2.402 b ) 7.863 c ) 5.989 VF ) 28.28 (VF0 ) 39.6) c/a ) 2.493 c/b ) 0.762

a Vector lengths are in Ångstroms; volumes VF are per Al2O3 formula unit in Å3. b kJ/mol. c For references on the experimental data of different polymorphs and on the models for γ and δ phases, see text. For the RH and CI phases the structural data are obtained at pressures of 113 and 163 GPa, respectively; VF0 for the CI phase refers to the VF value extrapolated to zero pressure from a fitted equation of state.

inadequate when modeling oxides and other highly polar inorganic solids. These coefficients can be derived rather, as

Relevance of Dispersion Interactions done here, from refractive indices (if these are available) through the calculation of molecular polarizabilities. The contributions of these terms to total energies can reach values in the order of 1 eV per formula unit, and suffice to correct the known failure that DFT and hybrid DFT methods have when trying to reproduce the correct result that rutile is the most stable of the TiO2 polymorphs. Furthermore the excess energies ∆E of anatase and brookite in respect to rutile computed with this DFT-D method approach reasonably and with the correct order those experimentally determined in the literature. An attempt to include DI in the calculation through the use 1 of a fully nonlocal functional (the vdW-DF method ) produces for both anatase and brookite less negative ∆E values, although the correct positive sign is not reached. It remains to be seen whether the result could be improved by modifying the exchange functional or including some semilocal component in the correlation functional, as suggested in the literature. Still, the magnitude of the effect found with the vdW-DF method is significant if one compares it with the experimental data, and together with the results obtained here with Grimme’s method (once its parameters are adapted to these oxides) for TiO2 and also Al2O3 polymorphs allows this author to conclude that any method that tries to reproduce the true relative energies of oxide phases (and probably also those of other ionic solids) should not consider the contribution of DI as negligible. Certainly the present level of theory (whether DFT-D or vdWDF), as applied here, is not accurate enough for these purposes. Other more reliable approaches will be needed for this. In fact, it is proposed here that the good reproduction of experimental ∆E values in oxide polymorphs be used as an additional validity test for any present or future theory of DI and nonlocal functionals. This would facilitate, in addition, a more reliable modeling of DI when studying weak adsorption on oxides. Of course this is also a reason for asking experimentalists to provide a higher amount of reference values of this kind. Further work to assess the efect of DI on other properties of these oxides may be of interest. For example, bulk moduli and phonon frequencies could be computed and compared with experimental data where these are available. The effect of DI on the energies of the oxide surfaces (and on the energies of adsorption of molecules on them) also remains to be proven. Work on these issues is in progress; although the surface-related task has additional difficulties, since the polarizabilities and Cij coefficients corresponding to the surface ions are likely to be significantly higher than those of the bulk ones, and there is no simple way to obtain their values from refractive indices or similar properties. Acknowledgment. Thanks are given to Tomas Bucko for kindly providing corrected routines of the VASP code that calculates DI effects through Grimme’s method and to Timo Thonhauser for kindly providing the development version of the routines in the Quantum Espresso code that computes the nonlocal vdW-DF functional. Financial aid from the Spanish Plan Nacional de Investigacio´n (Project Nos. CTQ2006-15600 and MAT2009-14625) and from COST actions D41 and 540 is also acknowledged. References and Notes (1) (a) Dion, M.; Rydberg, H.; Schro¨der, E.; Langreth, D. C.; Lundqvist, B. I. Phys. ReV. Lett. 2004, 92, 246401. (b) Dion, M.; Rydberg, H.; Schro¨der, E.; Langreth, D. C.; Lundqvist, B. I. Phys. ReV. Lett. 2005, 95, 109902. (c) Thonhauser, T.; Cooper, V. R.; Li, S.; Puzder, A.; Hyldgaard, P.; Langreth, D. C. Phys. ReV. B 2007, 76, 125112. (d) Vydrov, O. A.; Wu, Q.; Van Voorhis, T. J. Chem. Phys. 2008, 129, 014106.

J. Phys. Chem. C, Vol. 114, No. 51, 2010 22725 (2) (a) Vydrov, O. A.; Van Voorhis, T. J. Chem. Phys. 2009, 130, 104105. (b) Vydrov, O. A.; Van Voorhis, T. Phys. ReV. Lett. 2009, 103, 063004. (3) Harl, J.; Kresse, G. Phys. ReV. B 2008, 77, 045136. (4) Tavernelli, I.; Lin, I.-C.; Rothlisberger, U. Phys. ReV. B 2009, 79, 045106. (5) Langreth, D. C.; Lundqvist, B. I.; Chakarova-Kack, S. D.; Cooper, V. R.; Dion, M.; Hyldgaard, P.; Kelkkanen, A.; Kleis, J.; Kong, L.; Li, S.; Moses, P. G.; Murray, E.; Puzder, A.; Rydberg, H.; Schroder, E.; Thonhauser, T. J. Phys.: Condens. Matter 2009, 21, 084203. (6) (a) Langreth, D. C.; Lundqvist, B. I. Phys. ReV. Lett. 2010, 104, 099303. (b) Vydrov, O. A.; Van Voorhis, T. Phys. ReV. Lett. 2010, 104, 099304. (7) Lee, K.; Murray, E. D.; Kong, L.; Lundqvist, B. I.; Langreth, D. C. Phys. ReV. B 2010, 82, 081101(R). (8) (a) Klimesˇ, J.; Bowler, D. R.; Michaelides, A. J. Phys.: Condens. Matter 2010, 22, 022201. (b) Cooper, V. R. Phys. ReV. B 2010, 81, 161104. (9) Romaner, L; Nabok, D.; Puschnig, P.; Zojer, E.; Ambrosch-Draxl, C. New J. Phys. 2009, 11, 053010. (10) Grimme, S. J. Comput. Chem. 2004, 25, 1463. (11) Grimme, S. J. Comput. Chem. 2006, 27, 1787. (12) (a) Becke, A. D.; Johnson, E. R. J. Chem. Phys. 2005, 122, 154104. (b) Silvestrelli, P. L. Phys. ReV. Lett. 2008, 100, 053002. (c) Tkatchenko, A.; Scheffler, M. Phys. ReV. Lett. 2009, 102, 073005. (13) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. J. Chem. Phys. 2010, 132, 154104. (14) Catlow, C. R. Computational Materials Science; Catlow, R. C.; Kotomin, E., Eds.; IOS Press: Amsterdam, 2003; Vol. 187, 1-29. (15) Smith, S. J.; Stevens, R.; Liu, S.; Li, G.; Navrotsky, A.; BoerioGoates, J.; Woodfield, B. F. Am. Mineral. 2009, 94, 236. (16) Ranade, M. R.; Navrotsky, A.; Zhang, H. Z.; Banfield, J. F.; Elder, S. H.; Zaban, A.; Borse, P. H.; Kulkarni, S. K.; Doran, G. S.; Whitfield, H. J. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 6476. (17) Castro, R. H.; Ushakov, S. V.; Gengembre, L.; Gouveˆa, D.; Navrotsky, A. Chem. Mater. 2006, 18, 1867. (18) Cambi, R.; Cappelletti, D.; Liuti, G.; Pirani, F. J. Chem. Phys. 1991, 95, 1852. (19) Tessman, J. R.; Kahn, A. H.; Shockley, W. Phys. ReV. 1953, 92, 8902. (20) Dimitrov, V.; Sakka, S. J. Appl. Phys. 1996, 79, 173. (21) Dimitrov, V.; Komatsu, T. J. Solid State Chem. 2002, 163, 100. (22) Kordes, E. Z. Phys. Chem. B 1939, 44, 249. (23) Fawcett, W. Proc. Phys. Soc. 1963, 82, 33. (24) Wu, Q.; Yang, W. J. Chem. Phys. 2002, 116, 515. (25) Tang, K. T.; Toennies, P. J. Chem. Phys. 1984, 80, 3726. (26) Sato, T.; Nakai, H. J. Chem. Phys. 2009, 131, 224104. (27) Shannon, R. D. Acta Crystallogr. A 1976, 32, 751. (28) (a) Kresse, G.; Hafner, J. Phys. ReV. B 1993, 47, 558. (b) Kresse, G.; Hafner, J. J. Phys.: Condens. Matter 1994, 6, 8245. (c) Kresse, G.; Furthmu¨ller, J. Comput. Mater. Sci. 1996, 6, 15. (29) http://cms.mpi.univie.ac.at/marsweb/. (30) Kresse, G.; Joubert, D. Phys. ReV. B 1999, 59, 1758. (31) Blum, V.; Gehrke, R.; Hanke, F.; Havu, P.; Havu, V.; Ren, X.; Reuter, K.; Scheffler, M. Comput. Phys. Commun. 2009, 180, 2175. (32) http://www.fhi-berlin.mpg.de/aims/. (33) Roma´n-Pe´rez, G.; Soler, J. M. Phys. ReV. Lett. 2009, 103, 096102. (34) Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G. L.; Cococcioni, M.; Dabo, I.; Dal Corso, A.; Fabris, S.; Fratesi, G.; de Gironcoli, S.; Gebauer, R.; Gerstmann, U.; Gougoussis, C.; Kokalj, A.; Lazzeri, M.; Martin-Samos, L.; Marzari, N.; Mauri, F.; Mazzarello, R.; Paolini, S.; Pasquarello, V; Paulatto, L.; Sbraccia, C.; Scandolo, S.; Sclauzero, G.; Seitsonen, A. P.; Smogunov, A.; Umari, P.; Wentzcovitch, R. M. J. Phys.: Condens. Matter 2009, 395502. A development version of the code made by T. Thonhauser et al. following the principles described in refs 1 and 34 was used for this. (35) http://www.pwscf.org/. (36) Burdett, J. K.; Hughbanks, T.; Miller, G. J.; Richardson, J. W., Jr.; Smith, J. V. J. Am. Chem. Soc. 1987, 109, 3639. (37) Djerdj, I.; Tonejc, A. M. J. Alloys Compd. 2006, 413, 159. (38) Muscat, J.; Swamy, V.; Harrison, N. M. Phys. ReV. B 2002, 65, 224112. (39) Labat, F.; Baranek, P.; Domain, C.; Minot, C.; Adamo, C. J. Chem. Phys. 2007, 126, 154703. (40) Ma, X. G.; Liang, P.; Miao, L.; Bie, S. W.; Zhang, C. K.; Xu, L.; Jiang, J. J. Phys. Status Solidi B 2009, 246, 2132. (41) Vittadini, A. Unpublished result, personally communicated to this author concomitant to the editing process. (42) Marchand, R.; Brohan, L.; Tournoux, M. Mater. Res. Bull. 1980, 15, 1129. (43) Banfield, J. F.; Veblen, D. R.; Smith, D. J. Am. Mineral. 1991, 76, 343. (44) Akimoto, J.; Gotoh, Y.; Oosawa, Y.; Nonose, N.; Kumagai, T.; Aoki, Y.; Takei, H. J. Solid State Chem. 1994, 113, 27.

22726

J. Phys. Chem. C, Vol. 114, No. 51, 2010

(45) Sasaki, T.; Watanabe, M.; Fujiki, Y. Acta Cryst. B 1993, 49, 838. (46) Simons, P. Y.; Dachille, F. Acta Crystallogr. 1967, 23, 334. (47) Grey, I. E.; Li, C.; Madsen, I. C.; Braunshausen, G. Mater. Res. Bull. 1988, 23, 743. (48) Hwang, S.-L.; Shen, P.; Chu, H.-T.; Yui, T.-F. Science 2000, 288, 321. (49) Sato, H.; Endo, S.; Sugiyama, M.; Kikegawa, T.; Shimomura, O.; Kusaba, K. Science 1991, 251, 786. (50) Dubrovinskaia, N. A.; Dubrovinsky, L. S.; Ahuja, R.; Prokopenko, V. B.; Dmitriev, V.; Weber, H.-P.; Osorio-Guillen, J. M.; Johansson, B. Phys. ReV. Lett. 2001, 87, 275501. (51) Dubrovinsky, L. S.; Dubrovinskaia, N. A.; Swamy, V.; Muscat, J.; Harrison, N. M.; Ahuja, R.; Holm, B.; Johansson, B. Nature 2001, 410, 653. (52) El Goresy, A.; Dubrovinsky, L.; Gillet, P.; Graup, G.; Chen, M. Am. Mineral. 2010, 95, 892. (53) Mattesini, M.; de Almeida, J. S.; Dubrovinsky, L.; Dubrovinskaia, N.; Johansson, B.; Ahuja, R. Phys. ReV. B 2004, 70, 212101. (54) Paglia, G.; Buckley, C. E.; Rohl, A. L.; Hunter, B. A.; Hart, R. D.; Hanna, J. V.; Byrne, L. T. Phys. ReV. B 2003, 68, 144110. (55) Zhou, R.-S.; Snyder, R. L. Acta Crystallogr. B 1991, 47, 617. (56) Wolverton, C.; Hass, K. C. Phys. ReV. B 2000, 63, 024102. (57) Krokidis, X.; Raybaud, P.; Gobichon, A.-E.; Rebours, B.; Euzen, P.; Toulhoat, H. J. Phys. Chem. B 2001, 105, 5121.

Conesa (58) Digne, M.; Sautet, P.; Raybaud, P.; Euzen, P.; Toulhoat, H. J. Catal. 2004, 226, 54. (59) Mene´ndez-Proupin, E.; Gutie´rrez, G. Phys. ReV. B 2005, 72, 035116. (60) Generated ordered structure (before relaxation): space group C2/c (Nr. 15), a ) 7.9 Å, b ) 7.9 Å, c ) 12.49 Å, β ) 108.44°, Al1 ) (0.8334, 0.5000, 0.1668), Al2 ) (0.9166, 0.2500, 0.8332), Al3 ) (0.2500, 0.2500, 0.5000), Al4 ) (0.0000, 0.1250, 0.2500), Al5 ) (0.8334, 0.3750, 0.4168), O1 ) (0.9234, 0.510, 0.3266), O2 ) (0.8468, 0.260, 0.6734), O3 ) (0.5900, 0.510, 0.6598), O4 ) (0.1802, 0.260, 0.3402), O5 ) (0.2567, 0.510, 0.9932), O6 ) (0.0135, 0.760, 0.0068). (61) Repelin, Y.; Husson, E. Mater. Res. Bull. 1990, 25, 611. (62) Tsybulya, S. V.; Kryukova, G. N. Powder Diffr. 2003, 18, 309. (63) Paglia, G.; Buckley, C. E.; Rohl, A. L.; Hart, R. D.; Winter, K.; Studer, A. J.; Hunter, B. A.; Hanna, J. V. Chem. Mater. 2004, 16, 220. (64) Jorgensen, J. E.; Mosegaard, L.; Thomsen, L. E.; Jensen, T. R.; Hanson, J. C. J. Solid State Chem. 2007, 180, 180. (65) Levin, I.; Brandon, D. J. Am. Ceram. Soc. 1998, 81, 1995. (66) Devi, U. M. Ceram. Intern. 2004, 30, 555. (67) Lin, J. F.; Degtyareva, O.; Prewitt, C. T.; Dera, P.; Sata, N.; Gregoryanz, E.; Mao, H.-K.; Hemley, R. J. Nat. Mater. 2004, 3, 389. (68) Oganov, A. R.; Ono, S. Proc. Natl. Acad. Sci. 2005, 102, 10828. (69) Ono, S.; Oganov, A. R.; Koyama, T.; Shimizu, H. Earth Planet. Sci. Lett. 2006, 246, 326. (70) Murali, K. R.; Thirumoorthy, P. J. Alloys Compd. 2010, 500, 93.

JP109105G