Proton Conduction in Acceptor-Doped BaSnO3: The Impact of the Interaction between Ionic Defects and Acceptor Impurities

Barium stannate is known as a promising proton-conducting material for clean energy applications. In this work, we elucidate the effect of the interaction of protons and oxygen vacancies with acceptor impurities on proton conduction in acceptor-doped BaSnO3. The analysis relies on our theoretical developments in hydration and proton hopping in proton-conducting perovskites. The transport theory, based on the master equation and effective medium approximation, provides the analytical description of hopping conduction considering the effects of disorder and changes in the potential energy landscape for protons caused by acceptor impurities. Using the proposed approach, we establish the dependence of the proton mobility and conductivity on the energies of the acceptor-bound states of ionic defects and external conditions. It is shown that the considered interactions can substantially affect the effective activation energies and prefactors of these transport coefficients. We also demonstrate that the correlation between the ionic radius rA of an acceptor impurity and the energies of its interaction with ionic defects leads to a non-monotonic dependence of the proton conductivity on rA. The obtained results are in reasonable agreement with the experimental data on the bulk conductivity of BaSnO3 doped with different acceptors.


Introduction
Acceptor-doped proton-conducting oxides are garnering significant attention due to their potential use in clean energy applications such as protonic ceramic fuel cells and electrolyzers [1][2][3][4]. Acceptor impurities, required for the hydration and appearance of protonic charge carriers, can substantially affect various properties of proton-conducting oxides. Tuning the properties important for applications by optimal acceptor doping is one of the key issues in the development of these materials. The influence of acceptor impurities on hydration and proton transfer has been extensively studied in both experimental [5][6][7][8][9][10][11] and theoretical [12][13][14][15][16][17][18][19][20][21] works. Abundant evidence has been obtained showing that the interaction between acceptor impurities and ionic defects is of fundamental importance for these phenomena in proton-conducting oxides. The main theoretical results on the influence of such interactions on proton transport were obtained by computer simulationsin particular, by the Monte Carlo method [18][19][20][21]. In our recent study [22], we proposed an analytical theory of proton conduction in acceptor-doped perovskites accounting for the fundamental effects caused by acceptor impurities (disorder, acceptor-bound defect states, changes in the potential energy landscape for proton hopping, percolation effect). This approach, which relies on the master equation for proton hopping and effective medium approximation, allowed us to describe the experimental data on the proton conductivity of BaZr 1-x R x O 3-δ [22] (R here and below denotes an acceptor impurity).
In this work, based on the theoretical developments in hydration [8,16] and proton hopping conduction [22], we explore the impact of the interaction between impurities and ionic defects on proton transport in acceptor-doped BaSnO 3 . High proton conductivity and chemical stability make acceptor-doped BaSnO 3 potentially interesting as a proton-conducting material [7,23,24]. Due to its high electron mobility, n-type BaSnO 3 is a promising material for electronic applications [25][26][27]. We have recently demonstrated that barium stannate is a good model object to analyze the role of the interaction of protons and oxygen vacancies with acceptor impurities in the hydration of proton-conducting oxides [8]. Taking into account this interaction allowed us to explain the effect of dopant type and concentration on the observed hydration behavior of BaSn 1-x R x O 3-δ [8].
Here, we elucidate the effect of the considered interactions on the dependence of the proton mobility u H and conductivity σ H , as well as their effective activation energies and prefactors, on external conditions for acceptor-doped BaSnO 3 . The relationships between the studied transport characteristics and the trapping energies of ionic defects by acceptor impurities are established. It is shown that the low-temperature behavior of σ H and its effective activation energy E σ a is determined by the proton trapping energy, while, at higher temperatures, σ H and E σ a depend on the trapping energies of both protons and oxygen vacancies. This, in particular, can alter the order in which the values of σ H and E σ a corresponding to different impurities increase as temperature changes. To compare our findings with experiments, we exploit the values of the trapping energies for specific impurities determined by the DFT method [28]. The obtained theoretical results are in reasonable agreement with the experimental data on the bulk proton conductivity of BaSn 1-x R x O 3-δ [7,29].
In addition, we analyze the implications of the trapping effect for the dependence of the proton conductivity on the ionic radius r A of the acceptor impurity in BaSnO 3 . The calculated dependence σ H (r A ) is non-monotonic, in accordance with the experimental observations for BaSn 1-x R x O 3-δ and other acceptor-doped perovskites.

