Atomistic Simulation of Pure and Doped BaSnO - American

Oct 29, 2009 - Laboratoire Structures, Propriétés et Modélisation des Solides, Grande Voie des Vignes, Ecole Centrale Paris,. 92295 Châtenay-Malabry ...
1 downloads 0 Views 211KB Size
20486

J. Phys. Chem. C 2009, 113, 20486–20492

Atomistic Simulation of Pure and Doped BaSnO3 Y. Z. Wang,† E. Bevillon,† A. Chesnaud,‡ G. Geneste,† and G. Dezanneau*,† Laboratoire Structures, Proprie´te´s et Mode´lisation des Solides, Grande Voie des Vignes, Ecole Centrale Paris, 92295 Chaˆtenay-Malabry Cedex, France, and Centre des Mate´riaux, Mines Paris, Paristech, CNRS UMR 7633, BP 87, F-91003 EVry cedex, France ReceiVed: June 11, 2009; ReVised Manuscript ReceiVed: September 15, 2009

Atomistic simulations have been used to investigate energetic aspects of defects in pure and doped BaSnO3, in close comparison with experimental data. Semiempirical potentials were used to determine the solution energies of trivalent cation on either A or B site of the perovskite structure. Density functional theory calculations were also performed to gain further insight into the redox properties of this material. Both approaches give a complete view of the way doping atoms will be inserted inside the BaSnO3 structure. For instance, we calculated with good precision the evolution of cell parameters with dopant nature and concentration. Then, by taking into account dopant vacancy interactions, we also show that the oxygen vacancy preferentially locates close to the dopant in the case of small size dopant while stays at the second neighbor position in the case of bigger dopant. These calculations also show that In and Y are the dopants giving the most homogeneous energetic map for oxygen vacancy movements and might thus be anticipated as the best acceptor-like dopants for enhanced anion conduction properties. Finally, we calculated an oxygen ion migration activation energy of 0.65 eV. 1. Introduction Acceptor-doped BaSnO3 was presented few years ago as a good proton conductor1-3 with potential applications in fuel cells. Sn here was partially substituted by trivalent metals Y1and In2 and oxygen vacancies appearing by charge compensation were readily filled by water in wet atmospheres. More recent applications concern gas sensors4,5 or photoelectrochemical applications.6,7 At high hydrostatic pressure, an interesting dielectric behavior has also been predicted.8 Despite the potential applications of BaSnO3, its use stays limited for now and to our knowledge, there are only very few experimental and theoretical studies on BaSnO3 based compounds, especially for ion conduction properties. This is all the more surprising since BaSnO3, due to its similarities with BaZrO3 and BaCeO3, may present some fundamental interest for a better understanding of anion and proton perovskite conductors. In our research group, an extensive study of acceptor-doped barium stannate compounds was thus initiated from a theoretical and experimental approach. From one DFT study, we have already shown that the dopant-oxygen vacancy and dopant-proton interaction depends strongly on the dopant nature. In the case of big dopants, the next nearest neighbor oxygen position relative to the dopant position has been found to be more stable than the nearest neighbor position.9 An older study based on an atomistic study also gave some first results on the surface energies and redox reactions.10 Atomic simulations based on semiempirical interatomic potentials have been applied successfully for the study of ion conducting compounds. In particular, this helped to calculate defect energies and to propose ionic transport model at the atomic level,11-14 having a predictive role in the design and improvement of materials. * To whom correspondence should be addressed. Tel.: +33 (0) 141131420. Fax: +33 (0) 141131437. E-mail: [email protected]. † Ecole Centrale Paris. ‡ CNRS UMR 7633.

In this paper, we have explored by atomistic simulations the properties of pure and doped BaSnO3. These simulations were essentially based on semiempirical potentials but DFT calculations were also performed to give complementary information in particular on electronic properties. We will show that these calculations help not only to give local information on the defects energies but also on their macroscopic counterpart such as the evolution of cell parameters with doping. This article is voluntarily limited to dry samples, that is, the incorporation of protons inside the structure will not be studied here. 2. Methods 2.1. Simulation Techniques. Two kinds of approaches were used here, one based on semiempirical potentials, the other one based on DFT calculations. The first one is less CPU-costing but cannot give precise information on electronic states. The second one is in this sense more complete but limited to smaller systems. Atomistic semiempirical simulations start with the specification of an interatomic potential model, which describes the potential energy of the system as a function of the atomic coordinates. This allows modeling both perfect and defective lattices. In both cases, an energy minimization procedure allows finding the lowest energy configuration for the crystal structure. One of the most widely used models is the Born model that consists of the sum of a long-range Coulombic interaction between charged ions and of short-range pair potentials of the Buckingham form

( )

Vij(rij) ) A exp -

rij C - 6 F rij

(1)

where A, F, and C are adjustable parameters assigned to cation-anion or anion-anion interactions. In the present study, the values of these parameters were taken from a previous

10.1021/jp9054878 CCC: $40.75  2009 American Chemical Society Published on Web 10/29/2009

Atomistic Simulation of Pure and Doped BaSnO3

