Optimum strain state for oxygen diffusion in yttria-stabilised zirconia

Optimum strain state for oxygen diffusion in yttria-stabilised zirconia

Solid State Ionics 190 (2011) 75–81 Contents lists available at ScienceDirect Solid State Ionics j o u r n a l h o m e p a g e : w w w. e l s ev i e...

1MB Sizes 2 Downloads 60 Views

Solid State Ionics 190 (2011) 75–81

Contents lists available at ScienceDirect

Solid State Ionics j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / s s i

Optimum strain state for oxygen diffusion in yttria-stabilised zirconia Wakako Araki ⁎, Yoshio Arai Dept. of Mechanical Engineering, Saitama University, 255 Shimo-okubo, Sakura-ku, Saitama 338 8570, Japan

a r t i c l e

i n f o

Article history: Received 6 December 2010 Received in revised form 21 March 2011 Accepted 23 March 2011 Available online 17 April 2011 Keywords: Strain Ionic conductivity Oxygen diffusion Stabilised zirconia Solid oxide fuel cell

a b s t r a c t The strain effect of the oxygen diffusion in the yttria-stabilised zirconia in the present study was investigated by means of a molecular dynamics simulation. The simulation was conducted for various parameters such as the strain, temperature, yttria concentration, and potential parameter, and in addition, the biaxial and hydrostatic cases were considered. For a uniaxial strain, the oxygen diffusion was enhanced in the tensile direction whereas it was hindered in the compressive direction. The maximum improvement was achieved for a smaller strain at a lower temperature and also with a lower yttria concentration. For the biaxial and hydrostatic strains, the total diffusion coefficient was enhanced simply as a result of the enhancement when using the uniaxial strain in each direction. The detailed deformation analysis reveals that the optimum strain state for the highest oxygen diffusion in the tensile direction can be obtained when the oxygen ion is largely displaced in the fluorite lattice structure. © 2011 Elsevier B.V. All rights reserved.

1. Introduction Zirconia stabilised with trivalent rare earth elements such as yttria and scandia has become very attractive for use as an electrolyte material for solid oxide fuel cell (SOFC) because of its excellent properties, such as its high oxygen-ionic conductivity, great chemical stability and mechanical strength. There has been a lot of research over the past several decades devoted to improving the ionic conductivity of the electrolytes including the stabilized zirconia that ensure the higher power performance of the SOFC [1–3]. The effect of various dopants on the conductivity of stabilised zirconia has been discussed taking into account the strain field caused by the dopant [4]. Multilayered electrolyte systems have recently been proposed and the colossal improvements in the ionic conductivity, i.e., one or two orders, have been reported [5–8]. In a multilayered system, each electrolyte layer with a thickness of a few nano-meters is extensively strained in the biaxial directions, i.e., larger than 5%. The improvement in conductivity has been so far attributed to various causes such as the strain, defects, dislocations, interfaces, and the grain boundaries. However, the mechanism has not been fully explained. An investigation of uniaxial strain effect on conductivity should be essential in order to understand effects of biaxial and hydrostatic strains. The effect of uniaxial strain on the oxygen conductivity has been investigated by means of experiments and simulations [9–12]. A clear improvement in oxygen migration has been observed in the tensile direction while a decrease in conductivity has been obtained in the compressive direction. The improvement was considered to be related

⁎ Corresponding author. Tel.: + 81 48 858 3435; fax: + 81 48 856 2577. E-mail address: [email protected] (W. Araki). 0167-2738/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.ssi.2011.03.012