Hydration
To describe the hydration of barium stannate, we exploit the statistical approach recently proposed to elucidate the effect of acceptor-bound complexes of ionic defects on the hydration and defect thermodynamics of ABO 3 perovskites [16]. Herein, we consider hydration taking into account two-particle complexes formed by acceptor-bound protons and oxygen vacancies. In this model, there exist two types of states for ionic defects. These states, bound and free, correspond to oxygen sites located in the vicinity of acceptor impurities and away from them, respectively. In our previous studies [8,16], we demonstrated that the implemented approach allows one to correctly describe the hydration of acceptor-doped perovskites BaZrO 3 , BaCeO 3 and BaSnO 3 .
In the exploited model [16], the concentration of protons c H in cubic perovskites AB 1-x R x O 3-δ in the case of moderate dopant content can be written as Here, p H 2 O is the partial pressure of water vapor (in atm); c R is the concentration of acceptor impurities (per formula unit); K hydr is the equilibrium constant of hydration, given by: where ∆H 0 hydr and ∆S 0 hydr are the enthalpy and entropy of hydration in the absence of trapping; K trap hydr is the component of the equilibrium constant associated with trapping. Within the adopted approach, K trap hydr can be expressed as [16] K trap hydr = where ∆E H and ∆E V are the trapping energies of protons and oxygen vacancies, defined as the difference between the formation energies in free and bound states; 2 are the probabilities that an oxygen site is located near acceptor impurities and away from them, respectively. Figure 1 shows the dependence of K trap hydr on the proton trapping energy ∆E and the ratio ∆E V /∆E H . There are two regions of the ∆E H and ∆E V values, in which hydration is enhanced (K trap hydr > 1) or suppressed (K trap hydr < 1) by trapping. The boundary between these two regions is determined by K trap hydr = 1. A detailed thermodynamic analysis of hydrogen dissolution taking into account acceptor-bound states of ionic defects is given in Reference [16].
where ΔEH and ΔEV are the trapping energies of protons and oxygen vacancies, defined as the difference between the formation energies in free and bound states; pb = 1 − (1 − cR) 2 and pf = (1 − cR) 2 are the probabilities that an oxygen site is located near acceptor impurities and away from them, respectively. Figure 1 shows the dependence of hydr trap on the proton trapping energy ΔE and the ratio ΔEV/ΔEH. There are two regions of the ΔEH and ΔEV values, in which hydration is enhanced ( hydr trap > 1) or suppressed ( hydr trap < 1) by trapping. The boundary between these two regions is determined by hydr trap = 1. A detailed thermodynamic analysis of hydrogen dissolution taking into account acceptor-bound states of ionic defects is given in Reference [16]. Trapping-related component of the equilibrium constant of hydration as a function of the proton trapping energy ΔEH/kT and the ratio ΔEV/ΔEH for an acceptor-doped perovskite AB0.9R0.1O3−δ. The points on the surface indicate the boundary separating the regions, in which hydration is enhanced ( hydr trap > 1) or suppressed ( hydr trap < 1) by trapping.

Proton Transport
The consideration of proton transport in this study relies on our recently developed analytical description of proton-hopping conduction in proton-conducting oxides [22]. Let us outline the physical model and the main assumptions underlying this approach.
The proton migration at elevated temperatures in the studied cubic perovskites AB1-xRxO3-δ is considered to be the result of thermally activated hopping between neighboring oxygen sites [19,22]. To analyze the effect of acceptor impurities on conduction, we consider two models of the potential energy landscape for proton hopping [22]. The first model implies that acceptor impurities deepen potential wells for protons at the nearest oxygen sites (bound states), but have little effect on the saddle point energies for proton inter-site transitions (Figure 2a). In the second model, impurities considerably reduce both the proton energy at the nearest oxygen sites and the saddle point energy for transitions between neighboring bound sites ( Figure 2b).