J. Phys. Chem. C, Vol. 113, No. 47, 2009 20487

TABLE 1: Interatomic Potentials Parametersa M· · ·O O2Ba2+ Sn4+ La3+ Nd3+ Sm3+ Gd3+ Y3+ In3+ Sc3+ Fe3+

A (eV)

F (Å)

C (eV, Å6)

Y (e)

9547.96 873.829 1056.8 2088.89 1995.20 1944.44 1885.75 1766.40 1495.65 1575.85 1414.60

0.2192 0.3863 0.3683 0.3460 0.3430 0.3414 0.3399 0.33849 0.3327 0.3211 0.3128

32.00 0.0 0.0 23.25 22.59 21.49 20.34 19.43 4.33 0.0 0.0

-2.389 9.203 1.58 3.0 3.0 3.0 3.0 3.0 3.0 3.0 3.0

K (eV Å-2) ref 23.09 459.2 2037.8

15 15 15 16 16 16 16 16 17 17 17

a A, F, C for Buckingham potential and Y and K for shell charge and spring force constant.

TABLE 2: Calculated and Experimental Properties of BaSnO3 properties lattice parameter (Å) lattice energy (eV) static dielectric constant high frequency dielectric constant elastic constant (GPa) Young modulus (GPa) Poisson coefficient Shear modulus (GPa) compressibility (GPa-1)

semiempirical this study 4.1236 -144.32 23.63 4.24 C11)280.36 C44)137.64 189.72 0.329 100.4 5.4 × 10-3

DFT28

experiment29

4.156

4.116

22.3 4.89

207 0.233 84

244 99.9 6.86 × 10-3

study15-17 (Table 1). Those corresponding to Ba and Sn were shown to well reproduce the perfect lattice properties of simple oxides and also of BaSnO3. These properties are gathered in Table 2. Besides being able to simulate perfect lattices, an important feature of this simulation is to predict how the lattice accommodates various kinds of defects. For this purpose, the enthalpy of a defect, such as a vacancy, an interstitial ion, or a substitutional ion can be calculated via the Mott-Littleton approach.18 It consists in considering two regions around the studied defect. In the region closer to the defect, the ions are strongly perturbed by the defect and are thus allowed to relax till the forces applying on atoms tend to be zero. In the second region, the atoms are considered to be less perturbed. In the region 2a, the forces on the individual ions due to short-range interatomic potentials and the Coulomb term are explicitly calculated and the local displacements evaluated. In the region 2b, the relaxation energy is evaluated taking into account only the colombic (long-range) potential and is evaluated implicitly since this region goes to infinity. In the present calculations, the radii of regions 1 and 2a were 12 and 24 Å respectively. Region sizes were chosen so as to represent a compromise between calculations time and the evolution of the defect energy with the chosen radii values. In all calculations, the electronic polarizability of O2-, Ba2+, and Sn4+ ions are taken into account through the shell model. This model consists of a mass-less shell with charge Y|e| that is allowed to move with respect to a massive core of charge X|e|. The core and shell charges are connected by a harmonic spring force constant k. The use of the shell model, if not necessary for structure optimization, is mandatory if we consider other properties like the optical dielectric constant. It may also be useful when considering interstitial defects. In all case, calculations have been achieved using the Gulp code.19

To give complementary information to semiempirical simulations, we also performed a number of density functional calculations with the SIESTA code.20,21 We used the generalized gradient approximation in the treatment of Perdew, Burke, and Ernzerhof (GGA-PBE22 functional) to describe the exchangecorrelation energy. The core electrons are treated through Troullier-Martins pseudopotentials,23 that include as valence electrons the 5s, 5p, and 6s for Ba, the 4d, 5s, and 5p for Sn, and the 2s and 2p for O. This numerical scheme provides for cubic BaSnO3 an equilibrium lattice constant of 4.21 Å with a typical overestimation of the GGA (while a previous calculation performed with the ABINIT code and no semicore electron for Sn provided 4.156 Å). The defects studied by this method (neutral oxygen vacancy, neutral oxygen interstitial) are included in a 2 × 2 × 2 supercell (40 atoms if no defect). The configurations are optimized until all the components of the Hellman-Feynman forces are below 0.02 eV/Å. We mainly used density-functional calculations to study the energetics of some elementary redox reactions in barium stannate, and also to examine the charge distribution in specific cases. The generalized gradient approximation has been used rather than the local density approximation (LDA) due to its well-known ability to reproduce quite well the energetics of microscopic processes. It is known that the LDA largely overestimates all the cohesive energies, bond energies, and so forth while the GGA (at least the most recent ones like the PBE) provides in general with accurate energies. The reader is referred for instance to the very complete work of Bjo¨rketun et al.24,25 on barium zirconate, performed in the GGA and providing hydration energies in good agreement with experiments. In some cases anyway, we performed calculations using charged supercells, by using as usual a uniform compensating background to neutralize the system. Besides, adding one electron in the supercell would probably lead to electronic localization (polaronic state) on the Sn atoms (able to reduce from +4 to +2 in BaSnO3) but our GGA calculations systematically provide a uniform distribution of the additional electronic charge. Such a result can be questionable and related to the possible inability of the standard functionals to describe such localized states. The use of other approximations such as GGA+U might improve the description of such polaronic states. 2.3. Defect volume and lattice parameters. Also, the Mott-Littleton defect calculations have been used to predict unit cell volumes using an approximation of no interaction between defects. More concretely, the defect relaxation volume V can be calculated from17,26,27