to the strain energy [10,11], although the mechanism for the improvement caused by the strain has not yet been clarified. Therefore, the strain effect on the oxygen migration in the yttriastabilised zirconia was investigated for this study by conducting a molecular dynamics simulation. The simulation was conducted using various parameters such as the strain, temperature, yttria concentration, and potential parameter. The oxygen diffusion property as well as the stress–strain curve under various conditions has been analysed, and the deformation process such as the lattice deformation and the oxygen displacement has also been closely examined. The optimum strain state for the oxygen migration has been clarified based on these simulation results. 2. Molecular dynamics simulation In the present molecular dynamics (MD) simulation, the oxygen diffusions in the 4, 8, and 14 mol%-yttria-doped zirconia (4YSZ, 8YSZ, and 14YSZ) that were subjected to various stresses were analyzed. The fluorite crystal structures of 4YSZ, 8YSZ, and 14YSZ were modelled by using an MD program (Materials Explorer 5.0, Fujitsu). Each model consisted of 10 × 10× 10 unit cells, where Y3+ ions were randomly dispersed. The structure was first relaxed at 300 K for 10 ps, followed by another relaxation process at a simulated temperature for 10 ps, where the time step was 0.5 fs. The simulation temperatures were 973, 1273, and 2000 K. The stress was then instantly applied to the MD model and maintained for 100 ps at 1273 and 2000 K and for 500 ps at 973 K. The applied stress was mostly uniaxial in the [100], [110], or [111] directions and varied from 0 to the fracture stress level. The oxygen diffusion behaviour was investigated under various stress conditions. In addition, the oxygen diffusion was investigated under the biaxial and hydrostatic stress conditions.

76

W. Araki, Y. Arai / Solid State Ionics 190 (2011) 75–81

3. Simulation result of 8YSZ at 1273 K

Fig. 1. Stress–strain curve with apparent elastic modulus of 8YSZ at 1273 K.

The simulation was carried out by using the Gear's integration method. The Born–Mayer–Huggins interatomic potential was mainly used with the parameters given by Schelling et al. [13], which has been commonly used [14]. The temperature was controlled by using the velocity scaling method, and the Coulomb interaction was calculated by using the Ewald method. The periodic boundary condition was taken into consideration, and the shape of the MD cell was assumed to be cubic during the relaxation processes and deformable during the tensile stress loading. The number of atoms, the temperature, and the pressure were kept constant during the simulation (NTP ensemble). The movement of oxygen ions is evaluated by using the mean square displacement, LMSD, calculated by using [15] LMSD = ð1 = NÞΣ½r ðt Þ–rð0Þ

2

  i = 1e N ;

ð1Þ

where r is the displacement of each ion, N is a number of ions, and t is the time. The diffusion coefficient of the oxygen ion, D, is obtained from the LMSD of the oxygen ion according to the following equation: D = lim LMSD = 6t ðt→∞Þ;

ð2Þ

The activation energy for the oxygen diffusion, ΔH, is given by: D = B expð−ΔH = RT Þ;

ð3Þ

where B is the pre-exponential factor and R is the gas constant. In this study, the total diffusion coefficient and also the directional diffusion coefficients in the [100], [010], and [001] directions are calculated.