Proton Transport
The consideration of proton transport in this study relies on our recently developed analytical description of proton-hopping conduction in proton-conducting oxides [22]. Let us outline the physical model and the main assumptions underlying this approach.
The proton migration at elevated temperatures in the studied cubic perovskites AB 1-x R x O 3-δ is considered to be the result of thermally activated hopping between neighboring oxygen sites [19,22]. To analyze the effect of acceptor impurities on conduction, we consider two models of the potential energy landscape for proton hopping [22]. The first model implies that acceptor impurities deepen potential wells for protons at the nearest oxygen sites (bound states), but have little effect on the saddle point energies for proton inter-site transitions (Figure 2a). In the second model, impurities considerably reduce both the proton energy at the nearest oxygen sites and the saddle point energy for transitions between neighboring bound sites ( Figure 2b).  Under equilibrium, the rate of the thermally activated transition from the occupied site i to the empty site j is where the superscript "0" denotes the equilibrium value; Qij is the potential barrier for the jump i → j and νi is the prefactor (with dimension of frequency), which is assumed to be the same for all sites: νb = νf = ν. The barriers Qij for the different types of sites in the pairs (i, j) are: Qff = Qfb = Q, Qbf = Q + ΔEH and Qbb = Q + ΔEH − ΔQ. For the potential energy landscapes depicted in Figure 2a,b, the parameter ΔQ takes the values of 0 and ΔEH, respectively. The energy landscape with ΔQ = 0, corresponding to the known lattice gas model with site-energy disorder, was previously used in computer simulations and interpretation of the experimental data on proton dynamics in proton-conducting oxides (see, e.g., Refs. [5,19,20]). The model with ΔQ > 0 was also exploited in several works: in Monte Carlo simulations of proton transport and the interpretation of nuclear magnetic resonance data [9,19].
Proton hopping in our work [22] is considered to be governed by the standard master equation, which, in the mean field approximation, gives the rate equation for the occupation probability fi of site i. Under an external electric field, both the rate Wij of the transitions i → j and the occupation probability fi deviate from their equilibrium values 0 and 0 , resulting in a current of proton charge carriers. The calculation of this current is a complex problem due to the effects of disorder and different types of inter-site transitions.
Our approach [22] to the analysis of proton conduction is based on the mapping of the hopping problem onto the random resistor network of Miller and Abrahams [30], and treating it using effective medium theory (see, e.g., Ref. [31]). The local conductances gij between pairs of neighboring oxygen sites (i, j = f, b) corresponding to our problem can be written as follows [22]: where 0 and 0 are the equilibrium occupation probabilities for free and bound sites.
To obtain the expressions for the conductances (5) and (6), we exploited the detailed bal- Under equilibrium, the rate of the thermally activated transition from the occupied site i to the empty site j is where the superscript "0" denotes the equilibrium value; Q ij is the potential barrier for the jump i → j and ν i is the prefactor (with dimension of frequency), which is assumed to be the same for all sites: ν b = ν f = ν. The barriers Q ij for the different types of sites in the pairs (i, j) are: Q ff = Q fb = Q, Q bf = Q + ∆E H and Q bb = Q + ∆E H − ∆Q. For the potential energy landscapes depicted in Figure 2a,b, the parameter ∆Q takes the values of 0 and ∆E H , respectively. The energy landscape with ∆Q = 0, corresponding to the known lattice gas model with site-energy disorder, was previously used in computer simulations and interpretation of the experimental data on proton dynamics in proton-conducting oxides (see, e.g., Refs. [5,19,20]). The model with ∆Q > 0 was also exploited in several works: in Monte Carlo simulations of proton transport and the interpretation of nuclear magnetic resonance data [9,19].
Proton hopping in our work [22] is considered to be governed by the standard master equation, which, in the mean field approximation, gives the rate equation for the occupation probability f i of site i. Under an external electric field, both the rate W ij of the transitions i → j and the occupation probability f i deviate from their equilibrium values W 0 ij and f 0 i , resulting in a current of proton charge carriers. The calculation of this current is a complex problem due to the effects of disorder and different types of inter-site transitions.
Our approach [22] to the analysis of proton conduction is based on the mapping of the hopping problem onto the random resistor network of Miller and Abrahams [30], and treating it using effective medium theory (see, e.g., Ref. [31]). The local conductances g ij between pairs of neighboring oxygen sites (i, j = f, b) corresponding to our problem can be written as follows [22]: where f 0 f and f 0 b are the equilibrium occupation probabilities for free and bound sites. To obtain the expressions for the conductances (5) and (6), we exploited the detailed balance condition and Boltzmann statistics ( f 0 i 1) for protons. The latter is possible because, for the considered moderate dopant content, we can neglect the prohibition for several protons to occupy the same oxygen sites simultaneously [19,22].
In the effective medium approximation, the effective conductance g eff of the network of randomly distributed resistors is determined by the equation [31] where w(g) is the probability distribution function for g ij values; z is the coordination number for the network sites (z = 8 in our case). Within the considered model with bound and free sites, where p lm is the probability that two nearest neighboring oxygen sites are of types l and m (l, m = b, f ). For the adopted uniform distribution of acceptor impurities, Equations (5)-(8) have an exact analytical solution for g eff . The corresponding expression for the macroscopic conductivity is [22] Here, V 0 is the volume of the formula unit; u H is the proton mobility and u 0 H is the proton mobility in the absence of the interaction between protons and acceptor impurities: where A u is the prefactor. The component of the proton mobility associated with the interaction of ionic defects with acceptor impurities is defined as u where p bb = c R 1 + c R − c 2 R is the probability that two nearest neighboring oxygen sites are of type b. The above expression for σ H corresponds to the diagonal component of the conductivity tensor in a crystal with cubic symmetry. In Equation (9) and below, the tensor indices are omitted.
Consider the main features of the proton mobility behavior as a function of dopant content c R for the adopted models of the potential energy landscape. The results for the second type of the landscape are given for ∆Q = ∆E H ; however, the behavior of u H (c R , ∆Q/kT, ∆E H /kT) is quite general [22].
At ∆Q = 0, M trap = 1 and the expression for the proton mobility simplifies: In this case, u H decreases with increasing the concentration of acceptors c R and the proton trapping energy ∆E H (see Figure 3a).
In the model of the potential landscape with ∆Q = ∆E H , the barriers for the transitions b → b are significantly lower than in the first model: Q bb = Q ff (see Figure 2). Low barriers Q bb result in a non-monotonic dependence u H (c R ), as seen in Figure 3b. At low dopant concentrations, when the clusters of bound sites are isolated, u H decreases with increasing c R , as in the first model. Further increase in c R results in the overlapping of the regions of bound states and the formation of an infinite cluster at c R = c * R , where c * R is the percolation threshold. At c R > c * R , the mobility u H increases with increasing c R due to a growing contribution of proton transfer over the network of pair-connected bound sites. The c * R value in the considered problem can be found analytically ( c * R ≈ 0.21). For a more detailed discussion concerning the behavior of u H , see [22]. value in the considered problem can be found analytically ( R * ≈ 0.21). For a more detailed discussion concerning the behavior of uH, see [22].
Note that the behavior of uH, predicted within our analytical approach, agrees with the results of Monte Carlo simulations for similar potential energy landscapes [19][20][21]. It should also be noted that, in our consideration, we neglect the correlation effects caused by the interactions between charge carriers. These effects can be significant at high dopant content. However, as the Monte Carlo results showed [19,20], the influence of protonproton and proton-vacancy interactions on proton conduction is not too significant and can be neglected in many cases up to concentrations cR ~ 0.20-0.25. At moderate dopant content and reasonable values of the parameter ΔQ and trapping energies, the difference in the proton mobilities, corresponding to the considered models of the potential energy landscape, is not too significant for most problems [22]. This difference is illustrated in Figure 4 for a perovskite AB0.9R0.1O3-δ. Since altering the heights of the barriers for transitions between bound states has little effect on the results at the considered dopant content and relevant model parameters, further analysis is given for the potential energy landscape with ΔQ = 0.  Note that the behavior of u H , predicted within our analytical approach, agrees with the results of Monte Carlo simulations for similar potential energy landscapes [19][20][21]. It should also be noted that, in our consideration, we neglect the correlation effects caused by the interactions between charge carriers. These effects can be significant at high dopant content. However, as the Monte Carlo results showed [19,20], the influence of protonproton and proton-vacancy interactions on proton conduction is not too significant and can be neglected in many cases up to concentrations c R~0 .20-0.25.
At moderate dopant content and reasonable values of the parameter ∆Q and trapping energies, the difference in the proton mobilities, corresponding to the considered models of the potential energy landscape, is not too significant for most problems [22]. This difference is illustrated in Figure 4 for a perovskite AB 0.9 R 0.1 O 3-δ . Since altering the heights of the barriers for transitions between bound states has little effect on the results at the considered dopant content and relevant model parameters, further analysis is given for the potential energy landscape with ∆Q = 0. value in the considered problem can be found analytically ( R * ≈ 0.21). For a more detailed discussion concerning the behavior of uH, see [22].
Note that the behavior of uH, predicted within our analytical approach, agrees with the results of Monte Carlo simulations for similar potential energy landscapes [19][20][21]. It should also be noted that, in our consideration, we neglect the correlation effects caused by the interactions between charge carriers. These effects can be significant at high dopant content. However, as the Monte Carlo results showed [19,20], the influence of protonproton and proton-vacancy interactions on proton conduction is not too significant and can be neglected in many cases up to concentrations cR ~ 0.20-0.25. At moderate dopant content and reasonable values of the parameter ΔQ and trapping energies, the difference in the proton mobilities, corresponding to the considered models of the potential energy landscape, is not too significant for most problems [22]. This difference is illustrated in Figure 4 for a perovskite AB0.9R0.1O3-δ. Since altering the heights of the barriers for transitions between bound states has little effect on the results at the considered dopant content and relevant model parameters, further analysis is given for the potential energy landscape with ΔQ = 0.