V ) -kTVC

( ) ∂fv ∂V

T

(2)

where kT is the isothermal compressibility and VC is the unit cell volume. (∂fv/∂V)T is the change at constant temperature in the Helmholtz defect formation energy as a function of the change in the unit cell volume. Given the defect relaxation volume, the unit cell volume can then be determined from the following equation

Vcell ) where

∑ Vdefect × cdefect + Vperfect

(3)

20488

J. Phys. Chem. C, Vol. 113, No. 47, 2009

∑ Vdefect ) VM

 Sn

1 [M] + VV..o, cdefect ) ) ,V 2 [M] + [Sn] perfect (4.1236)3 ) 70.118 Å3

This approach will be used to determine both the influence of dopant size and dopant content on cell parameters. We will see that this approach works very well in almost all cases even if the approximation of noninteracting defects is obviously not correct for high doping levels. 2.4. Experimental Section. In order to give more weight to our calculations, we also prepared some doped BaSnO3 compounds and compared their properties to those obtained from calculations. The pure and doped BaSnO3 powders were prepared by a gel polymerization route.28 The precursor obtained was heat-treated at 1473 K during 4 h leading to the desired compounds. The powders were then ball-milled with yttria stabilized zirconia balls at 400 rpm for 2 h. Green pellets (10 mm in diameter, 1-2 mm thickness) were obtained by pressing ceramic powders uniaxially at 4 tons for 1 min, and then applying an isostatic pressure of 750 MPa during 10 min. At this stage, the typical density for the green pellets was about 60% of the theoretical one. The samples were then sintered at 1873 K for 12 h. The X-ray powder diffraction patterns of pellets were recorded on a high accuracy two-axis diffractometer in a Bragg-Brentano geometry using CuKR1,2 radiation as X-ray source which was generated by a 18 kW rotating anode. The data were collected, in the 2θ region between 20 and 110° with a 0.02° step size. The lattice parameters of the samples were determined by a combined use of the Topas and cellref softwares. 3. Results and Discussion 3.1. Basic Properties of BaSnO3. Before studying defects, it is necessary to simulate the properties of pure BaSnO3. Table 2 shows the values of simulated properties including lattice parameter, lattice energy, mechanical properties, and dielectric properties. The calculated a lattice parameter of the Pm3m perovskite structure is found to be in very good agreement with the experimental value. Indeed, the lattice parameter of BaSnO3 is reproduced to within 0.17% of the experimental value. The calculated mechanical and dielectrics properties are also consistent with the experimental values29 and with those obtained by means of first-principles density functional calculations.28 To evaluate thermal expansion coefficients, the zero static internal stress approximation method (ZSISA) was used in this work to determine the cell parameter at different temperatures. In the ZSISA method, only the strain derivatives with respect to the free energy are used while the internal derivatives neglect the free energy contribution.30 This approach has been found to be more stable at high temperature than total free energy minimization. The temperature evolution of the lattice parameter (not shown) is linear and the thermal expansion coefficient calculated is 1.06 × 10-5 K-1 (300-800 K), close to the experimental data 9.6 × 10-6 K-128,29 for the same range of temperatures. Therefore, the potentials used allow to reasonably well reproducing the structural, elastic, and dielectric properties of pure BaSnO3 and are thus validated for the following calculations. 3.2. Intrinsic Atomic Defects. We first calculated the main point defects formation energy. Solutions energies were then evaluated considering isolated, that is, noninteracting point defects (vacancies and interstitials). To calculate the energy of a vacancy, an ion was extracted from its lattice site and sent to infinity. An interstitial was treated similarly by taking an ion

Wang et al. TABLE 3: Isolated Point Defects

defect nature

defect equilibrium

2+

Ba vacancy Sn4+ vacancy O2- vacancy O2- interstitial Schottky

′′ ′′ BaBa× + SnSn× + 3OO× ) VBa ′′ + VSn + 3VO· · + BaSnO3 O Frenkel OO× ) VO· · + Oi′′ Ba Schottky-type BaBa× + OO× ) VBa ′′ + VO· · + BaO ′′ ′′ + 2VO· · + SnO2 Sn Schottky-type SnSn× + 2OO× ) VSn

formation energy (per defect /eV) 19.1 83.99 18.64 -9.05 2.94 4.795 2.615 6.757

