- Email: [email protected]
Influence of elastic properties on superplasticity in doped yttria-stabilized zirconia
ARTICLE IN PRESS Journal of Physics and Chemistry of Solids 70 (2009) 15–19
Contents lists available at ScienceDirect
Journal of Physics and Chemistry of Solids journal homepage: www.elsevier.com/locate/jpcs
Influence of elastic properties on superplasticity in doped yttria-stabilized zirconia Yuriy Natanzon a, Marek Boniecki b, Zbigniew Łodziana a,c, a b c
´ ski Institute of Nuclear Physics, Polish Academy of Sciences, ul. Radzikowskiego 152, 31-342 Krako ´w, Poland Henryk Niewodniczan ´lczyn ´ ska 133, 01-919 Warsaw, Poland Institute of Electronic Materials Technology, ul. Wo ¨ bendorf, Switzerland Department of Environment, Energy and Mobility, EMPA, 8600 Du
a r t i c l e in fo
abstract
Article history: Received 14 July 2008 Received in revised form 12 August 2008 Accepted 17 August 2008
We report the first principles calculations of elastic and electronic properties of yttria-stabilized tetragonal zirconium dioxide (YZP) doped with GeO2 , TiO2 and SiO2 . Electronic structure and isotropic elastic properties of YZP do not change upon addition of dopants. Addition of dopants affects the shear C 66 elastic constant that decreases with the increasing dopant concentration. A simple model that connects elastic softening to enhancement of superplasticity in doped fine-grained zirconia ceramics is proposed. & 2008 Elsevier Ltd. All rights reserved.
Keywords: A. Oxides C. Ab initio calculations D. Elastic properties D. Electronic structure
1. Introduction For the first time superplasticity in tetragonal yttria-stabilized zirconia (YZP) was observed over than 20 years ago by Wakai et al. [1,2]. Since that time superplastic properties of YZP have been an exciting field of research (see review of Jimenez-Melendo et al. [3]). At present there is no commonly accepted explanation of this phenomenon. One of the most accepted explanation of superplasticity in YZP is grain boundary sliding related to deformation. It was confirmed by number of microstructural observations [2,4]. Superplasticity is driven by the diffusion of Zr IV ions at the grain boundary (or so-called diffusion creep) [5]. Yttrium dopants segregate toward the surface on grain boundary and that enhances Zr diffusion [6–8]. Such diffusion contributes to the increase of superplastic strain rate and facilitates intergranular sliding during deformation [3]. Go´mez-Garcı´a et al. [8] have shown that segregation of Y cations cause electrostatic potential difference at the grain boundary. Intergranular sliding is then supported by the electrostatic interaction between the grains. It is well known that superplastic properties of ZrO2 are drastically changed by addition of small amount of metal oxides [9]. Some dopants can significantly increase the superplastic strain rate of zirconia to allow the elongation of up to 1000% of its
initial length [9,10]. This behavior is confirmed for such dopants as GeO2 [9–11], TiO2 [9–12], SiO2 [4,13,14], MgO [11,12] and Al2 O3 [13,14] as well as their combinations. It was suggested that changes in superplasticity are due to the formation of glassy thin film of dopant, which acts as a lubricant for intergranular sliding [2]. However, the effect is still observed when the intermediate glass layer is not present but instead the dopant segregation near the grain boundary is observed [15]. Another reason for the effect could be the related changes in ionic bonding strength of cation–oxygen bonds in the vicinity of grain boundary [9,10]. Unfortunately the existing models do not provide clear microscopic picture of the superplasticity in ZrO2 related to metal oxide doping. Especially there is lack of atomic level insight in the literature [9,10]. In order to gain more detailed microscopic insight into phenomenon of superplasticity in doped YZP we have carried out a number of first principles calculations of elastic properties of yttria-stabilized tetragonal zirconia with GeO2 , TiO2 and SiO2 dopant metal oxides. Our calculations indicate the correlation between changes in shear elastic constants and superplasticity. We show that some of the shear elastic properties depend on Y/dopant concentration while the isotropic properties, like bulk modulus do not.
2. Model and method of calculation Corresponding author at: Henryk Niewodniczan´ski Institute of Nuclear Physics,
Polish Academy of Sciences, ul. Radzikowskiego 152, 31-342 Krako´w, Poland. E-mail addresses: [email protected] (Y. Natanzon), Marek.Boniecki@ itme.edu.pl (M. Boniecki), [email protected] (Z. Łodziana). 0022-3697/$ - see front matter & 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.jpcs.2008.08.014
Depending on the external thermodynamic conditions three phases of ZrO2 can be stable. Phase transition sequence of ZrO2 includes high temperature cubic phase ðFm3mÞ which exists
ARTICLE IN PRESS 16
Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
between T ¼ 2650–2990 K and tetragonal phase (P42 =nmc) at T ¼ 1400–2650 K. Monoclinic phase (P21 =c) exists at the temperatures below 1400 K [16–18]. However, high temperature phases can be stabilized with various metal oxides like MgO, CaO, CeO2 and Y2 O3 and exist at lower temperatures. As seen from phase diagram of YZP in Ref. [19], tetragonal phase of zirconia exists for yttria concentrations below 6 mol% and temperature range of about 600–2000 K. In this paper we consider tetragonal phase of zirconia with the concentration of yttria below 6.7 mol%. Our present calculations were performed within densityfunctional theory (DFT) approach [20] with generalized-gradient approximation (GGA) [21] for exchange-correlation functional. The linear combination of atomic orbitals (LCAO) basis set was used, as implemented in SIESTA program [22]. The implementation of DFT used in the present approach uses specific techniques to speed up the calculations and balance between the desired accuracy and calculation time. Only the valence electrons which take part in bond formation are considered. The other electrons (also called ‘‘core electrons’’) are replaced by effective norm-conserving potentials with Troullier and Martins parametrization [23]. In this calculation we have used GGA pseudopotentials in Perdew–Burke–Ernzerhof parametrization [24]. In order to achieve better accuracy we have also treated 4s and 4p core electrons of Zr and Y atoms as valence 2 states. Thus valence electron configuration was 4s2 4p6 5s2 4d for 1 2 2 6 2 2 6 2 2 2 2 Zr, 4s 4p 5s 4d for Y, 3s 3p 4s 3d for Ti, 3s 3p for Si, 4s 4p2 for Ge, 3s2 3p1 for Al and 3s2 for Mg. Secondly, all the integrals in SIESTA are calculated on the finite real space grid. Plane-wave cutoff of 476 eV (350 Ry) was used in this calculation and a k point grid with a step of 0.03 A˚. The wave functions are expanded into strictly localized atomic orbitals, which consist of the product of numerical radial function and spherical harmonics [22,25]. Double-z with polarization (DZP) basis set was used for the valence electrons of all atoms and double-z (DZ) for the semicore 4s and 4p electrons of Zr and Y. All the pseudopotentials and basis sets were properly tested on the corresponding oxides to reproduce good experimental values of lattice constants and bulk modulus. Additionally, we have calculated total energies for different phases of ZrO2 in order to check if the energy sequence is correct and states Emonoclinic oEtetragonal oEcubic [19]. In the present studies the supercell technique was used to vary the concentration of yttrium and other dopants [26]. At the first step experimental unit cell parameters of ZrO2 were optimized and bulk modulus was calculated in order to assess the accuracy of calculations. The calculated lattice constants and bulk modulus are shown in Table 1. Calculated bulk modulus agree within 5% with this measured by XRD [18] and the lattice parameters are in good agreement with previous theoretical [17,27] and experimental results [18,28]. To study stabilized ZrO2 two formula units of zirconia were replaced with one formula unit of yttria. In order to keep the
Table 1 Comparison of calculated lattice constants and elastic properties of pure ZrO2 with experimental measurements and other calculations Lattice constants (A˚)
This work (GGA) Ref. [17] (LDA) Ref. [27] (GGA) Ref. [28] (neutron diffraction) Ref. [18] (XRD)
a¼b
c
3.61 3.66 3.645 3.60 3.59
5.20 5.30 5.29 5.18 5.18
Bulk modulus (GPa)
209.5 200 – – 198
LDA stands for local density approximation, GGA for generalized gradients approximation and X-ray diffraction is abbreviated as XRD.
correct stoichiometry for each two Y dopant atoms one oxygen should be removed. Relative position of Y atoms and oxygen vacancy was chosen according to results of Bogicevic et al. [19,29] and Eichler [27]. Additional tests confirm that Y atoms located in ða; bÞ plane at the distance of 3.5 A˚ and the oxygen being the first neighbor of dopants provide the lowest energy configuration. The optimized unit cell of zirconia was multiplied in selected directions to obtain desired size of the supercell. The concentrations of Y and dopant oxides were the same depending on the size of the supercell as each supercell contained two Y dopants and one metal oxide dopant (the smaller the cell the larger the concentration and vice versa). In fact, the concentration of 6.7 mol% goes actually beyond the stability of the tetragonal phase and corresponds already to the cubic ZrO2 . However, the concentration of dopants increases towards the surface [6] and such a large concentration was used for the comparison purposes. For each dopant concentration and for each supercell size we have performed full structural relaxation (the supercell size and atomic positions within the supercell) to obtain the ground state configuration of the system. The conjugate-gradients method was used for such an optimization [30].
3. Results 3.1. Electronic structure Tetragonal zirconium dioxide is an insulator. Its measured band gap is Eg ¼ 4:2 eV [31], the band gap calculated by Eichler in GGA is Eg ¼ 4:21 eV [27]. In the present calculation the band gap is Eg ¼ 3:7 eV. Slight underestimation is typical for the approximations used here. Stabilization of zirconia with yttrium oxide as well as further aliovalent doping with metal oxides does not significantly influence the electronic structure of ZrO2 [19]. In the present calculation isovalent dopants were studied. Calculated electronic density of states of pure and doped zirconia is shown in Fig. 1 and does not depend on the type of the dopant. Only for titanium oxide dopant some localized unoccupied donor states are observed. The valence band mainly consists of oxygen 2p orbitals and the unoccupied conduction band include empty states of metal cation. The DOS for yttria concentration of 2.8% and 6.7% is presented in the bottom part of Fig. 1. The band gap shrinks by 0.16 eV (5%) for high dopant concentration. Results presented above indicate that macroscopic electronic properties of ZrO2 do not change after the dopant addition. 3.2. Elastic properties The phenomena of superplasticity is related to changes of the macroscopic elastic properties of fine-grained ZrO2 . As it is generally attributed to the grain boundary sliding, fewer attention is given to the elastic properties of the single crystal grains, which also affect the plastic properties of YZP. In this section we present the results that show what is the influence of doping on the elasticity of bulk ZrO2 crystal. The basic elastic properties include volume compressibility (bulk modulus) and tensor of the elastic constants. Below we present how these properties depend on concentration and the type of dopant. 3.2.1. Bulk modulus Bulk moduli were obtained by series of the calculations where uniform lattice deformation in the range of 3% to 3% of its initial volume with the step of 0.5% was applied. For each calculation
ARTICLE IN PRESS Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
17
215
Bulk modulus (GPa)
210
205
200 No dopant Ti Ge Si
195
190 1.5
2
2.5
3
4 3.5 4.5 5 5.5 Yttria concentration (mol%)
6
6.5
7
Fig. 2. Y/dopant concentration dependence of bulk modulus of YZP. Filled circles correspond to undoped YZP and empty symbols to YZP with Ti, Ge or Si dopants. Solid line is a linear fit of the average values for each concentration.
Table 2 Calculated elastic constants of pure ZrO2 in comparison with experimental measurements and other calculations C 11 C 12 C 13 C 33 C 44 C 66 This work (GGA) Ref. [35] (LDA) Ref. [36] (neutron diffraction, 12 mol% Ce-doped ZrO2 ) Ref. [37] (Brillouin scattering, clear monoclinic ZrO2 )
293 248 111 382 221 72 327 100 62 340 33 160
385 346 264 325
51 187 42 167 59 64 66 95
All the values are given in GPa.
analysis of the components of the elastic constants tensor can reveal more detailed information about system response to nonuniform strains that arise during grain sliding. In the next section we analyze the shear properties of YZP. Fig. 1. Total density of states (DOS) for 2.8 mol% doped tetragonal Y-stabilized zirconia (YZP): DOS for pure ZrO2 is shown at the bottom panel, the undoped YZP is second from the bottom (thick line is for Y concentration of 6.7mol%). Upper panels are for YZP with Si, Ge and Ti dopants.
structural relaxation of atomic positions was performed. To obtain the value of bulk modulus the calculated dependence of the total energy E of the system on the supercell volume V was fitted by Murnaghan equation of state [32]. The dependence of the bulk modulus on Y/dopant concentration is presented in Fig. 2. From the linear fit of the average values of bulk modulus one can see that no systematic changes of the bulk modulus with Y/dopant concentration are observed. The values of the bulk modulus are slightly lower than that of pure ZrO2 (see Table 1). The fluctuations are related to the accuracy of calculations (relaxation of large supercells with more than 100 atoms) and as well as to the errors in the Murnaghan equation fit. Significant change of bulk modulus for the concentration of 6.7 mol% is related to the phase transition. For the concentrations above 6 mol% the tetragonal phase is unstable with respect to the cubic one [19]. The observed fluctuations (see Fig. 2) reflect this instability. Weak concentration dependence of the bulk modulus is in good agreement with models of the superplastic behavior of doped ZrO2 which indicate that the macroscopic grain deformation is not related to the increase of superplasticity [3]. However, elastic response of the system is the tensor quantity. Detailed
3.2.2. Shear elastic constants Elastic constants of pure and doped ZrO2 were calculated by series of finite deformations of the supercell as described in the work of Nielsen and Martin [33] and Hongzhi Yao et al. [34]. For each constant and each dopant concentration the appropriate deformation of 1% was applied and the components of the stress tensor were calculated as the derivatives of Kohn–Sham total energy with respect to the strain tensor. Thus elastic constants were calculated by dividing the components of the stress tensor by the components of strain tensor: C nn ¼ sn =n , where n ¼ 4; 5; 6 and i ¼ 0:01 only for the component for i ¼ n. To minimize the numerical errors, each elastic constant was obtained as the average of two calculations with the values of stress of þ1% and 1%. For each deformation a full relaxation of the atomic positions in the supercell was performed. The elastic constant tensor for pure ZrO2 is presented in Table 2. The calculated values of elastic constants are in good agreement to available theoretical results. Some of the components show deviation from the experimental data; however, the measurements are usually performed on stabilized zirconia [36] or on the stable monoclinic phase with the approximation of elastic constants temperature dependence to the tetragonal phase [37]. As will be shown below significant difference between theoretical and experimental values of C 66 elastic constant is due to the dopant or stabilizer addition to ZrO2 .
ARTICLE IN PRESS Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
65
190
60
180 C66 (GPa)
C44 (GPa)
18
55 50 45
No dopant Ti Ge Si
170 160 No dopant Ti Ge Si
150 140
40 1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 Yttria concentration (mol%)
1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 Yttria concentration (mol%)
Fig. 3. Y/dopant concentration dependence of C 44 elastic constant (a) and C 66 elastic constant (b). Filled circles correspond to undoped YZP and empty symbols are for YZP with Ti, Ge or Si dopant. Solid line is a linear fit of the average values for each concentration.
Calculated dopant concentration dependence of C 44 and C 66 elastic constants for doped Y-stabilized ZrO2 are presented in Fig. 3. As one can see from the linear fit of the average values there are no systematic changes in C 44 constant (see Fig. 3(a)) and systematic reduction of the C 66 elastic constant with increasing dopant concentration is observed (see Fig. 3(b)). The largest reduction is observed for Si dopant where for 4.8 mol% increase of dopant concentration C 66 drops by 20%. For yttria-stabilized ZrO2 reduction of shear elastic constants with increasing concentration of Y is also observed, and this change is below 10%.
4. Discussion It is well established that the concentration of dopant increases towards the grain surface in doped YZP [6–8]. According to results presented above this phenomena would affect shear modulus in the grain boundary region. The strain rate of YZP dependence on the inverse of the shear modulus is already predicted from a simple phenomenological model [3,8]. The empirical equation for the strain rate is given by [3]: 2 2 s , G
q ADGb b ¼ kT d qt
(1)
where A is the dimensionless constant, D the grain diffusion coefficient, b the Burgers vector, k the Boltzmann constant, T the temperature, d the grain size, s the applied stress and G the shear modulus. The shear modulus, which is the ratio of shear stress and the shear strain by definition equals the combination shear elastic constants C 66 or C 44 (see Table 2) depending on the crystalline direction. The present calculations are in a good agreement with this relation and provide microscopic insight into the phenomenon. Below we propose a simple model that relates dopant concentration with the shear elastic constants and flow stress. Intergrain sliding in YZP occurs along crystalline facets that are exposed in the crystallites. The Wulff construction of the tetragonal nano-crystallite of ZrO2 is presented in Fig. 4(a). The construction is based on the surface energy calculations by Christensen and Carter [38]. It is clear that the (111) facet will dominate the surface, two other facets that are present on the surface are (1 0 0) and (0 0 1). The (0 0 1) facet is directly related to C 66 elastic constant, while shear deformation of the (111) surface is linear combination of shear related to C 44 (23 of total shear) and C 66 (13 of total shear). As the C 44 shear modulus does not change upon doping, decrease of C 66 will contribute to the decrease of the apparent shear modulus G increasing the strain rate as can be seen from Eq. (1). Below we focus on effects solely related to C 66 shear
Fig. 4. (a) Wulff construction of the zirconia grain shape: the (111) surface is marked as dark gray (blue) and it is dominating one, the (0 0 1) surface at the top (and bottom) is marked as gray (red) and (1 0 0) facet is light gray; (b) a schematic model showing the sliding of the ZrO2 grains. The d denotes the size of the layer where the dopant concentration varies, arrows show the direction of the applied force F and the direction of sliding. The grain surface layers undergo shear deformation during sliding.
modulus, such considerations provide the upper boundary of the influence of elastic softening on the superplastic flow. For simplicity we consider zirconia grains as cuboidal objects sliding along (0 0 1) contact facets following the idea of Go´mezGarcı´a et al. [8]. We assume that in the layer of thickness d from the surface the concentration of dopant increases linearly towards the surface. This assumption together with the result presented in Fig. 3(b) provide the expression for the concentration dependence of C 66 elastic constant: C 66 ¼ C 066 ð1 ayÞ,
(2)
C 066
is the elastic constant for the bulk, a is the coefficient where which depends on the type of dopant. The y is the distance from the bottom of layer d to the contact surface plane; 0pypd, and ayp1. Schematic picture of the model is presented in Fig. 4(b). During sliding grains undergo a shear deformation on the contact surfaces. The work applied on the contact boundary is dW ¼ 2F dy, where F is a force and dy, shear deformation angle. In the linear regime the deformation is proportional to the applied 2 force, thus dW ¼ C 66 dy , and the expression for the deformation angle reads
dy ¼
2F C 066 ð1 ayÞ
.
(3)
With assumptions of the linear change of dopant concentration the deformation of ZrO2 grain in the vicinity of the crystalline contact surface is presented in Fig. 5. With our secure assumption of linear change of concentration and for the linear elastic regime the deformation angle is already doubled at the surface. Such a
ARTICLE IN PRESS Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
19
elastic constants induces increase of the strain rate at the grain surfaces and promotes intergranular sliding by the decrease of intergranular friction. The increase of the strain rate in doped tetragonal phase of Y-stabilized ZrO2 is proportional to the inverse of shear modulus that decrease upon doping. This effect contributes to overall superplasticity besides the electrostatic and kinetic effects.
Acknowledgment The authors acknowledge CPU time allocation at Interdisciplinary Center for Mathematical and Computer Modelling of Warsaw University (Grant no. G28-22). This work was supported by Ministry of Science and Higher Education of Poland under Grant no. 3 T08D 013 29 during the years 2005–2008. References
Fig. 5. Change of the shear deformation angle (y) with respect to depth inside the grain caused by the increase of dopant concentration. y0 is the deformation for constant bulk concentration of dopant. Origin of the x-axis corresponds to deformation in the bulk, y=d ¼ 1 corresponds to the grain surface. The insert shows schematic deformation of the single grain.
large change of deformation brings single-crystalline grains closer to nonlinear elastic regime. In the limiting case of a 1=d the deformation angle at the surface goes to infinity (and C 66 goes to zero). This case corresponds to the liquid-like phase formed at the grain contact plane. Such liquid-like layer resembles a lubricant that significantly decreases friction between the grains. This fact allows us to suggest that in other physical conditions the increase of the deformation angle caused by the decrease of C 66 will also cause the friction decrease and promote grain boundary sliding.
5. Conclusion We have presented density functional calculations of elastic and electronic properties of doped tetragonal yttria-stabilized zirconium dioxide. The electronic properties of the system does not change significantly with doping and the band gap of 3.7 eV is not affected by addition of metals. Only for titanium oxide dopant some localized unoccupied donor states are observed. Increase of yttrium/dopant concentration from 2.8 to 6.7 mol% reduces the band gap by 5%. Calculated yttrium/dopant concentration dependence of the bulk modulus does not show any systematic changes; thus isotropic compressibility of ZrO2 single crystalline grains does not change upon doping. In contrast, C 66 shear elastic constants are reduced with increasing dopant concentration. Changes of the shear properties in doped ZrO2 may explain discrepancies between calculated and measured tensor of elastic constants of pure ZrO2 . Based on the Wulff construction of the ZrO2 crystallite, we have proposed the model that relates changes of the shear properties and the strain rate. We suggest that the decrease of C 66
[1] F. Wakai, S. Sakaguchi, Y. Matsumoto, Adv. Ceram. Mater. 1 (1986) 259. [2] F. Wakai, N. Kondo, H. Ogawa, T. Nagano, S. Tsukerawa, Mater. Characterization 37 (1996) 331. [3] M. Jimenez-Melendo, A. Dominiguez-Rodriguez, A. Bravo-Leo´n, J. Am. Ceram. Soc. 81 (1998) 2761. [4] A. Nakamura, N. Shibata, N. Morishige, K. Matsunaga, Y. Ikuhara, J. Am. Ceram. Soc. 88 (2005) 938. [5] M.Z. Berbon, T.G. Langdon, Acta Mater. 47 (1999) 2485. [6] D. Go´mez-Garcı´a, J. Martı´nez-Ferna´dez, A. Domı´nguez-Rodrı´guez, J. Castaing, J. Am. Ceram. Soc. 80 (1997) 1668. [7] Y. Ikuhara, Y. Nagai, T. Yamamoto, T. Sakuma, Interface Sci. 7 (1999) 77. [8] D. Go´mez-Garcı´a, C. Lorenzo-Martı´n, A. Munoz-Bernabe´, A. Domı´nguezRodrı´guez, Philos. Mag. 83 (2003) 93. [9] A. Kuwabara, M. Nakano, H. Yoshida, Yu. Ikuhara, T. Sakuma, Acta Mater. 52 (2004) 5563. [10] H. Yoshida, J. Ceram. Soc. Jpn. 114 (2006) 155. [11] J. Mimurada, M. Nakano, K. Sasaki, Y. Ikuhara, T. Sakuma, J. Am. Ceram. Soc. 84 (2001) 1817. [12] Y. Sakka, T.S. Suzuki, T. Matsumoto, K. Morita, K. Hiraga, Y. Moriyoshi, Solid State Ionics 172 (2004) 499. [13] K. Hiraga, J. Ceram. Soc. Jpn. 115 (2007) 395. [14] K. Hiraga, B.-N. Kim, K. Morita, T.S. Suzuki, Y. Sakka, J. Ceram. Soc. Jpn. 113 (2005) 191. [15] F. Wakai, N. Kondo, Y. Shinoda, Curr. Opin. Solid State Mater. Sci. 4 (1999) 461. [16] G. Teufer, Acta Cryst. 15 (1962) 1187. [17] J.K. Dewhurst, J.E. Lowther, Phys. Rev. B 57 (1998) 741. [18] P. Bouvier, E. Djurado, G. Lucazeau, T. Le Bihan, Phys. Rev. B 62 (2000) 8731. [19] A. Bogicevic, C. Wolverton, G.M. Crosbie, E.B. Stechel, Phys. Rev. B 64 (2001) 014106. [20] W. Kohn, L.J. Sham, Phys. Rev. 140 (1965) 1133. [21] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Phys. Rev. B 46 (1992) 6671. [22] J. Soler, E. Artacho, J.D. Gale, A. Garcı´a, J. Junquera, P. Ordejo´n, D. Sa´nchezPortal, J. Phys. Condens. Matter 14 (2002) 2745. [23] N. Troullier, J.L. Martins, Phys. Rev. B 43 (1991) 1993. [24] J.P. Perdew, K. Burke, M. Ernzerhof, Phys. Rev. Lett. 77 (1996) 3865. [25] O.F. Sankey, D.J. Niklewski, Phys. Rev. B 40 (1989) 3979. [26] J. Carrasco, N. Lopez, F. Illas, Phys. Rev. Lett. 93 (2004) 225502. [27] A. Eichler, Phys. Rev. B 64 (2001) 174103. [28] C.J. Howard, R.J. Hill, B.E. Reichert, Acta Crystallogr. B Struct. Sci. 44 (1988) 116. [29] A. Bogicevic, C. Wolverton, Phys. Rev. B 67 (2003) 024106. [30] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flanneryet, Numerical Recipes: The Art of Scientific Computing, third edition, Cambridge University Press, New York, 2007. [31] D.W. McComb, Phys. Rev. B 54 (1996) 7094. [32] F.D. Murnaghan, Proc. Natl. Acad. Sci. 30 (1944) 244. [33] O.H. Nielsen, R.M. Martin, Phys. Rev. Lett. 50 (1983) 697. [34] H. Yao, L. Ouyang, W.-Y. Ching, J. Am. Ceram. Soc. 90 (2007) 3194. [35] J.E. Lowther, Phys. Rev. B 173 (2006) 134110. [36] E.H. Kisi, Ch.J. Howard, J. Am. Ceram. Soc. 81 (1998) 1682. [37] S.-K. Chen, Y. Fang, M. Grimsdithch, M.Z. Li, V. Nevitt, W.M. Robertson, E.S. Zouboulis, J. Am. Ceram. Soc. 74 (1991) 1742. [38] A. Christensen, E.A. Carter, Phys. Rev. B 58 (1998) 8050.
Contents lists available at ScienceDirect
Journal of Physics and Chemistry of Solids journal homepage: www.elsevier.com/locate/jpcs
Influence of elastic properties on superplasticity in doped yttria-stabilized zirconia Yuriy Natanzon a, Marek Boniecki b, Zbigniew Łodziana a,c, a b c
´ ski Institute of Nuclear Physics, Polish Academy of Sciences, ul. Radzikowskiego 152, 31-342 Krako ´w, Poland Henryk Niewodniczan ´lczyn ´ ska 133, 01-919 Warsaw, Poland Institute of Electronic Materials Technology, ul. Wo ¨ bendorf, Switzerland Department of Environment, Energy and Mobility, EMPA, 8600 Du
a r t i c l e in fo
abstract
Article history: Received 14 July 2008 Received in revised form 12 August 2008 Accepted 17 August 2008
We report the first principles calculations of elastic and electronic properties of yttria-stabilized tetragonal zirconium dioxide (YZP) doped with GeO2 , TiO2 and SiO2 . Electronic structure and isotropic elastic properties of YZP do not change upon addition of dopants. Addition of dopants affects the shear C 66 elastic constant that decreases with the increasing dopant concentration. A simple model that connects elastic softening to enhancement of superplasticity in doped fine-grained zirconia ceramics is proposed. & 2008 Elsevier Ltd. All rights reserved.
Keywords: A. Oxides C. Ab initio calculations D. Elastic properties D. Electronic structure
1. Introduction For the first time superplasticity in tetragonal yttria-stabilized zirconia (YZP) was observed over than 20 years ago by Wakai et al. [1,2]. Since that time superplastic properties of YZP have been an exciting field of research (see review of Jimenez-Melendo et al. [3]). At present there is no commonly accepted explanation of this phenomenon. One of the most accepted explanation of superplasticity in YZP is grain boundary sliding related to deformation. It was confirmed by number of microstructural observations [2,4]. Superplasticity is driven by the diffusion of Zr IV ions at the grain boundary (or so-called diffusion creep) [5]. Yttrium dopants segregate toward the surface on grain boundary and that enhances Zr diffusion [6–8]. Such diffusion contributes to the increase of superplastic strain rate and facilitates intergranular sliding during deformation [3]. Go´mez-Garcı´a et al. [8] have shown that segregation of Y cations cause electrostatic potential difference at the grain boundary. Intergranular sliding is then supported by the electrostatic interaction between the grains. It is well known that superplastic properties of ZrO2 are drastically changed by addition of small amount of metal oxides [9]. Some dopants can significantly increase the superplastic strain rate of zirconia to allow the elongation of up to 1000% of its
initial length [9,10]. This behavior is confirmed for such dopants as GeO2 [9–11], TiO2 [9–12], SiO2 [4,13,14], MgO [11,12] and Al2 O3 [13,14] as well as their combinations. It was suggested that changes in superplasticity are due to the formation of glassy thin film of dopant, which acts as a lubricant for intergranular sliding [2]. However, the effect is still observed when the intermediate glass layer is not present but instead the dopant segregation near the grain boundary is observed [15]. Another reason for the effect could be the related changes in ionic bonding strength of cation–oxygen bonds in the vicinity of grain boundary [9,10]. Unfortunately the existing models do not provide clear microscopic picture of the superplasticity in ZrO2 related to metal oxide doping. Especially there is lack of atomic level insight in the literature [9,10]. In order to gain more detailed microscopic insight into phenomenon of superplasticity in doped YZP we have carried out a number of first principles calculations of elastic properties of yttria-stabilized tetragonal zirconia with GeO2 , TiO2 and SiO2 dopant metal oxides. Our calculations indicate the correlation between changes in shear elastic constants and superplasticity. We show that some of the shear elastic properties depend on Y/dopant concentration while the isotropic properties, like bulk modulus do not.
2. Model and method of calculation Corresponding author at: Henryk Niewodniczan´ski Institute of Nuclear Physics,
Polish Academy of Sciences, ul. Radzikowskiego 152, 31-342 Krako´w, Poland. E-mail addresses: [email protected] (Y. Natanzon), Marek.Boniecki@ itme.edu.pl (M. Boniecki), [email protected] (Z. Łodziana). 0022-3697/$ - see front matter & 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.jpcs.2008.08.014
Depending on the external thermodynamic conditions three phases of ZrO2 can be stable. Phase transition sequence of ZrO2 includes high temperature cubic phase ðFm3mÞ which exists
ARTICLE IN PRESS 16
Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
between T ¼ 2650–2990 K and tetragonal phase (P42 =nmc) at T ¼ 1400–2650 K. Monoclinic phase (P21 =c) exists at the temperatures below 1400 K [16–18]. However, high temperature phases can be stabilized with various metal oxides like MgO, CaO, CeO2 and Y2 O3 and exist at lower temperatures. As seen from phase diagram of YZP in Ref. [19], tetragonal phase of zirconia exists for yttria concentrations below 6 mol% and temperature range of about 600–2000 K. In this paper we consider tetragonal phase of zirconia with the concentration of yttria below 6.7 mol%. Our present calculations were performed within densityfunctional theory (DFT) approach [20] with generalized-gradient approximation (GGA) [21] for exchange-correlation functional. The linear combination of atomic orbitals (LCAO) basis set was used, as implemented in SIESTA program [22]. The implementation of DFT used in the present approach uses specific techniques to speed up the calculations and balance between the desired accuracy and calculation time. Only the valence electrons which take part in bond formation are considered. The other electrons (also called ‘‘core electrons’’) are replaced by effective norm-conserving potentials with Troullier and Martins parametrization [23]. In this calculation we have used GGA pseudopotentials in Perdew–Burke–Ernzerhof parametrization [24]. In order to achieve better accuracy we have also treated 4s and 4p core electrons of Zr and Y atoms as valence 2 states. Thus valence electron configuration was 4s2 4p6 5s2 4d for 1 2 2 6 2 2 6 2 2 2 2 Zr, 4s 4p 5s 4d for Y, 3s 3p 4s 3d for Ti, 3s 3p for Si, 4s 4p2 for Ge, 3s2 3p1 for Al and 3s2 for Mg. Secondly, all the integrals in SIESTA are calculated on the finite real space grid. Plane-wave cutoff of 476 eV (350 Ry) was used in this calculation and a k point grid with a step of 0.03 A˚. The wave functions are expanded into strictly localized atomic orbitals, which consist of the product of numerical radial function and spherical harmonics [22,25]. Double-z with polarization (DZP) basis set was used for the valence electrons of all atoms and double-z (DZ) for the semicore 4s and 4p electrons of Zr and Y. All the pseudopotentials and basis sets were properly tested on the corresponding oxides to reproduce good experimental values of lattice constants and bulk modulus. Additionally, we have calculated total energies for different phases of ZrO2 in order to check if the energy sequence is correct and states Emonoclinic oEtetragonal oEcubic [19]. In the present studies the supercell technique was used to vary the concentration of yttrium and other dopants [26]. At the first step experimental unit cell parameters of ZrO2 were optimized and bulk modulus was calculated in order to assess the accuracy of calculations. The calculated lattice constants and bulk modulus are shown in Table 1. Calculated bulk modulus agree within 5% with this measured by XRD [18] and the lattice parameters are in good agreement with previous theoretical [17,27] and experimental results [18,28]. To study stabilized ZrO2 two formula units of zirconia were replaced with one formula unit of yttria. In order to keep the
Table 1 Comparison of calculated lattice constants and elastic properties of pure ZrO2 with experimental measurements and other calculations Lattice constants (A˚)
This work (GGA) Ref. [17] (LDA) Ref. [27] (GGA) Ref. [28] (neutron diffraction) Ref. [18] (XRD)
a¼b
c
3.61 3.66 3.645 3.60 3.59
5.20 5.30 5.29 5.18 5.18
Bulk modulus (GPa)
209.5 200 – – 198
LDA stands for local density approximation, GGA for generalized gradients approximation and X-ray diffraction is abbreviated as XRD.
correct stoichiometry for each two Y dopant atoms one oxygen should be removed. Relative position of Y atoms and oxygen vacancy was chosen according to results of Bogicevic et al. [19,29] and Eichler [27]. Additional tests confirm that Y atoms located in ða; bÞ plane at the distance of 3.5 A˚ and the oxygen being the first neighbor of dopants provide the lowest energy configuration. The optimized unit cell of zirconia was multiplied in selected directions to obtain desired size of the supercell. The concentrations of Y and dopant oxides were the same depending on the size of the supercell as each supercell contained two Y dopants and one metal oxide dopant (the smaller the cell the larger the concentration and vice versa). In fact, the concentration of 6.7 mol% goes actually beyond the stability of the tetragonal phase and corresponds already to the cubic ZrO2 . However, the concentration of dopants increases towards the surface [6] and such a large concentration was used for the comparison purposes. For each dopant concentration and for each supercell size we have performed full structural relaxation (the supercell size and atomic positions within the supercell) to obtain the ground state configuration of the system. The conjugate-gradients method was used for such an optimization [30].
3. Results 3.1. Electronic structure Tetragonal zirconium dioxide is an insulator. Its measured band gap is Eg ¼ 4:2 eV [31], the band gap calculated by Eichler in GGA is Eg ¼ 4:21 eV [27]. In the present calculation the band gap is Eg ¼ 3:7 eV. Slight underestimation is typical for the approximations used here. Stabilization of zirconia with yttrium oxide as well as further aliovalent doping with metal oxides does not significantly influence the electronic structure of ZrO2 [19]. In the present calculation isovalent dopants were studied. Calculated electronic density of states of pure and doped zirconia is shown in Fig. 1 and does not depend on the type of the dopant. Only for titanium oxide dopant some localized unoccupied donor states are observed. The valence band mainly consists of oxygen 2p orbitals and the unoccupied conduction band include empty states of metal cation. The DOS for yttria concentration of 2.8% and 6.7% is presented in the bottom part of Fig. 1. The band gap shrinks by 0.16 eV (5%) for high dopant concentration. Results presented above indicate that macroscopic electronic properties of ZrO2 do not change after the dopant addition. 3.2. Elastic properties The phenomena of superplasticity is related to changes of the macroscopic elastic properties of fine-grained ZrO2 . As it is generally attributed to the grain boundary sliding, fewer attention is given to the elastic properties of the single crystal grains, which also affect the plastic properties of YZP. In this section we present the results that show what is the influence of doping on the elasticity of bulk ZrO2 crystal. The basic elastic properties include volume compressibility (bulk modulus) and tensor of the elastic constants. Below we present how these properties depend on concentration and the type of dopant. 3.2.1. Bulk modulus Bulk moduli were obtained by series of the calculations where uniform lattice deformation in the range of 3% to 3% of its initial volume with the step of 0.5% was applied. For each calculation
ARTICLE IN PRESS Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
17
215
Bulk modulus (GPa)
210
205
200 No dopant Ti Ge Si
195
190 1.5
2
2.5
3
4 3.5 4.5 5 5.5 Yttria concentration (mol%)
6
6.5
7
Fig. 2. Y/dopant concentration dependence of bulk modulus of YZP. Filled circles correspond to undoped YZP and empty symbols to YZP with Ti, Ge or Si dopants. Solid line is a linear fit of the average values for each concentration.
Table 2 Calculated elastic constants of pure ZrO2 in comparison with experimental measurements and other calculations C 11 C 12 C 13 C 33 C 44 C 66 This work (GGA) Ref. [35] (LDA) Ref. [36] (neutron diffraction, 12 mol% Ce-doped ZrO2 ) Ref. [37] (Brillouin scattering, clear monoclinic ZrO2 )
293 248 111 382 221 72 327 100 62 340 33 160
385 346 264 325
51 187 42 167 59 64 66 95
All the values are given in GPa.
analysis of the components of the elastic constants tensor can reveal more detailed information about system response to nonuniform strains that arise during grain sliding. In the next section we analyze the shear properties of YZP. Fig. 1. Total density of states (DOS) for 2.8 mol% doped tetragonal Y-stabilized zirconia (YZP): DOS for pure ZrO2 is shown at the bottom panel, the undoped YZP is second from the bottom (thick line is for Y concentration of 6.7mol%). Upper panels are for YZP with Si, Ge and Ti dopants.
structural relaxation of atomic positions was performed. To obtain the value of bulk modulus the calculated dependence of the total energy E of the system on the supercell volume V was fitted by Murnaghan equation of state [32]. The dependence of the bulk modulus on Y/dopant concentration is presented in Fig. 2. From the linear fit of the average values of bulk modulus one can see that no systematic changes of the bulk modulus with Y/dopant concentration are observed. The values of the bulk modulus are slightly lower than that of pure ZrO2 (see Table 1). The fluctuations are related to the accuracy of calculations (relaxation of large supercells with more than 100 atoms) and as well as to the errors in the Murnaghan equation fit. Significant change of bulk modulus for the concentration of 6.7 mol% is related to the phase transition. For the concentrations above 6 mol% the tetragonal phase is unstable with respect to the cubic one [19]. The observed fluctuations (see Fig. 2) reflect this instability. Weak concentration dependence of the bulk modulus is in good agreement with models of the superplastic behavior of doped ZrO2 which indicate that the macroscopic grain deformation is not related to the increase of superplasticity [3]. However, elastic response of the system is the tensor quantity. Detailed
3.2.2. Shear elastic constants Elastic constants of pure and doped ZrO2 were calculated by series of finite deformations of the supercell as described in the work of Nielsen and Martin [33] and Hongzhi Yao et al. [34]. For each constant and each dopant concentration the appropriate deformation of 1% was applied and the components of the stress tensor were calculated as the derivatives of Kohn–Sham total energy with respect to the strain tensor. Thus elastic constants were calculated by dividing the components of the stress tensor by the components of strain tensor: C nn ¼ sn =n , where n ¼ 4; 5; 6 and i ¼ 0:01 only for the component for i ¼ n. To minimize the numerical errors, each elastic constant was obtained as the average of two calculations with the values of stress of þ1% and 1%. For each deformation a full relaxation of the atomic positions in the supercell was performed. The elastic constant tensor for pure ZrO2 is presented in Table 2. The calculated values of elastic constants are in good agreement to available theoretical results. Some of the components show deviation from the experimental data; however, the measurements are usually performed on stabilized zirconia [36] or on the stable monoclinic phase with the approximation of elastic constants temperature dependence to the tetragonal phase [37]. As will be shown below significant difference between theoretical and experimental values of C 66 elastic constant is due to the dopant or stabilizer addition to ZrO2 .
ARTICLE IN PRESS Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
65
190
60
180 C66 (GPa)
C44 (GPa)
18
55 50 45
No dopant Ti Ge Si
170 160 No dopant Ti Ge Si
150 140
40 1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 Yttria concentration (mol%)
1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 7 Yttria concentration (mol%)
Fig. 3. Y/dopant concentration dependence of C 44 elastic constant (a) and C 66 elastic constant (b). Filled circles correspond to undoped YZP and empty symbols are for YZP with Ti, Ge or Si dopant. Solid line is a linear fit of the average values for each concentration.
Calculated dopant concentration dependence of C 44 and C 66 elastic constants for doped Y-stabilized ZrO2 are presented in Fig. 3. As one can see from the linear fit of the average values there are no systematic changes in C 44 constant (see Fig. 3(a)) and systematic reduction of the C 66 elastic constant with increasing dopant concentration is observed (see Fig. 3(b)). The largest reduction is observed for Si dopant where for 4.8 mol% increase of dopant concentration C 66 drops by 20%. For yttria-stabilized ZrO2 reduction of shear elastic constants with increasing concentration of Y is also observed, and this change is below 10%.
4. Discussion It is well established that the concentration of dopant increases towards the grain surface in doped YZP [6–8]. According to results presented above this phenomena would affect shear modulus in the grain boundary region. The strain rate of YZP dependence on the inverse of the shear modulus is already predicted from a simple phenomenological model [3,8]. The empirical equation for the strain rate is given by [3]: 2 2 s , G
q ADGb b ¼ kT d qt
(1)
where A is the dimensionless constant, D the grain diffusion coefficient, b the Burgers vector, k the Boltzmann constant, T the temperature, d the grain size, s the applied stress and G the shear modulus. The shear modulus, which is the ratio of shear stress and the shear strain by definition equals the combination shear elastic constants C 66 or C 44 (see Table 2) depending on the crystalline direction. The present calculations are in a good agreement with this relation and provide microscopic insight into the phenomenon. Below we propose a simple model that relates dopant concentration with the shear elastic constants and flow stress. Intergrain sliding in YZP occurs along crystalline facets that are exposed in the crystallites. The Wulff construction of the tetragonal nano-crystallite of ZrO2 is presented in Fig. 4(a). The construction is based on the surface energy calculations by Christensen and Carter [38]. It is clear that the (111) facet will dominate the surface, two other facets that are present on the surface are (1 0 0) and (0 0 1). The (0 0 1) facet is directly related to C 66 elastic constant, while shear deformation of the (111) surface is linear combination of shear related to C 44 (23 of total shear) and C 66 (13 of total shear). As the C 44 shear modulus does not change upon doping, decrease of C 66 will contribute to the decrease of the apparent shear modulus G increasing the strain rate as can be seen from Eq. (1). Below we focus on effects solely related to C 66 shear
Fig. 4. (a) Wulff construction of the zirconia grain shape: the (111) surface is marked as dark gray (blue) and it is dominating one, the (0 0 1) surface at the top (and bottom) is marked as gray (red) and (1 0 0) facet is light gray; (b) a schematic model showing the sliding of the ZrO2 grains. The d denotes the size of the layer where the dopant concentration varies, arrows show the direction of the applied force F and the direction of sliding. The grain surface layers undergo shear deformation during sliding.
modulus, such considerations provide the upper boundary of the influence of elastic softening on the superplastic flow. For simplicity we consider zirconia grains as cuboidal objects sliding along (0 0 1) contact facets following the idea of Go´mezGarcı´a et al. [8]. We assume that in the layer of thickness d from the surface the concentration of dopant increases linearly towards the surface. This assumption together with the result presented in Fig. 3(b) provide the expression for the concentration dependence of C 66 elastic constant: C 66 ¼ C 066 ð1 ayÞ,
(2)
C 066
is the elastic constant for the bulk, a is the coefficient where which depends on the type of dopant. The y is the distance from the bottom of layer d to the contact surface plane; 0pypd, and ayp1. Schematic picture of the model is presented in Fig. 4(b). During sliding grains undergo a shear deformation on the contact surfaces. The work applied on the contact boundary is dW ¼ 2F dy, where F is a force and dy, shear deformation angle. In the linear regime the deformation is proportional to the applied 2 force, thus dW ¼ C 66 dy , and the expression for the deformation angle reads
dy ¼
2F C 066 ð1 ayÞ
.
(3)
With assumptions of the linear change of dopant concentration the deformation of ZrO2 grain in the vicinity of the crystalline contact surface is presented in Fig. 5. With our secure assumption of linear change of concentration and for the linear elastic regime the deformation angle is already doubled at the surface. Such a
ARTICLE IN PRESS Y. Natanzon et al. / Journal of Physics and Chemistry of Solids 70 (2009) 15–19
19
elastic constants induces increase of the strain rate at the grain surfaces and promotes intergranular sliding by the decrease of intergranular friction. The increase of the strain rate in doped tetragonal phase of Y-stabilized ZrO2 is proportional to the inverse of shear modulus that decrease upon doping. This effect contributes to overall superplasticity besides the electrostatic and kinetic effects.
Acknowledgment The authors acknowledge CPU time allocation at Interdisciplinary Center for Mathematical and Computer Modelling of Warsaw University (Grant no. G28-22). This work was supported by Ministry of Science and Higher Education of Poland under Grant no. 3 T08D 013 29 during the years 2005–2008. References
Fig. 5. Change of the shear deformation angle (y) with respect to depth inside the grain caused by the increase of dopant concentration. y0 is the deformation for constant bulk concentration of dopant. Origin of the x-axis corresponds to deformation in the bulk, y=d ¼ 1 corresponds to the grain surface. The insert shows schematic deformation of the single grain.
large change of deformation brings single-crystalline grains closer to nonlinear elastic regime. In the limiting case of a 1=d the deformation angle at the surface goes to infinity (and C 66 goes to zero). This case corresponds to the liquid-like phase formed at the grain contact plane. Such liquid-like layer resembles a lubricant that significantly decreases friction between the grains. This fact allows us to suggest that in other physical conditions the increase of the deformation angle caused by the decrease of C 66 will also cause the friction decrease and promote grain boundary sliding.
5. Conclusion We have presented density functional calculations of elastic and electronic properties of doped tetragonal yttria-stabilized zirconium dioxide. The electronic properties of the system does not change significantly with doping and the band gap of 3.7 eV is not affected by addition of metals. Only for titanium oxide dopant some localized unoccupied donor states are observed. Increase of yttrium/dopant concentration from 2.8 to 6.7 mol% reduces the band gap by 5%. Calculated yttrium/dopant concentration dependence of the bulk modulus does not show any systematic changes; thus isotropic compressibility of ZrO2 single crystalline grains does not change upon doping. In contrast, C 66 shear elastic constants are reduced with increasing dopant concentration. Changes of the shear properties in doped ZrO2 may explain discrepancies between calculated and measured tensor of elastic constants of pure ZrO2 . Based on the Wulff construction of the ZrO2 crystallite, we have proposed the model that relates changes of the shear properties and the strain rate. We suggest that the decrease of C 66
[1] F. Wakai, S. Sakaguchi, Y. Matsumoto, Adv. Ceram. Mater. 1 (1986) 259. [2] F. Wakai, N. Kondo, H. Ogawa, T. Nagano, S. Tsukerawa, Mater. Characterization 37 (1996) 331. [3] M. Jimenez-Melendo, A. Dominiguez-Rodriguez, A. Bravo-Leo´n, J. Am. Ceram. Soc. 81 (1998) 2761. [4] A. Nakamura, N. Shibata, N. Morishige, K. Matsunaga, Y. Ikuhara, J. Am. Ceram. Soc. 88 (2005) 938. [5] M.Z. Berbon, T.G. Langdon, Acta Mater. 47 (1999) 2485. [6] D. Go´mez-Garcı´a, J. Martı´nez-Ferna´dez, A. Domı´nguez-Rodrı´guez, J. Castaing, J. Am. Ceram. Soc. 80 (1997) 1668. [7] Y. Ikuhara, Y. Nagai, T. Yamamoto, T. Sakuma, Interface Sci. 7 (1999) 77. [8] D. Go´mez-Garcı´a, C. Lorenzo-Martı´n, A. Munoz-Bernabe´, A. Domı´nguezRodrı´guez, Philos. Mag. 83 (2003) 93. [9] A. Kuwabara, M. Nakano, H. Yoshida, Yu. Ikuhara, T. Sakuma, Acta Mater. 52 (2004) 5563. [10] H. Yoshida, J. Ceram. Soc. Jpn. 114 (2006) 155. [11] J. Mimurada, M. Nakano, K. Sasaki, Y. Ikuhara, T. Sakuma, J. Am. Ceram. Soc. 84 (2001) 1817. [12] Y. Sakka, T.S. Suzuki, T. Matsumoto, K. Morita, K. Hiraga, Y. Moriyoshi, Solid State Ionics 172 (2004) 499. [13] K. Hiraga, J. Ceram. Soc. Jpn. 115 (2007) 395. [14] K. Hiraga, B.-N. Kim, K. Morita, T.S. Suzuki, Y. Sakka, J. Ceram. Soc. Jpn. 113 (2005) 191. [15] F. Wakai, N. Kondo, Y. Shinoda, Curr. Opin. Solid State Mater. Sci. 4 (1999) 461. [16] G. Teufer, Acta Cryst. 15 (1962) 1187. [17] J.K. Dewhurst, J.E. Lowther, Phys. Rev. B 57 (1998) 741. [18] P. Bouvier, E. Djurado, G. Lucazeau, T. Le Bihan, Phys. Rev. B 62 (2000) 8731. [19] A. Bogicevic, C. Wolverton, G.M. Crosbie, E.B. Stechel, Phys. Rev. B 64 (2001) 014106. [20] W. Kohn, L.J. Sham, Phys. Rev. 140 (1965) 1133. [21] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Phys. Rev. B 46 (1992) 6671. [22] J. Soler, E. Artacho, J.D. Gale, A. Garcı´a, J. Junquera, P. Ordejo´n, D. Sa´nchezPortal, J. Phys. Condens. Matter 14 (2002) 2745. [23] N. Troullier, J.L. Martins, Phys. Rev. B 43 (1991) 1993. [24] J.P. Perdew, K. Burke, M. Ernzerhof, Phys. Rev. Lett. 77 (1996) 3865. [25] O.F. Sankey, D.J. Niklewski, Phys. Rev. B 40 (1989) 3979. [26] J. Carrasco, N. Lopez, F. Illas, Phys. Rev. Lett. 93 (2004) 225502. [27] A. Eichler, Phys. Rev. B 64 (2001) 174103. [28] C.J. Howard, R.J. Hill, B.E. Reichert, Acta Crystallogr. B Struct. Sci. 44 (1988) 116. [29] A. Bogicevic, C. Wolverton, Phys. Rev. B 67 (2003) 024106. [30] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flanneryet, Numerical Recipes: The Art of Scientific Computing, third edition, Cambridge University Press, New York, 2007. [31] D.W. McComb, Phys. Rev. B 54 (1996) 7094. [32] F.D. Murnaghan, Proc. Natl. Acad. Sci. 30 (1944) 244. [33] O.H. Nielsen, R.M. Martin, Phys. Rev. Lett. 50 (1983) 697. [34] H. Yao, L. Ouyang, W.-Y. Ching, J. Am. Ceram. Soc. 90 (2007) 3194. [35] J.E. Lowther, Phys. Rev. B 173 (2006) 134110. [36] E.H. Kisi, Ch.J. Howard, J. Am. Ceram. Soc. 81 (1998) 1682. [37] S.-K. Chen, Y. Fang, M. Grimsdithch, M.Z. Li, V. Nevitt, W.M. Robertson, E.S. Zouboulis, J. Am. Ceram. Soc. 74 (1991) 1742. [38] A. Christensen, E.A. Carter, Phys. Rev. B 58 (1998) 8050.