Model Parameters for Barium Stannate
To determine the hydration and proton transfer parameters, which are independent of the interaction between defects and impurities, we used the experimental data on hydrogen dissolution and bulk conductivity for BaSn 0.875 Sc 0.125 O 3-δ [7]. The trapping energies of protons ∆E H and oxygen vacancies ∆E V required for the calculation of hydration and conductivity are taken from the DFT study [28] (see Table 1). Note that we use the same set of trapping energies as in our previous work on the hydration of acceptor-doped BaSnO 3 [8]. Trapping energies of protons and oxygen vacancies for BaSnO 3 doped with different acceptor impurities [28] ∆E Ga The trapping-independent components of the hydration enthalpy ∆H 0 hydr and entropy ∆S 0 hydr (Equation (2)) obtained by the fitting of the experimental hydration isobars for BaSn 0.875 Sc 0.125 O 3-δ are presented in Table 1. The validation of the model by comparison of the theoretical predictions, obtained using the determined parameters, with the hydration data for Y-, In-and Gd-doped BaSnO 3 is given in our previous study [8]. In general, the exploited theory provides good agreement with the experimental isobars for BaSnO 3 doped with various acceptors [8]. It is noteworthy that, according to the thermogravimetry measurements [7], the effective and nominal dopant content of the considered oxides differs. The possible reasons for such difference are discussed in more detail elsewhere [14,16]. Henceforth, we use the effective dopant content [7] for the calculation of the hydration and transport properties of barium stannate with nominal composition BaSn 0.875 Sc 0.125 O 3-δ (see Table 1).
The activation energy Q and prefactor A u of the proton mobility in the absence of the interaction with acceptor impurities (Equation (10)) were determined by fitting of the conductivity data for Sc-doped BaSnO 3 [7]. Since protons provide the dominant contribution to charge transfer at low temperatures [7], only the bulk conductivity data below~700 K were used for fitting. The parameters Q and A u were determined for the potential energy landscape for proton hopping with ∆Q = 0. As mentioned in Section 2.2, two models of the potential landscape ( Figure 2) yield similar results at the moderate dopant content and energy parameters considered in this study. Therefore, here and below, we consider only one model.