Fig. 1 shows the stress–strain curve in the [100] direction of 8YSZ at 1273 K from −4 to 12 GPa. The relationship between the stress and the strain is almost linear from 0 to 1.8% of the strain, and then it clearly shows the nonlinearity, followed by a fracture at a strain of 3.7%, i.e., the applied stress of 12 GPa. Fig. 1 also shows the apparent elastic modulus in the [100] direction, where the modulus is simply derived from the ratio of the increases in the stress to the strain. The modulus is about 450 GPa and gradually decreases, and then it starts to show a significant decrease when the strain is around 1.5%. The diffusion property is then calculated for each strain state. Fig. 2(a) shows the total mean square displacement for 100 ps, calculated at different strains. The displacement slightly increases with a strain up to 1.8%, but it significantly decreases with a strain of 2.6%. Fig. 2(b) shows the total diffusion coefficient at different strains calculated from the slope of the mean square displacement. The diffusion coefficient is about 2.8 × 108 m2s− 1 without any strain and increases with a strain up to about 1.5%, but then is followed by a plummet, and has a low coefficient of about 2.0 m2s− 1 until the fracture. Fig. 3(a)–(c) shows the directional mean square displacement in the [100], [010], and [001] directions of 8YSZ at 1273 K for 100 ps, calculated for different strains. The displacement in the [100] direction, i.e., the tensile direction, significantly increases with a strain up to 1.8%, and then it plummets with a strain of 2.6% as shown in Fig. 3(a). In other vertical directions, i.e., in the [010] and [001] directions, however, the displacements are almost constant up to a strain of 1.8% followed by a huge decrease with the strain of 2.6%, as shown in Fig. 3(b) and (c). Fig. 4 summarizes the directional diffusion coefficients in the [100], [010], and [001] directions of 8YSZ at 1273 K calculated for different strains. The diffusion coefficient in the [100] direction clearly increases with a strain up to 1.8%, while the coefficients in the other vertical directions are almost constant or slightly decrease. The maximum improvement for 8YSZ at 1273 K is about 44%, where the maximum improvement rate is evaluated from the strain state which provides the largest directional mean square displacement in the tensile direction. With a further increase in the strain, all the directional coefficients show a huge decrease. In the uniaxial compression case, the diffusion coefficient in the [100] direction, i.e., in the compressive direction, decreases whereas the ones in the other vertical directions increase. The diffusion property in each direction contributes to the total diffusion shown in Fig. 2, so that the total diffusion coefficient has a maximum value with a strain of about 1.5% but the improvement is relatively small. In addition, the comparison of the directional diffusion coefficient in Fig. 4 with the stress–strain curve in Fig. 1 reveals that the diffusion coefficient in the tensile direction increases while the elastic modulus gradually decreases, but all the coefficients greatly decrease at around the beginning of the significant decrease of the modulus. The correlation

Fig. 2. Strain effect on total oxygen diffusion property of 8YSZ at 1273 K: (a) total mean square displacement, (b) total diffusion coefficient.

W. Araki, Y. Arai / Solid State Ionics 190 (2011) 75–81

77

Fig. 3. Strain effect on directional mean square displacement of 8YSZ at 1273 K: (a) [100], (b) [010], (c) [001].

between the diffusion enhancement and the decreasing behaviour of the modulus will be discussed in chapter 5, allowing for the deformation process. The present simulation result has demonstrated, therefore, that the oxygen diffusion can be enhanced in the tensile direction while it can be deteriorated in the compressive direction. It also indicates that there could be an optimum strain state for the oxygen diffusion in the yttria-stabilised zirconia, as it has been suggested in the reference [16]. 4. Consideration with various parameters The strain effect on the oxygen diffusion has been investigated under various conditions in order to clarify the improvement mechanism as well as to obtain the optimum strain state for the oxygen diffusion. 4.1. Temperature The temperature effect on the improvement of the directional diffusion coefficient has been taken into consideration. Figs. 5 and 6 show the stress–strain curve and the directional diffusion coefficients of 8YSZ at 973 K and 2000 K, respectively. For the case at 973 K, the diffusion coefficient cannot be calculated when the strain is larger than 3.4%, because the phase transformation from the cubic phase to the orthogonal I phase occurs just after the deformation. (It should be noted here that, in the molecular dynamics simulation with the present potential parameter set, the fluorite phase transforms to the orthogonal I phase, which has a similar structure to the monoclinic phase [14]). The apparent elastic modulus is generally lower and the diffusion coefficient is higher at a higher temperature, as is well known. At both temperatures, the oxygen diffusion is clearly enhanced in the tensile

Fig. 4. Strain effect on directional diffusion coefficient of 8YSZ at 1273 K: (a) [100], (b) [010], (c) [001].

direction and it is at the maximum value at around the beginning of the huge decrease of the modulus. The improvement of the diffusion coefficient in the tensile direction occurs at a smaller strain at a lower temperature. The maximum 50% improvement is achieved at about 1.3% of the strain at 973 K while the 36% one is achieved at about 2.6% of the strain at 2000 K. The improvement is larger at the lower temperature, but they are not that different.