from infinity and placing it at an interstitial position. The calculated energies are reported in Table 3 with the exception of cation interstitial and Frenkel energies, since they led to calculations instabilities. Anyway, they are undoubtly higher in energy and thus will not be considered further. The lattice energy used in the calculation of Schottky-type defects were obtained from atomistic simulation of BaSnO3, SnO2, and BaO. Examination of Table 3 reveals that the creation of Frenkel and Schottky defects are associated to significant energies. The most favorable energy corresponds to the formation of Ba Schottky-type disorder, resulting possibly in BaO loss at high temperature. This behavior is very well-known in Ba-containing oxides and a recent study on cation nonstoichiometry in BaCeO3 confirms that barium loss does occur at high temperatures.31 However, as a general rule, the magnitude of the calculated energies suggests a low concentration of intrinsic atomic defects. 3.3. Redox Reactions. It is essential to understand the stability and behavior of electronic defects with respect to the lattice defects since most applications, that is, fuel cells or gas sensors, require the materials to be used in highly reducing and/ or highly oxidizing atmospheres. Allan et al.10 investigated defect and surface structure of the ASnO3 (Ca, Sr, Ba) stannates as oxidation catalysts by atomistic simulations. The formation energies of cation vacancies and interstitials were calculated to be very high compared to that of anions. Therefore, in the present study we only considered the possibility of oxygen vacancy for reduction and oxygen interstitial or oxygen vacant site filling for oxidation. For the undoped system, the following reactions for oxidation and reduction were considered. For reduction

1 OOx ) VO· · + O2(g) + 2e 2

(4)

1 O (g) ) Oi′′ + 2h · 2 2

(5)

For oxidation:

Taking into account the fact that charge-compensating oxygen vacancies are created through acceptor-doping, we also considered the following equation for acceptor-doped samples

1 VO· · + O2(g) ) OOx + 2h · 2

(6)

These calculations were performed by means of DFT. For the reduction reaction, a neutral oxygen vacancy is placed in the supercell described above and the system is optimized. The

Atomistic Simulation of Pure and Doped BaSnO3

Figure 1. Structure of BaSnO3 with one oxygen interstitial atom (barium, tin, and oxygen atoms are represented respectively with blue, gray, and red spheres). The peroxide species is visible on the left part of the drawing.

question is whether the two electrons left by the missing oxygen atom occupy an extended state (like those of the conduction band in a semiconductor) or a spatially localized state (polaronlike). After optimization, and although the oxygen vacancy would be a priori surrounded by two equivalent Sn atoms, we observed that only one Sn atom has a Mulliken charge that differs significantly from the others; seven among the eight Sn atoms of the supercells have equal Mulliken charge within 0.03 e, while the eighth (one of the two surrounding the vacancy) has a Mulliken charge 0.18 e lower than the seven others. Ba atoms have all their charges equal within 0.01 e, while O atoms also have charges equal within 0.04 e. This result indicates that the charge left by the missing O is essentially localized on one single Sn (small bound polaron) and does not occupy any extended state. It is noteworthy that the charge does not equally localize on the two Sn surrounding the vacancy but prefers to localize on one single Sn. This is compatible with the known oxidation degrees of Sn (+4 and +2). As a consequence, the distance from the reduced Sn to the neighboring O is increased (between 2.16 and 2.20 Å) with respect to the Sn-O distance in pure barium stannate (2.105 Å in GGA). This would justify the possible semiempirical modeling (not considered here) of reduced barium stannate by the replacement of one Sn4+ ion by one with lower charge. It is noteworthy that adding (respectively removing) an additional electron in the supercell (without introducing any point defect) results in a uniform distribution of the additional (respectively missing) charge on the Sn atoms rather than the localization on a single Sn. It means that in this latter case, the additional (respectively missing) charge corresponds to the density of probability of an extended conduction band state (respectively valence band state). The reduction of Sn from +4 to +2 is thus either intrinsically related to the existence of the point defect (O vacancy), or would require corrections to the exchange-correlation effects beyond GGA (such as GGA+U) to allow such an electronic localization without point defect in the structure. For the oxidation reaction, we first examined pure BaSnO3 in which one extra oxygen atom is inserted as interstitial in the 2 × 2 × 2 supercell. It is clear that neither Ba nor Sn can be more oxidized. In such systems, one usually observes the formation of suboxide species (peroxide O22- or superoxide ions O2-).32,33 This is precisely what happens in this case; the extra O atom forms a bond with one oxygen atom of the structure (see Figure 1). The two atoms share the charge (equal Mulliken charges within 0.01, which are higher than the charges of all

J. Phys. Chem. C, Vol. 113, No. 47, 2009 20489 the other O atoms by approximately 0.21 per oxygen atom), and their bond is 1.49 Å, which is typically the bond length of peroxide ions O22-. The charges of the Ba and Sn are almost not affected. In acceptor-doped compounds, no oxygen interstitial is necessary since oxygen vacancies are already present. We have thus studied the case of In-doped and Y-doped barium stannate (25% doping rate) without oxygen vacancies, that is, under strong oxidizing atmosphere. It is usually admitted that filling the oxygen vacancies leads to the appearance of holes in the electronic structure. Indeed, we find the systems to be metallic with empty states close to the valence band top (holes). In this case, the hole is most often taken into account by considering a polaron on oxygen atoms. We also found experimentally that acceptor-doped barium stannate compounds present a p-type behavior at high temperature for high oxygen partial pressure. Concerning energies associated to these reactions, the reduction equation (4) is associated to an energy value of 4.0 eV. This value is high but significantly smaller than the one of 6.6 eV obtained for BaZrO3 from similar DFT calculations.25 This would signify that BaSnO3 would reduce more easily than BaZrO3. Concerning oxidation reaction, the insertion of a neutral oxygen atom is associated with a very low positive energy of 1.1 eV. This means that this reaction is not favorable. Besides, this reaction becomes less and less favored at high temperature due to the entropic contribution of gaseous oxygen to the reaction total free energy. The oxidation reaction (6) is associated with a negative energy of respectively -0.95 eV in agreement with the calculated value of +1.21 eV for the inverse reaction of charged (+2) oxygen vacancy formation energy in BaZrO3 by Sundell et al.25 This value does not account properly for the stabilizing interaction between dopants and charged oxygen vacancies in the initial state of the reaction. Thus, a more precise treatment would lead to a much less negative value since the interaction energies between the Y (respectively In) dopants and charged O vacancies are stabilizing for the left member of the reaction and have been calculated in ref 9 to be -0.75 (respectively -0.63) eV. A negative value of enthalpy would mean that the system tends toward the incorporation of oxygen and holes in oxidative atmosphere at least at low temperature and that charged oxygen vacancies would be stable only at high temperature. In reality, if the literature on this subject is scarce, no experimental results come to confirm the incorporation of oxygen at low temperature. Besides, if a p-type behavior has indeed been measured in barium stannate in air at temperatures between 600 and 900 °C, the conduction values follow a 1/4 power law with oxygen partial pressure meaning that the main charge carrier in terms of concentration are still oxygen vacancies. 3.4. Dopant-Ion Substitution. The substitution of a range of trivalent cations was investigated on both Ba and Sn sites in order to gain some insight into substitution energy cost and subsequent preferential sites of dopants in the host lattice. The substitution of dopants on the Sn site was analyzed according to the standard “acceptor-doping” defect reaction, with the creation of oxygen vacancies as charge-compensating defects