Dependence of the Proton Conductivity on the Trapping Energies of Ionic Defects
The dependence of the proton conductivity and its activation energy on the trapping energies of protons and oxygen vacancies is shown in Figure 5. The effective activation energy is calculated as ing on the relationship between the energies ΔEH and ΔEV (Figure 1), the proton conductivity σH is always reduced by the trapping effect at dopant concentrations below the percolation threshold, as illustrated in Figure 5a. At low temperatures, when the oxide is fully hydrated, an increase in the proton trapping energy ΔEH results in the reduction of σH due to a decrease in the proton mobility. However, at higher temperatures, when the oxide is partially hydrated, the dependence σH(ΔEH) can be non-monotonic (for more details, see [22]). Increasing the trapping energy of oxygen vacancies ΔEV shifts the hydration isobars to the low-temperature region [16]. Correspondingly, at a certain ΔEV value, when the contribution of oxygen vacancies to the charge neutrality condition becomes noticeable, σH begins to decrease with increasing ΔEV, as can be seen in Figure 5a. Such behavior of the proton conductivity results in a significant change in the effective activation energy (Figure 5b). The points on the surfaces in Figure 5 indicate the theoretical values of σH and a σ calculated using the trapping energies corresponding to particular acceptor impurities.  Table 1).
The effect of specific dopants on the temperature dependence of the components of the equilibrium constant hydr trap and proton mobility H trap related to the interaction of defects with acceptor impurities is illustrated in Figure 6. hydr trap depends on the trapping energies of both protons ΔEH and oxygen vacancies ΔEV. The hydration properties of an oxide can be improved by choosing an acceptor dopant with maximum and minimum values of ΔEH and ΔEV, respectively. Figure 6a shows that maximum equilibrium constant is expected for Sc and Lu, while La provides the worst hydration among the considered dopants. In contrast to hydr trap , the proton mobility under the considered conditions is determined only by the trapping of protons and decreases with increasing ΔEH (Figure 6b). Accordingly, the highest proton mobility is expected for dopants with the lowest values of ΔEH. Note that for dopants with large ionic radii, there may be an additional effect of trapping on proton conduction; this will be discussed further in Section 3.5.
In the low-temperature region, when the oxide is fully hydrated, the change in the proton conductivity σH upon replacement of an acceptor impurity is predominantly determined by the change in the ΔEH value. At higher temperatures, when the contribution  Table 1).
In contrast to hydration, which can be enhanced or suppressed by trapping depending on the relationship between the energies ∆E H and ∆E V (Figure 1), the proton conductivity σ H is always reduced by the trapping effect at dopant concentrations below the percolation threshold, as illustrated in Figure 5a. At low temperatures, when the oxide is fully hydrated, an increase in the proton trapping energy ∆E H results in the reduction of σ H due to a decrease in the proton mobility. However, at higher temperatures, when the oxide is partially hydrated, the dependence σ H (∆E H ) can be non-monotonic (for more details, see [22]). Increasing the trapping energy of oxygen vacancies ∆E V shifts the hydration isobars to the low-temperature region [16]. Correspondingly, at a certain ∆E V value, when the contribution of oxygen vacancies to the charge neutrality condition becomes noticeable, σ H begins to decrease with increasing ∆E V , as can be seen in Figure 5a. Such behavior of the proton conductivity results in a significant change in the effective activation energy (Figure 5b). The points on the surfaces in Figure 5 indicate the theoretical values of σ H and E σ a calculated using the trapping energies corresponding to particular acceptor impurities. The effect of specific dopants on the temperature dependence of the components of the equilibrium constant K trap hydr and proton mobility u trap H related to the interaction of defects with acceptor impurities is illustrated in Figure 6.  Figure 6a shows that maximum equilibrium constant is expected for Sc and Lu, while La provides the worst hydration among the considered dopants. In contrast to K trap hydr , the proton mobility under the considered conditions is determined only by the trapping of protons and decreases with increasing ∆E H (Figure 6b). Accordingly, the highest proton mobility is expected for dopants with the lowest values of ∆E H . Note that for dopants with large ionic radii, there may be an additional effect of trapping on proton conduction; this will be discussed further in Section 3.5. of oxygen vacancies to the charge neutrality condition is significant, σH depends on both energies ΔEH and ΔEV. Therefore, the order in which the conductivity value of an oxide doped with different impurities changes can differ at high and low temperatures. The calculations for each acceptor impurity were performed using the corresponding trapping energies (see Table 1).

Proton Conductivity as a Function of the Ionic Radius of the Acceptor Dopant
To further elucidate the effect of acceptor impurities on the transport properties of barium stannate, we consider the dependence of the proton conductivity σH and its effective activation energy a σ on the ionic radius rA of the dopant. In our model, this dependence results from the correlation between the radius rA and the trapping energies ΔEH and ΔEV. For barium stannate, such correlation was established by the DFT simulation [28]. Figure 7 shows the results of the calculations of σH and a σ for BaSn0.9R0.1O3-δ with different dopants. As can be seen, the conductivity increases with increasing ionic radius for small dopants (with In being an outlier at elevated temperatures) and decreases with increasing rA for large dopants (Figure 7a). The calculated activation energies alter substantially upon varying rA, and the dependence a σ (rA) is also non-monotonic (Figure 7b). . The values of σH and a σ were calculated using the trapping energies for the acceptor dopants from the DFT study [28].  Table 1).

Effect of Temperature on the Activation Energies and Prefactors of Proton Conductivity and Mobility
In the low-temperature region, when the oxide is fully hydrated, the change in the proton conductivity σ H upon replacement of an acceptor impurity is predominantly determined by the change in the ∆E H value. At higher temperatures, when the contribution of oxygen vacancies to the charge neutrality condition is significant, σ H depends on both energies ∆E H and ∆E V . Therefore, the order in which the conductivity value of an oxide doped with different impurities changes can differ at high and low temperatures.

Proton Conductivity as a Function of the Ionic Radius of the Acceptor Dopant
To further elucidate the effect of acceptor impurities on the transport properties of barium stannate, we consider the dependence of the proton conductivity σ H and its effective activation energy E σ a on the ionic radius r A of the dopant. In our model, this dependence results from the correlation between the radius r A and the trapping energies ∆E H and ∆E V . For barium stannate, such correlation was established by the DFT simulation [28]. Figure 7 shows the results of the calculations of σ H and E σ a for BaSn 0.9 R 0.1 O 3-δ with different dopants. As can be seen, the conductivity increases with increasing ionic radius for small dopants (with In being an outlier at elevated temperatures) and decreases with increasing r A for large dopants (Figure 7a). The calculated activation energies alter substantially upon varying r A , and the dependence E σ a (r A ) is also non-monotonic (Figure 7b). of oxygen vacancies to the charge neutrality condition is significant, σH depends on both energies ΔEH and ΔEV. Therefore, the order in which the conductivity value of an oxide doped with different impurities changes can differ at high and low temperatures.  Table 1).