4.2. Potential parameter The simulation has also been conducted with a different potential parameter set. Fig. 7 shows the stress–strain curve and the diffusion coefficient of 8YSZ at 1273 K calculated using the potential parameter set given by Brinkman et al. [17]. The Brinkman's parameter generally gives a slightly larger lattice constant for 8YSZ, higher elastic moduli, higher diffusion coefficients with a lower activation energy, and more cubic-statbilized structures than the Schelling's one, as has been reported in the references [11]. As for the diffusion coefficient, the

Fig. 5. Simulation results of 8YSZ at 973 K: (a) stress–strain curve, (b) directional diffusion coefficient.

78

W. Araki, Y. Arai / Solid State Ionics 190 (2011) 75–81

Fig. 6. Simulation results of 8YSZ at 2000 K: (a) stress–strain curve, (b) directional diffusion coefficient.

Fig. 8. Simulation results of 4YSZ at 2000 K: (a) stress–strain curve, (b) directional diffusion coefficient.

Schelling's parameter seems to give a better agreement with the experimental data, i.e. about 3 × 10− 7 cm2/s at 1273 K [18,19]. Although the above properties are different from the results for the Schelling's parameter shown in the previous chapter, a similar strain effect on the diffusion was observed. Compared to the results for the Schelling's parameter, the improvement in the tensile direction is relatively small, which is achieved at a smaller strain, followed by larger decreases in the elastic modulus and the diffusion coefficient.

4.3. Yttria concentration The effect of the yttria concentration has also been considered. Figs. 8 and 9 show the stress–strain curve and the directional diffusion coefficient of 4YSZ and 14YSZ at 2000 K. The apparent elastic modulus is generally higher with a lower yttria concentration because there are fewer oxygen vacancies. On the other hand, the diffusion coefficient of 8YSZ is the highest at 2000 K among these three, which results from the

Fig. 7. Simulation results of 8YSZ at 1273 K calculated using with Brinkman's parameter set: (a) stress–strain curve, (b) directional diffusion coefficient.

Fig. 9. Simulation results of 14YSZ at 2000 K: (a) stress–strain curve, (b) directional diffusion coefficient.

W. Araki, Y. Arai / Solid State Ionics 190 (2011) 75–81

79

Fig. 10. Strain effect on directional mean square displacement of 8YSZ at 1273 K subjected to uniaxial, biaxial, and hydrostatic stress: (a) [100], (b) [010], (c) [001].

lower activation energy and the higher pre-exponential factor of the stabilised zirconia for the higher yttria concentration [11]. The oxygen displacement is obviously enhanced in the tensile direction in both 4YSZ and 14YSZ. The improvement of the diffusion coefficient in the tensile direction is achieved when using the smaller strain for the lower yttria concentration. The maximum 44% improvement is achieved with a strain of 2.0% for 4YSZ while the one of 38% is achieved with a strain of 4.7% for 14YSZ. The improvement is slightly larger for the lower concentration, but again they are not significantly different.

4.4. Biaxial and hydrostatic strains The biaxial and hydrostatic cases have been also considered. Fig. 10 shows the directional mean square displacements of 8YSZ at 1273 K that is subjected to the uniaxial, biaxial, and hydrostatic strains. In the simulation, the same stress of 7.7 GPa, which is equivalent to a strain of 1.8% for the uniaxial case, is simply applied in both the [100] and [010] directions for the biaxial case and all in the [100], [010], and [001] directions for the hydrostatic case. As already explained above, in the uniaxial case, the diffusion in the tensile [100] direction is enhanced whereas the ones in the other directions are almost constant. In the biaxial case, the diffusions in both the [100] and [010] directions are clearly enhanced. In the hydrostatic case, the diffusions in all three directions are enhanced. Accordingly, the total diffusion coefficients are more clearly enhanced for the biaxial and hydrostatic cases, which are 29 and 70%, respectively. Thus, the improvement in oxygen diffusion for the biaxial and hydrostatic cases could be the accumulative effect of the strain in each direction. 5. Discussion