1 1 1 ×  + O×o f MSn + V..o + SnO2(s) M O (s) + SnSn 2 2 3 2 2

(7) For substitution on the Ba site, dopant incorporation was described by

20490

J. Phys. Chem. C, Vol. 113, No. 47, 2009

1 3 × 1  3 · f MBa + VBa + BaO(s) M O (s) + BaBa 2 2 3 2 2 2

Wang et al.

(8)

which leads to BaO elimination. The energies associated to these reactions were evaluated by combining corresponding defect and lattice energy terms. Such an approach provides a useful systematic guide to the relative energies for the doping of different atoms at the same site. The resulting solution energies for substitution on the Ba and Sn sites as a function of ionic Shannon radius34 are presented in Figure 2. Six-fold radius is used here in all cases for reasons of simplicity. Examination of the results reveals an important correlation between solution energy and ion size. This was already the case for other perovskite compounds.12,13,35,36 The solution energy seems to lower continuously with dopant size on A-site, while it passes through a minimum value for dopant on B-site. Indeed, from the values of solution energy we predict In3+ and Y3+ to be among the most favorable dopants on B-site for BaSnO3. This agrees with the experimental work of Schober1,2 where these cations were used as acceptor dopant for effective proton conduction. Nevertheless, the solution energy on B-site is almost the same for Sc3+, In3+, and Y3+ cations even if their ionic radius differs by more than 0.15 Å. Besides, Y solution energy on B-site is one of the lowest in spite of the big radius difference between Y3+ (0.9 Å34) and Sn4+ (0.69 Å34). This situation is rather different from the one observed in fluorite compounds for which the lowest solution energy generally corresponds to a dopant atom with a radius very close to that of the host atom.17,27,37 Nevertheless, such behavior has already been observed in perovskite compounds. Indeed, Islam et al.35 showed that in barium zirconate the favorable substitution on B-site corresponds to Y3+ and Gd3+ whose radii differ also by more than 0.148 Å with that of Zr4+. This can be interpreted by the fact that perovskite structure is known to accommodate much deviation from the ideal structure. Besides, in the case of BaSnO3 the tolerance factor is calculated to be 1.02, which means that the radius of the tin ion on B-site is small if compared to the available volume. This might explain why bigger atoms can be incorporated in the structure without an important energy cost. It is interesting to note that Sm3+ cation sits at the crossover point in Figure 3 which suggests possible substitution at either Ba or Sn sites. This also exists in other perovskite oxide, such as Nd in BaZrO3 and Y in CaZrO3.35 Our simulations suggest that similar effect might occur for larger dopants (e.g., Nd and La), which have a lower calculated solution energy on Ba site. In fact, the real situation is somewhat different since cations are incorporated in stoichiometric proportion corresponding to substitution on either site. Nevertheless, this means that for larger dopant some partial substitution on A-site should not be discarded. A similar phenomenon was reported in barium cerate compounds.13 Indeed, Haile et al.38 showed that B-site doped BaxCe0.85M0.15O3-d (x ) 0.85-1.20) can exhibit significant Nd dopant partitioning over both Ba and Ce sites. 3.5. Defect Volume and Lattice Parameters. In this section, the variation of lattice parameter with dopant radius and solute concentration was calculated. In this case, we only considered the substitution of atoms on B-site. We then compared these calculations to experimental data. Table 4 shows the volume associated to the punctual defect according to the method described above. We observe first that oxygen vacancy has a negligible influence on cell volume evolution. Actually, from these results, 33% of oxygen vacancies in BaSnO3 would result in a less than 1% contraction on the cell parameter. This was already shown by Mogensen et al.39 who demonstrated that