Proton Conductivity as a Function of the Ionic Radius of the Acceptor Dopant
To further elucidate the effect of acceptor impurities on the transport properties of barium stannate, we consider the dependence of the proton conductivity σH and its effective activation energy a σ on the ionic radius rA of the dopant. In our model, this dependence results from the correlation between the radius rA and the trapping energies ΔEH and ΔEV. For barium stannate, such correlation was established by the DFT simulation [28]. Figure 7 shows the results of the calculations of σH and a σ for BaSn0.9R0.1O3-δ with different dopants. As can be seen, the conductivity increases with increasing ionic radius for small dopants (with In being an outlier at elevated temperatures) and decreases with increasing rA for large dopants (Figure 7a). The calculated activation energies alter substantially upon varying rA, and the dependence a σ (rA) is also non-monotonic (Figure 7b).  0.021 atm). The values of σH and a σ were calculated using the trapping energies for the acceptor dopants from the DFT study [28]. . The values of σ H and E σ a were calculated using the trapping energies for the acceptor dopants from the DFT study [28]. Figure 8 reports the temperature dependence of the effective activation energy E σ a and prefactor σ 0 H of the proton conductivity for Sc-, Y-, Gd-and In-doped BaSnO 3 . The prefactor is determined by (13) and the expression for E σ a is given by Equation (12).