Fig. 11. Radial correlation function of 8YSZ at 1273 K: (a) Zr–Zr, (b) Zr–O, (c) Y–O, (d) O–O.

All the results from the present simulation have demonstrated that (1) the stress–strain curve has a non-linear relationship so that the apparent elastic modulus has a clear decrease when the strain is large, (2) the directional diffusion can be enhanced in the tensile direction and hindered in the compressive direction, and (3) the directional diffusion coefficient in the tensile direction reaches maximum value with a particular strain, where the apparent elastic modulus starts to show a clear decrease. In addition, the maximum improvement is achieved with a smaller strain at a lower temperature and also for the lower yttria concentration. In this chapter, thus, the optimum strain state for the highest improvement in the diffusion coefficient is discussed, allowing for the detailed deformation process in the uniaxial strain case, such as the lattice structure and the oxygen displacement within the lattice. Fig. 11 shows the radial correlation function g(r) of the Zr4+–Zr4+, 4+ Zr –O2−, Y3+–O2–, and O2−–O2− ions of 8YSZ at 1273 K, calculated using different uniaxial strains. The distance of each couple without strain is almost agreed with the experimental data [20]. The change in distance of Zr4+–Zr4+ shown in Fig. 11(a) is almost equivalent to the strain, so the first and second neighbouring distances slightly increase with the strain. The change is quite small because the applied strain is small, where the maximum strain at the fracture is 3.7% for 8YSZ at 1273 K. The changes in the distances of Zr4+–O2− and Y3+–O2− shown in Fig. 11(b) and (c) are also quite small. The first neighbouring distance slightly decreases whereas the second neighbouring distance shows a slight increase as the strain increases. Unlike these three just mentioned, the function O2−–O2− has a much clearer change with the strain. The peak the value of the function g(r) at the second neighbouring distance distinctively decreases and the peak is broader when the strain is larger, which means that the oxygen ions no longer exist on its ideal sites as will

80

W. Araki, Y. Arai / Solid State Ionics 190 (2011) 75–81

Fig. 12. Intensity of the second neighbouring distance of O–O.

be explained in detail later. Fig. 12 shows the peak value of the function g(r) at the second neighbouring distance of O2−–O2− with the strain of 8YSZ at 1273 K. The peak slightly decreases when the strain is smaller than 1.8%, but then it decreases significantly until the fracture, which is similar to the decreasing behaviour of the elastic modulus and, accordingly, could be related to the enhancement of the diffusion coefficient. This similarity was also confirmed in the simulations using different temperatures as shown in Fig. 12 and also for the various parameters discussed in the previous chapter. The deformation process has been closely investigated. Fig. 13 illustrates the observed deformation process of the unit cell. First, without any strain, the lattice structure is cubic-fluorite with the oxygen ions in the ideal site (Fig. 13(a)). When the tensile strain is applied in the [100] direction, the oxygen starts to move parallel to the tensile direction (Fig. 13(b)), followed by another oxygen displacement in the vertical direction (Fig. 13(c)), which are temporally equivalent conformations to the t″- and t′-phases, respectively [21]. As the strain is further increased, the lattice deformation and oxygen displacement simply proceed (Fig. 13(d)). When the oxygen displacement becomes too large for this conformation (Fig. 13(e)), the lattice starts to transform into another structure (Fig. 13(f)), i.e., the orthogonal I phase (Fig. 13(g)) here, but for most cases in the present simulation, a fracture simply occurs without succeeding in the phase transformation, except for 4YSZ and 8YSZ when subjected to large strains. This deformation process also explains the broadening (and decrease) of the function g(r) at the second neighbouring distance of O2−–O2− with a small change in the first neighbouring distance as shown in Figs. 11 and 12.