Figure 2. Calculated solution energy as a function of dopant radius.

Figure 3. Comparison of experimental and predicted lattice parameter. (a) Lattice parameter as a function of dopant radius in BaSn0.75M0.25O2.875 (b) Lattice parameter as a function of dopant concentration in BaSn1-xYxO3-x/2 (acceptor doping is compensated by oxygen vacancy creation).

oxygen vacancy (RV..o ) 1.3771∼1.4045 Å) has a volume similar to that of an oxygen ion. Then, we observe that only Fe3+ substitution on B-site cation leads to a negative defect volume, which corresponds to the fact that this atom has a smaller radius than tin. Figure 3 shows the predicted lattice parameters as a function of dopant radius and concentration for BaSn0.75M0.25O3-δ (M ) Fe, Sc, In, Y, Gd, Sm, Nd, and La) and BaSn1-xYxO3-δ (x ) 0-0.375). The examination of results reveals that lattice parameter predictions reproduce well the experimental data for the BaSn1-xMxO3-δ system, except for In3+, Sm3+, Nd3+, and La3+. For In3+ substitution, it is probable that the actual experimental composition is different from the nominal formula due to indium evaporation during thermal treatment. This would explain why the cell parameter of indium substituted compound is smaller than expected. For Sm3+, Nd3+, and La3+, the lattice

Atomistic Simulation of Pure and Doped BaSnO3

J. Phys. Chem. C, Vol. 113, No. 47, 2009 20491

TABLE 4: Defect Relaxation Volume As Calculated from Semi-Empirical Calculations defect

defect relaxation volume (Å3)/defect

Vo· · FeSn ′ ScSn ′ InSn ′ YSn ′ GdSn ′ SmSn ′ ′ NdSn LaSn ′

-1.35 -3.22 3.54 7.34 13.97 17.36 19.31 21.19 24.9

parameters predicted are higher than the experimental ones possibly due to partial substitution on A-site. This would confirm the tendency foreseen by solution energy calculations presented above. Figure 3 also shows the evolution of cell parameter for Y doped compounds with doping level ranging from 0 to 37.5%. In this case, the agreement is very good in spite of the fact that for these calculations defects are considered as noninteracting, which is evidently far from the reality for high doping levels. 3.6. Defect Association. It is well-known that the chargecompensating oxygen vacancies, which form upon acceptor doping, tend to interact attractively with dopant, leading to defect clustering. The formation of such clusters adds an additional binding (or association) term to the activation energy for oxideion conductivity. The binding energy is defined as the difference between the formation enthalpy of the appropriate defect cluster and the sum of the formation enthalpies of the isolated point defects

Ebind ) Ecluster -

∑ Eisolated-defect

(9)