Effect of Temperature on the Activation Energies and Prefactors of Proton Conductivity and Mobility
erials 2022, 15, x FOR PEER REVIEW 10 o Figure 8 reports the temperature dependence of the effective activation energy and prefactor σ H 0 of the proton conductivity for Sc-, Y-, Gd-and In-doped BaSnO3. T prefactor is determined by (σ H 0 ) = ( (σ H )) ( and the expression for a σ is given by Equation (12). Both parameters a σ and σ H 0 increase with decreasing temperature and attain a s uration limit in the region of complete oxide hydration (cH ≈ cR). Such behavior of a σ and σ H 0 (T) is mainly related to the variation in the proton concentration cH with tempe ture. To illustrate this relation, we calculated the temperature dependence of the acti tion energy a u and prefactor H 0 of the proton mobility using equations similar to ( and (13). As seen in Figure 8, a u and H 0 weakly depend on T, in contrast to a p nounced decline in a σ and σ H 0 with increasing temperature. As temperature decreases and cH approaches the saturation value cR, a σ (12) and (13), in the region of moderate dopant concentrations, tend to the limits (1 − 4pbb) -1 for the potential energy landscapes with ΔQ = 0 and ΔQ = ΔEH, respectively At high temperatures, when the proton concentration is low (cH << cR), the activat energy can be approximated by a σ,high T = + 0.5Δ hydr The last term in Equation (16) equals zero for ΔQ = 0 and attains a constant value, depen Both parameters E σ a and σ 0 H increase with decreasing temperature and attain a saturation limit in the region of complete oxide hydration (c H ≈ c R ). Such behavior of E σ a (T) and σ 0 H (T) is mainly related to the variation in the proton concentration c H with temperature. To illustrate this relation, we calculated the temperature dependence of the activation energy E u a and prefactor u 0 H of the proton mobility using equations similar to (12) and (13). As seen in Figure 8, E u a and u 0 H weakly depend on T, in contrast to a pronounced decline in E σ a and σ 0 H with increasing temperature.
As temperature decreases and c H approaches the saturation value c R , E σ a (12) and σ 0 H (13), in the region of moderate dopant concentrations, tend to the limits where M 0 trap is the low-temperature limit of the function M trap (11). M 0 trap = 1 and M 0 trap = (1 − 4p bb ) −1 for the potential energy landscapes with ∆Q = 0 and ∆Q = ∆E H , respectively.
At high temperatures, when the proton concentration is low (c H << c R ), the activation energy can be approximated by The last term in Equation (16) equals zero for ∆Q = 0 and attains a constant value, depending on ∆E H and c R , at high temperatures for ∆Q = ∆E H . Thus, the observed weak temperature dependence of E σ a at high T in Figure 8a is determined by the third term in Equation (16). In the case of negligible trapping, the high-temperature limit of E σ a is constant and equals Q + 0.5∆H 0 hydr . Equations (14) and (15) show that the low-temperature limit of E σ a linearly depends on the proton trapping energy ∆E H (Figure 8a), while the saturation limit of σ 0 H is the same for all dopants (Figure 8b). As a result, at low temperatures, the activation energy decreases in the order Sc > In > Gd > Y, in accordance with the ∆E H values (see Table 1). However, at higher temperatures, the trend is different since the activation energy and prefactor depend on both energies ∆E H and ∆E V . It should be noted that, outside the region of small dopant concentrations, varying c R within reasonable limits has virtually no effect on the calculated values of E σ a and σ 0 H . Figure 8a demonstrates that in order to obtain the low-temperature limit of the activation energy of the proton conductivity from the experimental dependence σ H (T), the temperature range should be chosen in the region of complete oxide hydration. However, this can be complicated providing that the oxide hydration is poor and/or the conductivity measurements are performed at elevated temperatures, when c H < c R .
It is important to note that the activation energy is usually determined within the temperature range that exceeds the region of complete oxide hydration. In this case, the obtained temperature-averaged E σ a value would be lower than the low-temperature limit determined by Equation (14). For example, the low-temperature limit of the activation energy for Sc-doped BaSnO 3 equals Q + ∆E Sc H ≈ 0.63 eV (see Figure 8a). At the same time, the theoretical value of E σ a averaged over the temperature range 500-700 K is approximately 0.48 eV, which is close to the result of Wang et al. [7]. Figure 9 shows the proton conductivity of Sc-, Y-, Gd-and In-doped BaSnO 3 calculated within our model, along with the experimental data [7,29]. The conductivity of Sc-doped BaSnO 3 is a result of the fitting procedure (see Section 3.1). The conductivity of BaSnO 3 doped with other acceptor impurities is obtained without any fitting, using the determined model parameters (Table 1) and the trapping energies from the DFT study [28]. The calculations for BaSn 0.875 R 0.125 O 3-δ were performed using the effective values c eff R of the dopant content [7]. For BaSn 0.75 R 0.25 O 3-δ , we used the nominal value c R = 0.25 since there are no hydration data for these samples in [29].

Comparison with Experimental Data
It should be noted that the conductivities of BaSn 0.875 R 0 . 125 O 3-δ and BaSn 0.75 R 0.25 O 3-δ do not differ significantly under the considered conditions. Such a weak dependence σ H (c R ) outside the regions of small and large c R values is not unusual; it was experimentally observed for other proton-conducting perovskites (see, e.g., Ref. [32]). The Monte Carlo simulations also showed that, under certain conditions, the dependence σ H (c R ) can be relatively weak at intermediate dopant concentrations [19,20]. This effect can be explained by the mutually compensating influence on the conductivity of two factors-an increase in the proton concentration c H and a decrease in the proton mobility u H with increasing c R . Such a decreasing dependence u H (c R ) at the considered moderate dopant content is caused by the trapping effect (see Figure 3 and Refs. [19,20]).
It can be seen from Figure 9 that the shifts in the calculated proton conductivity for In and Y relative to Sc generally follow the experimental data, although the slopes of the conductivity curves for In are somewhat different in the low-temperature region. The agreement between theory and experiment for Gd-doped BaSnO 3 is worse than for the oxide with other dopants. The reasons for this discrepancy in the case of Y and Gd can be partially explained by their large ionic radius, as will be seen further below. Another important factor that can lead to lower conductivity values, as compared to the theoretical results, is slow kinetics caused by the high density of the samples. In particular, it can hinder the attainment of the theoretically expected degree of hydration, especially at low temperatures. For example, the relative density of BaSn 0.875 R 0.125 O 3-δ (R = In, Gd) [7] and BaSn 0.75 In 0.25 O 3-δ [29] was above 98%. erials 2022, 15, x FOR PEER REVIEW 12 of Figure 9. Temperature dependence of the proton conductivity of (a) BaSn0.875R0.125O3−δ and BaSn0.75R0.25O3−δ in a humidified atmosphere (pH2O = 0.021 atm). The symbols correspond to the perimental data on bulk conductivity in wet Ar [7,29]. The parameters Au and Q for the proton m bility were determined by fitting of the conductivity data for Sc-doped BaSnO3 (Table 1). Solid bl lines are the fitting curves. Red, blue and green lines are the theoretical conductivities calcula using the determined parameters and DFT results for the trapping energies [28]. According to a number of DFT studies for BaSnO3 [13,28], in the case of acceptor d pants with large ionic radii, the trapping energies of protons and oxygen vacancies in first and second neighbor positions can be comparable. In order to roughly estimate t implications of this effect, we extend the trapping regions around such impurities up the second neighbors, accordingly redefining the probabilities pf and pb. The total numb of proton positions in such trapping regions is large, and we are beyond the applicabil of our theory (especially if the effect of impurities on the saddle point energies for int site transitions is significant, as in the potential landscape model with ΔQ = ΔEH). Nev theless, to demonstrate a possible trend, we provide these estimates for the potential lan scape with ΔQ = 0. As shown in Figure 9, the extension of the trapping region leads good agreement between theory and experimental data for Y-doped BaSnO3. Howev the results for Gd-doped BaSnO3 still do not agree quite well with the experiment, a hough the calculated conductivity values become closer to the experimental data.
In another study, Li and Nino [33] measured the bulk conductivities for BaSn0.9R0.1 δ (R = In, Lu, Er, Y, Gd) under oxidizing and reducing conditions. The results indicate th the order in which the conductivity corresponding to different acceptor impurit changes is quite different from that obtained by Wang et al. [7,29]. However, the exter conditions of the conductivity measurements in Reference [33] differed from those in experiments [7,29], which were carried out in a humidified Ar atmosphere. According the EMF measurements [33], a significant contribution to the total conductivity in oxid ing and reducing atmospheres is provided by electronic charge carriers. Since this con bution can also differ for samples with different dopants, a comparison of our theoreti The symbols correspond to the experimental data on bulk conductivity in wet Ar [7,29]. The parameters A u and Q for the proton mobility were determined by fitting of the conductivity data for Sc-doped BaSnO 3 (Table 1). Solid black lines are the fitting curves. Red, blue and green lines are the theoretical conductivities calculated using the determined parameters and DFT results for the trapping energies [28]. Dopant concentrations c R are taken to be equal to the (a) effective [7] and (b) nominal, c R = 0. 25 According to a number of DFT studies for BaSnO 3 [13,28], in the case of acceptor dopants with large ionic radii, the trapping energies of protons and oxygen vacancies in the first and second neighbor positions can be comparable. In order to roughly estimate the implications of this effect, we extend the trapping regions around such impurities up to the second neighbors, accordingly redefining the probabilities p f and p b . The total number of proton positions in such trapping regions is large, and we are beyond the applicability of our theory (especially if the effect of impurities on the saddle point energies for inter-site transitions is significant, as in the potential landscape model with ∆Q = ∆E H ). Nevertheless, to demonstrate a possible trend, we provide these estimates for the potential landscape with ∆Q = 0. As shown in Figure 9, the extension of the trapping region leads to good agreement between theory and experimental data for Y-doped BaSnO 3 . However, the results for Gd-doped BaSnO 3 still do not agree quite well with the experiment, although the calculated conductivity values become closer to the experimental data.
In another study, Li and Nino [33] measured the bulk conductivities for BaSn 0.9 R 0.1 O 3-δ (R = In, Lu, Er, Y, Gd) under oxidizing and reducing conditions. The results indicate that the order in which the conductivity corresponding to different acceptor impurities changes is quite different from that obtained by Wang et al. [7,29]. However, the external conditions of the conductivity measurements in Reference [33] differed from those in the experiments [7,29], which were carried out in a humidified Ar atmosphere. According to the EMF measurements [33], a significant contribution to the total conductivity in oxidizing and reducing atmospheres is provided by electronic charge carriers. Since this contribution can also differ for samples with different dopants, a comparison of our theoretical results with the bulk conductivity data [33] would be incorrect.
We now turn to the dependence of the proton conductivity on the ionic radius of the dopant r A . To compare the results with experiments, the bulk conductivity data for Sc-, In-, Y-and Gd-doped BaSnO 3 [7,29] and the calculated conductivities are plotted as a function of r A in Figure 10. It is seen that the behavior of the theoretical conductivities correlates well with the experimental data, including the downward shift for In at elevated temperatures. The estimates of σ H (r A ) for the oxide with large dopants Y and Gd were also made using the extended trapping regions around acceptor impurities, see above. The obtained values of σ H (blue points in Figure 10) are shifted downwards and closer to the experimental data.  A non-monotonic dependence σH(rA), similar to that predicted by our model for acceptor-doped BaSnO3 (see Figures 7a and 10), was also experimentally observed for perovskites BaZrO3 [10,17] and BaCeO3 [34] doped with different acceptor impurities.

Conclusions
We have applied our recently developed statistical theory of proton-hopping conduction in oxide perovskites to reveal the role of the interaction between ionic defects and acceptor impurities in proton transport in acceptor-doped barium stannate. Accounting for this interaction within the proposed approach allowed us to explain the observed behavior of the bulk proton conductivity σH of BaSn1-xRxO3-δ. The experimental dependences of σH on temperature, type of acceptor impurity and its ionic radius are described reasonably well. A number of results concerning the influence of impurities on proton conductivity and mobility are quite general for perovskites with moderate dopant content. For example, in the low-temperature region of complete oxide hydration, the main effect of the interaction between acceptor impurities and ionic defects on the behavior of σH is due to the proton trapping. In the low-temperature limit, the effective activation energy a σ of σH increases linearly with increasing the proton trapping energy. At higher temperatures, a σ depends on the trapping energies of both protons and oxygen vacancies and decreases with increasing temperature. The predicted non-monotonic dependence of σH on the dopant ionic radius is observed not only for BaSnO3, but also for other acceptor-doped perovskites. Our findings contribute to the understanding of the role of acceptor impurities in proton transport in oxides and can be useful for selecting optimal acceptor doping for proton-conducting materials.
Funding: This research received no external funding. A non-monotonic dependence σ H (r A ), similar to that predicted by our model for acceptor-doped BaSnO 3 (see Figures 7a and 10), was also experimentally observed for perovskites BaZrO 3 [10,17] and BaCeO 3 [34] doped with different acceptor impurities.

Conclusions
We have applied our recently developed statistical theory of proton-hopping conduction in oxide perovskites to reveal the role of the interaction between ionic defects and acceptor impurities in proton transport in acceptor-doped barium stannate. Accounting for this interaction within the proposed approach allowed us to explain the observed behavior of the bulk proton conductivity σ H of BaSn 1-x R x O 3-δ . The experimental dependences of σ H on temperature, type of acceptor impurity and its ionic radius are described reasonably well. A number of results concerning the influence of impurities on proton conductivity and mobility are quite general for perovskites with moderate dopant content. For example, in the low-temperature region of complete oxide hydration, the main effect of the interaction between acceptor impurities and ionic defects on the behavior of σ H is due to the proton trapping. In the low-temperature limit, the effective activation energy E σ a of σ H increases linearly with increasing the proton trapping energy. At higher temperatures, E σ a depends on the trapping energies of both protons and oxygen vacancies and decreases with increasing temperature. The predicted non-monotonic dependence of σ H on the dopant ionic radius is observed not only for BaSnO 3 , but also for other acceptor-doped perovskites. Our findings contribute to the understanding of the role of acceptor impurities in proton transport in oxides and can be useful for selecting optimal acceptor doping for proton-conducting materials.