The deformation process in Fig. 13(a)–(d) corresponds exactly to the slight decrease in intensity of the second neighbouring distance of O2−–O2−, which contributes to the gradual decrease in the elastic modulus due to the slight oxygen displacement. In addition, the process in Fig. 13(e) (before Fig. 13(f)) corresponds to the distinctive decrease and the peak-broadening of the second neighbouring distance of O2−–O2−, which results in the clear decreasing behaviour of the modulus due to the huge oxygen displacement and the lattice deformation. Since the diffusion coefficient in the tensile direction maximises around at the beginning of the decrease of the elastic modulus, the conformation shown in Fig. 13(d) for the largest oxygen displacement is the critical strain state, i.e., the optimum strain state for the maximum diffusion enhancement. For the biaxial and hydrostatic strains, the oxygen ions are similarly displaced in all the strained directions, resulting in the large improvements in total diffusion coefficient. It should be also noted that the oxygen displacement and the improvement in the diffusion coefficient disappear when the applied strain is removed: so, they are reversible. Fig. 14 shows an Arrhenius plot of the directional diffusion coefficient in the [100] direction of 8YSZ when there was no strain and for the optimum strain that produces the maximum level of improvement at each temperature. As the improvements are not significantly different at all temperatures, the activation energy is almost constant at 0.64 eV regardless of the strain. This could be because the diffusion path is generally unchanged, as has been reported in other studies [10,11], but rather, the pre-exponential factor seems to be changed due to the strain. The critical strain state can be achieved at smaller strains when the temperature is lower or when the yttria concentration is also lower, because the cubic-fluorite phase is less stable for those cases, which can help to more clearly explain the present results. In addition, in the results from other experiments [9,12], a higher improvement for a smaller strain has been obtained, which could be attributed to the fact that the crystal structure in the MD simulations are generally much more stabilised in the cubic phase than the actual ones [13] (as is discussed in more detail below). It can be concluded, therefore, that the diffusion can be enhanced by the oxygen displacement caused by the strain and the highest improvement can be obtained when as much of the oxygen as possible is displaced in the fluorite lattice structure. It must be necessary to compare the present simulation and the experimental data on the multi-layered system previously reported [5–8]: however, at present, it is still difficult mainly due to the potential parameters for the zirconia. In the present study, two different potential

Fig. 13. Deformation process of unit cell. The large circle represents the oxygen ion while the smaller circle represents the cations.

W. Araki, Y. Arai / Solid State Ionics 190 (2011) 75–81

81

(probably up to 100% at the highest), which could be only a partial contribution to the improvement in the ionic conductivity of the multi-layered electrolyte system, as it is much lower than the huge improvement in experimentally-obtained data.

6. Conclusion

Fig. 14. Arrhenius plot of the directional diffusion coefficient in the [100] direction.

parameter sets, which both have been frequently used in the literatures, have been used. With the parameter set given by Brinkman et al., the diffusion coefficient and the mechanical properties of YSZ can be simulated relatively well whereas the slightly better simulation of the activation energy can be obtained with the one by Schelling et al. [11]. In addition, as for the crystal structure, the Brinkman's parameter set gives a completely cubic structure regardless of the yttria concentration whereas the Schelling's one gives a slightly tetragonal phase but only with 0 to 2 mol% of yttria concentration (at room temperature) [13], although both parameter sets give a good lattice constant. Thus, there is still a difficulty in obtaining the accurate model which has exactly the same properties and crystal structures as the real one. Despite that, both parameters give a qualitatively same result regarding the strain effect on the oxygen diffusion in the present study. As far as the deformation process explained in Fig. 13, the simulation with a lower yttria concentration or at a lower temperature would be similar to the experimental data (for instance, the simulation of 4YSZ at 1273 K or 8YSZ at 973 K would be similar to the experiment of 8YSZ at 1273 K), because the simulation gives a more cubic-stabilized structure than the real structure. Also, it should be noted that the present simulation is mainly focused on the uniaxial strain effect on the diffusion coefficient, which is essential when considering effects of biaxial or hydrostatic strain. On the other hand, the electrolyte in the multi-layered system is generally subjected to the biaxial strains, on which the present study doesn't include much result. Still, as it might be useful to give a possible quantitative value regarding the improvement in the diffusion coefficient, a rough estimation of the improvement is given here. According to our simulation with the uniaxial strain, the improvement rate in the tensile direction is about 50% (at maximum). With the biaxial strain, which can be considered as the accumulated effect of the uniaxial strains, the in-plane improvement rate could be 60%. (Note that, in the hydrostatic case, the total improvement ratio was about 70%, though it is uncertain if it is the maximum value.) Therefore, the improvement of the diffusion coefficient caused by the elastic strain would be about 60%