We have considered the energies of simple defects config′ /Vo..) · at an adjacent uration consisting of a pair cluster (MSn (nearest neighbor, NN) site and at a second coordination shell (next-nearest neighbor, NNN) site. The cluster binding energies associated to these two configurations are reported in Figure 4. Figure 4 reveals that the most stable configuration of the defect cluster depends on the size of the dopant cation. For small dopant cations (such as Fe3+, Sc3+ and In3+), it is energetically more favorable for the vacancy to locate at a NN site, whereas larger dopant cations (such as Gd3+, Sm3+, Nd3+ and La3+) favor the oxygen vacancy on the NNN site. For Y3+, the difference of binding energy is small between NN and NNN configuration. Thus as the oxygen vacancy moves away from the Y3+ ion,

Figure 4. Binding energy of doping atom and oxygen vacancies as a function of dopant radius.

Figure 5. Contour plot of potential energy surface for oxygen vacancy migration (showing the weakly curved path between adjacent anion sites).

there is very little change in the binding energy. This indicates that Y doping would probably lead to the best conduction properties for this system. This inversion between NN and NNN site stability for the O vacancy has also been observed from ab initio calculations9 at large dopant concentration. Although the larger concentration leads to larger binding energies, the largest dopants are also found in most cases by these authors to stabilize the O vacancy in the NNN site. 3.7. Oxygen Ion Migration. Finally, the simulations can give some information about oxygen migration. It is generally accepted that oxygen diffusion is based on a conventional hopping mechanism.40 Atomistic simulations have been successfully used to investigate oxygen vacancy migration in compounds like AZrO3,12 ACeO3,13,41 and LaBO3.40,42 To find the activation energy for oxygen ion migration, it is necessary to locate the lowest saddle point on the energy surface for the various possible pathways. The saddle point is the highest point along the lowest energy pathway for an oxygen ion moving between two unoccupied sites. The activation energy for O2migration is then just the difference in lattice energy between the two configurations consisting of the moving O2- positioned on the saddle point and of the moving oxygen ion positioned on a crystallographic oxygen site. In this study, the energy profile was mapped out by calculating the defect energy of an oxygen ion migrating along the diffusion path between adjacent oxygen sites followed the method of Cherry et al.40 Each energy point was obtained after relaxation of the lattice after fixing the coordinates of the migrating ion. The contour plot of potential energy surface for oxygen vacancy migration is presented in Figure 5. It clearly shows that the migration “channel” follows a small deviation from the direct path between two oxygen sites. Migrating ion within the perovskite lattice is commonly assumed to take a direct linear path along the edge of the BO6 octahedron for which the saddle point is located midway between the two anion sites.40 However, Figure 5 shows that oxygen follows a curved route with the saddle point bowed away from the adjacent Sn cation. This curved pathway had also been foreseen in other perovskite compounds like BaCeO313 and CaZrO3.11 In the case of BaSnO3, the saddle point obtained is displaced at a distance of 0.12 Å from the linear path and the migration energy is calculated to be 0.65 eV. It has been reported that the migration path in BaCeO3 with cubic structure is observed at a distance of 0.35-0.40 Å from the linear path and the migration energy is 0.84 eV.13 For undoped CaZrO3 and SrZrO3, the same distances are found to be 0.42 and 0.20 Å and the activation energies are 0.42 and 0.58 eV, respectively.12 These results tend

20492

J. Phys. Chem. C, Vol. 113, No. 47, 2009

to show that the curvature is more pronounced for larger B-site atoms. The small value of such distance for BaSnO3 is also coherent with the fact previously mentioned that tin is small on its site in this compound. According to experimental results, oxygen migration energies lies in the range from 0.7 to 1.2 eV in acceptor-doped BaSnO3 compounds. According to DFT calculations, the oxygen migration energy lies in the range from 0.65 and 0.79 eV also in acceptor-doped compounds. Thus the value obtained here might reflect only part of the vacancy hopping mechanism. In particular, it is important to keep in mind that this calculated energy relate to intrinsic migration of oxygen vacancies and does not include possible association energies or dopant effects. These additional energy terms may account for the larger activation energies that are found experimentally. 4. Conclusions In this paper atomic simulation techniques have been used to get further insight into the defect properties of pure and doped BaSnO3. The conclusions are drawn from our results as follows. (1) The potentials used reproduce well the observed properties of cubic BaSnO3 and might be used for other stannate compounds. Actually, pyrochlore structures (not presented here) are also very well reproduced with the tin-oxygen potential used here. (2) Calculated energies for the creation of intrinsic defects are relatively high with the most favorable energy for the formation of BaO vacancy pairs. (3) Oxidation and reduction behavior is similar to that observed in other barium-based perovskite compounds though a higher reducibility might be observed if compared to barium zirconate. (3) The trivalent dopants, Sc, In, and Y have relatively favorable solution energies at the Sn site. For cations greater than Sm3+, some partial substitution on A-site cannot be discarded and even is in part proved by experimental values. (4) The predicted lattice parameters, deduced from defect volume calculations, as a function of dopant radius and concentration are in good agreement with experimental values. ′ /Vo..) · association with dopant (5) The binding energy of (MSn size shows that small size dopants prefer a cluster geometry in which the oxygen vacancy resides in a nearest neighbor site, while for the big size dopants, the oxygen vacancy is located in a next nearest neighbor site. The yttrium dopant is the one that probably leads to the smallest association effects. (6) The oxygen migration energy is calculated to be around 0.65 eV, in partial agreement with the experimental values. Oxygen ions were found to follow a slightly curved trajectory on hopping between adjacent vacancy sites. Acknowledgment. Part of the study has been supported by the Agence Nationale de la Recherche through the NANOMATSOFC “jeunes chercheurs” project with number JC05_42231. We would like to thank Professor Julian D. Gale for letting us using the potentials he refined on barium stannate compounds. References and Notes (1) Murugaraj, P.; Kreuer, K. D.; He, T.; Schober, T.; Maier, J. Solid State Ionics 1997, 98, 1.

Wang et al. (2) Schober, T. Solid State Ionics 1998, 109, 1. (3) Kreuer, K. D. Solid State Ionics 1999, 125, 285. (4) Cerda, J.; Arbiol, J.; Dezanneau, G.; Diaz, R.; Morante, J. R. Sens. Actuators, B 2002, 84, 21. (5) Manorama, S. V.; Reddy, C. V. G.; Rao, V. J. Appl. Surf. Sci. 2001, 174, 93. (6) Borse, P. H; Lee, J. S; G.Kim, H J. Appl. Phys. 2006, 100, 124915. (7) Yuan, Y.; Lv, J.; Jiang, X.; Li, Z.; Yu, T.; Zou, Z.; Ye, J. Appl. Phys. Lett. 2007, 91, 094107. (8) Bevillon, E.; Geneste, G. Phys. ReV. B 2007, 75, 214106. (9) Bevillon, E.; Gregory, G. Phys. ReV. B 2008, 77, 184113. (10) (a) Hines, R. I.; Allan, N. L.; Flavell, W. R J. Chem. Soc., Faraday Trans. 1996, 92 (12), 2057–2063. (b) Hines, R. I.; Allan, N. L.; Flavell, W. R J. Chem. Soc., Faraday Trans. 1996, 92, 7. (11) Islam, M. S. J. Mater. Chem. 2000, 10, 1027. (12) Davies, R. A.; Islam, M. S.; Gale, J. D. Solid State Ionics 1999, 126, 323. (13) Glockner, R.; Islam, M. S.; Norby, T. Solid State Ionics 1999, 122, 145. (14) Khan, M. S.; Islam, M. S.; Bates, D. R. J. Phys. Chem. B 1998, 102, 3099. (15) Berry, F. J.; Moore, E.; Mortimer, M.; Ren, X.; Heap, R.; Slater, P.; Thomas, M. F. J. Solid State Chem. 2008, 181, 2185. (16) Pirzada, M.; Grimes, R. W.; Minervini, L.; Maguire, J. F.; Sickafus, K. E. Solid State Ionics 2001, 140, 201. (17) Minervini, L.; Zacate, M. O.; Grimes, R. W. Solid State Ionics 1999, 116, 339. (18) Lewis, G. V.; Catlow, C. R. A. J. Phys. Chem. Solids 1986, 47, 89. (19) Gale, J. D. J. Chem. Soc., Faraday Trans. 1997, 93, 629. (20) Ordejσn, P.; Artacho, E.; Soler, J. M. Phys. ReV. B 1996, 53, R10441. (21) Soler, J. M.; Artacho, E.; Gale, J. D; Garcı´a, A.; Junquera, J.; Ordejo´n, P.; Sa´nchez-Portal, D. J. Phys.: Condens. Matter 2002, 14, 2745. (22) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865. (23) Troullier, N.; Martins, J. L. Phys. ReV. B 1991, 43, 1993. (24) Bjo¨rketun, M. E.; Sundell, P. G.; Wahnstrom, G. Phys. ReV. B 2007, 76, 054307. (25) Sundell, P. G.; Bjo¨rketun, M. E.; Wahnstrom, G. Phys. ReV. B 2006, 73, 104112. (26) Catlow, C. R. A.; Corish, J.; Jacobs, P. W. M.; Lidiard, A. B. J. Phys. C: Solid State Phys. 1981, 14, L121. (27) Zacate, M. O.; Minervini, L.; Bradfield, D. J.; Grimes, R. W.; Sickafus, K. E. Solid State Ionics 2000, 128, 243. (28) Bevillon, E.; Chesnaud, A.; Wang, Y.; Dezanneau, G.; Geneste, G. J. Phys.: Condens. Matter 2008, 145217. (29) Maekawa, T.; Kurosaki, K.; Yamanaka, S. J. Alloys Compd. 2006, 416, 214. (30) Allan, N. L.; Barron, T. H. K.; Bruno, J. A. O. J. Chem. Phys. 1996, 105, 8300. (31) Shima, D.; Haile, S. M. Solid State Ionics 1997, 97, 443. (32) Geneste, G.; Morillo, J.; Finocchi, F. J. Chem. Phys. 2005, 122, 174707. (33) Kantorovich, L. N.; Gillan, M. J. Surf. Sci. 1997, 374, 373. (34) Shannon, R. D. Acta Crystallogr. 1976, A32, 751. (35) Islam, M. S.; Slater, P. R.; Tolchard, J. R.; Dinges, T. Dalton Transactions 2004, 3061. (36) Wu, J.; Davies, R. A.; Islam, M. S.; Haile, S. M. Chem. Mater. 2005, 17, 846. (37) Khan, M. S.; Islam, M. S.; Bates, D. R. J. Mater. Chem. 1998, 8, 2299. (38) Wu, J.; Li, L. P.; Espinosa, W. T. P.; Haile, S. M. J. Mater. Res. 2004, 2366. (39) Mogensen, M.; Lybye, D, Bonanos, N.; Hendriksen, P. V.; Poulsen, F. W. In International conference on ionic and mixed conducting ceramics; Ramanarayanan, T. A., Worrell, W. L., Mogensen, M., Eds.; Electrochemical Society Proceedings: San Francisco, CA, 2001; Vol 28, p 15. (40) Cherry, M.; Islam, M. S.; Gale, J. D.; Catlow, C. R. A. Solid State Ionics 1995, 77, 207. (41) Mather, G. C.; Islam, M. S. Chem. Mater. 2005, 17, 1736. (42) Islam, M. S. Solid State Ionics 2002, 154, 75.

JP9054878