The strain effect of the oxygen diffusion in the yttria-stabilised zirconia in the present study was investigated by means of a molecular dynamics simulation. The simulation was conducted for various parameters such as the strain, temperature, yttria concentration, and potential parameter, and in addition, the biaxial and hydrostatic cases were considered. For a uniaxial strain, the oxygen diffusion was enhanced in the tensile direction whereas it was hindered in the compressive direction. For 8YSZ at 1273 K, for instance, the maximum improvement was 44% with a strain of about 1.8%. The maximum improvement was achieved for a smaller strain at a lower temperature and also with a lower yttria concentration. For the biaxial and hydrostatic strains, the total diffusion coefficient was enhanced simply as a result of the enhancement when using the uniaxial strain in each direction. The possible improvement of the in-plane diffusion coefficient caused by the bi-axial elastic strain could be about 60%, which might partially contribute to the huge improvement reported in the studies on the multi-layered system. The detailed deformation analysis reveals that the optimum strain state for the highest oxygen diffusion in the tensile direction can be obtained when the oxygen ion is largely displaced in the fluorite lattice structure.

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21]

S.P.S. Badwal, F.T. Ciacchi, Ionics 6 (2000) 1. V.V. Kharton, F.M.B. Marques, A. Atkinson, Solid State Ionics 174 (2004) 135. J.W. Fergus, J. Power Sources 162 (2006) 30. J.A. Kilner, Faraday Discuss. 134 (2007) 9. N. Sata, K. Eberman, K. Eberl, J. Maier, Nature 408 (2000) 946. A. Peters, C. Korte, D. Hesse, N. Zakharov, J. Janek, Solid State Ionics 178 (2007) 67. N. Schichtel, C. Korte, D. Hesse, J. Janek, Phys. Chem. Chem. Phys. 11 (2009) 3043. S. Sanna, V. Esposito, A. Tebano, S. Licoccia, E. Traversa, G. Balestrino, Small 6 (2010) 1863. W. Araki, Y. Imai, T. Adachi, J. Eur. Ceram. Soc. 2 (2009) 2275. W. Araki, Y. Arai, Solid State Ionics 181 (2010) 441. W. Araki, Y. Arai, Solid State Ionics 181 (2010) 1534. T. Izumi, R. Narumi, K. Sato, K. Suzuki, K. Yashiro, T. Hashida, J. Mizusaki, ECS Trans. 838 (2008) 7. P.K. Schelling, S.R. Phillpot, D. Wolt, J. Am. Ceram. Soc. 84 (2001) 1609. M. Fevrè, A. Finel, R. Caudron, Phys. Rev. B 72 (2005) 104117. M.J.L. Sangster, M. Dixon, Adv. Phys. 25 (1976) 247. A. Kushima, B. Yildiz, J. Mater. Chem. 20 (2010) 4809. H.W. Brinkman, W.J. Briels, H. Verweij, Chem. Phys. Lett. 247 (1995) 386. P.S. Manning, J.D. Sirman, R.A. De Souza, J.A. Kilner, Solid State Ionics 100 (1997) 1. M. Kilo, C. Argirusis, G. Borchardt, R.A. Jackson, Phys. Chem. Chem. Phys. 5 (2003) 2219. P. Li, I.W. Chen, J.E. Penner-Hahn, Phys. Rev. B 48 (14) (1993) 10074. M. Yashima, K. Ohtake, M. Kakihana, H. Arashi, M. Yoshimura, J. Phys. Chem. Solids 57 (1996) 17.