Molecular dynamics studies of oxide ion transport in Sr-doped LaFeO₃: role of cationic environments and cooperativity

Abstract

Classical molecular dynamics simulation is employed to study La_1−xSr_xFeO_3−x/2 at 1200 K over a range of dopant concentrations, x = 0.1 to 0.6. Sr ions tend to improve the energetics of the oxide ion sites, resulting in a higher fraction of vacancies in the vicinity of La. The oxide ion migration in the system involves predominantly intra-octahedral hops along the edges of the FeO₆ octahedra. This ion migration is controlled by a triangular bottleneck of cations, formed by two La/Sr ions and one Fe ion, appearing midway between two neighboring oxide ion sites. It is noted that these bottlenecks pose higher barriers as more Sr ions are introduced. This increase in the microscopic energy barriers for oxide ion migration in the system corroborates the observed slowdown of ion transport with dopant concentration. The study also elucidates a dynamic correlation between the mobile species and the cationic framework, wherein as an oxide ion approaches the bottleneck, the cations move apart, increasing the cross-sectional area of the bottleneck.

---

1 Introduction

One of the key challenges to our society in this century is to reduce our dependency on fossil fuels. Carbon-neutral power generation can be achieved by resorting to energy sources like nuclear, solar, wind, and hydroelectric generators. However, when it comes to portable ‘green’ energy devices, we are limited to the choice of batteries and fuel cells. Fuel cells are primary energy devices, as opposed to batteries, that convert chemical energy in the fuel to electricity. Solid oxide fuel cells (SOFCs) are particularly attractive among their peers, owing to their scalability, longer life, high fuel efficiency, low emission, and fuel flexibility (permitting hydrogen as well as hydrocarbons).^1–10

The main components of a SOFC are the cathode, anode, and electrolyte. An ideal electrolyte must be a good ionic conductor but highly resistive to electronic conduction, whereas electrodes, i.e. cathode and anode, can be electronic as well as mixed ionic electronic conductors (MIECs). One of the shortcomings of pure electronic conducting cathode materials is that the oxygen reduction reaction is limited only to the triple phase boundary (TPB), the region where the cathode, electrolyte, and air meet, which reduces the efficiency of the fuel cell. MIEC cathodes have the advantage that the oxygen reduction reaction is not limited to TPB but over the entire surface of the cathode exposed to fuel.^11–15 Typically, SOFCs suffer from their requirement of a high operating temperature, in the range of 1100–1300 K, which causes material degradation and shortening of their lifespan.^14,16–18 Thus, in recent years, more attention has been paid to lowering the operating temperature of SOFCs to the intermediate regime (900–1100 K). One of the key challenges in lowering the operating temperature is the reduced electrocatalytic rate at the cathode.^{12,13,15–17,19,20} However, this deficiency can be better managed if the ionic conductivity of the cathode can be improved (that is, by overcoming the limitations due to TPB by MIEC). This strategy offers cost reduction and efficiency of SOFCs, which improves their market competitiveness essential for scaling up commercialization.

Over the past two decades, ABO₃ based perovskites, where large cations (such as La, Pr, Sm, Nd, Gd, etc.) occupy the A sites and transition metals (such as Mn, Fe, Co, etc.) occupy the B sites, have been extensively studied for their potential applications in SOFCs as electrolytes and cathode materials.^8,13,21–31 In general, ion transport in solids occurs primarily through vacancy migration and interstitial mechanisms.³² Oxide ion transport in perovskite materials occurs through the vacancy migration mechanism, where an oxide ion hops to one of its adjacent vacant sites, leaving the original site available for other ions to occupy. These vacancies are typically generated by aliovalent doping of the parent material.^8,33 Lanthanum ferrites (LaFeO₃), when doped with Sr, generate the necessary oxide ion vacancies for ion migration, while promoting favorable p-type electrical conductivity.^34,35 Reportedly, Sr-doping also improves the chemical stability of these ferrites;³⁶ however, their oxide ion (O²⁻) conductivity is rather low for SOFC applications at intermediate temperatures.

Computational studies in the past have offered valuable insights into ion conduction mechanisms in superionic solids.^37–46 In one of the earliest computational studies on lanthanum ferrites, Jones and Islam⁴⁷ examined the ion migration channels, the nature of bottlenecks for ion migration, and the effect of dopants on the microscopic barriers, using empirical-potential based calculations for the orthorhombic phase of LaFeO₃, with the lattice treated static. Since then many density functional theory (DFT) based studies have investigated the role of the Sr dopant in oxygen vacancy formation and migration energy in LaFeO₃ giving useful insights.^48–51 Very recently, the influence of grain boundaries on oxide ion diffusion was examined using a molecular dynamics study of orthorhombic and cubic phases of polycrystalline Sr substituted LaFeO₃ by Bonkowski et al.⁵²

Fig. 1 The rhombohedral structure (R3̄c) of LaFeO₃ showing corner shared FeO₆ octahedra. The red balls represent oxygen, with Fe occupying the center of the octahedra.

Over the typical operating temperatures (∼1200 K) these solids generally stabilize in the rhombohedral phase for a wide range of Sr concentrations (x = 0.1 to 0.6).⁵³ However, there are very few theoretical studies that examined the rhombohedral phases. Furthermore, most of the theoretical works available in the literature focus on the energetics of vacancy creation and migration, and hardly any dynamical properties were reported. In fact, except for the study of Bonkowski et al.,⁵² computational studies have been performed within a frozen cationic sub-lattice approximation. However, at the typical operating temperature of a fuel cell, the dynamics of the framework cations can significantly impact the oxide ion diffusion. In a recent review on the appropriateness of different computational techniques in the study of material properties, Schwarz et al. emphasized the need to invoke a dynamical approach as opposed to static calculations.⁵⁴

Thus, the present study explores the rhombohedral phase (Fig. 1) of La_1−xSr_xFeO_3−x/2 for Sr concentrations over x = 0.1 to 0.6 employing classical molecular dynamics simulations at 1200 K. The choice of this rather high-temperature simulation is to enhance the ion dynamics for better statistics. The focus of the article is on the microscopic structural and dynamical properties for a better understanding of the mechanism governing the oxide ion transport in these systems and the various factors that influence it.

2 Computational details

Classical molecular dynamics (MD) simulations of La_1−xSr_xFeO_3−x/2 are carried out in an isothermal–isobaric (NPT) ensemble at 1200 K for various Sr dopant concentrations (x = 0.1–0.6). Interactions between the ions are modeled using short-range Buckingham potential coupled with long-range Coulombic interactions,where q_i and q_j are respectively the formal charges of the ith and jth ions, and r is the separation between ion pairs. A_ij and ρ_ij are, respectively, the strength and range of the short-range repulsion between ions, and C_ij represents the van der Waals attraction between the pairs. The parameters for the short-range inter-atomic interactions (see Table 1) are taken from a previous study,⁴⁷ which was originally proposed as a core–shell model, but is replaced by a simpler core–core model in the present study to reduce the computational cost. The cation–cation interactions are considered to be purely Coulombic.

Table 1 The parameters of the interaction potential employed in the present study⁴⁷

Pairs	A ij (eV)	ρ ij (Å)	C ij (eV Å^6)
La^3+⋯O^2−	1545.21	0.3590	0
Sr^2+⋯O^2−	1400.00	0.3500	0
Fe^3+⋯O^2−	1156.36	0.3299	0
O ^2−⋯O^2−	22764.30	0.1490	43

---

NPT-MD simulations are carried out for 20 ns, after dedicating an initial 4 ns for equilibration. A time step of 2 fs is used for integration, and trajectories are stored at intervals of 1 ps. The simulation cell consists of 5 × 5 × 2 rhombohedral unit cells, each containing six formula units of La_1−xSr_xFeO_3−x/2, such that the undoped system has 1500 number of atoms (300-La, 300-Fe, 900-O). The number of oxide ions in the simulation cell varies from 885 to 810, with the corresponding oxygen vacancies varying from 15 to 90, respectively for x = 0.1 and 0.6. The NPT-MD simulations were started from the experimental structure,⁵⁵ having unit-cell parameters, a = b = 5.51107 Å and c = 13.41578 Å. Periodic boundary conditions and Ewald summation techniques are employed to mimic the bulk system. The Sr doping is performed at the La sites, picked randomly, and charge neutrality is maintained by removing one oxide ion for every two Sr inserted. Simulations are carried out using the LAMMPS package⁵⁶ and the trajectories are analyzed using in-house software.

3 Results and discussion

3.1 Structure

In order to establish the reliability of the potential parameters, the structural properties of the system are investigated. The typical time evolution of the lattice parameters of the rhombohedral unit cell of La_1−xSr_xFeO_3−x/2 during the NPT-MD simulations at 1200 K is presented in Fig. S1 in the ESI,† for select compositions. As shown in Fig. 2, the average lattice parameters of the systems increase marginally with Sr doping, and the values are within 3% of the experimental values.

Fig. 2 The average lattice parameters of a rhombohedral cell of La_1−xSr_xFeO_3−x/2, where x = 0.1 to 0.6, from NPT-MD simulations at 1200 K.

The spatial distributions of the atoms in the system are examined using the radial distribution function (RDF), calculated as

where V is the volume of the simulation box, N_α and N_β are the total number of atoms of type α and type β respectively, and n_i is the number of atoms of type β in a spherical shell of thickness Δr at a distance r from the reference atom α.

Fig. 3 shows RDFs between La–O, Sr–O, Fe–O, and O–O pairs for different compositions x = 0.1 to 0.6 at 1200 K. The vertical bars correspond to RDFs calculated from the experimentally reported structure for the undoped system.⁵⁵ As seen the dynamic RDFs of different pairs are in good agreement with their experimental ones, confirming the structural integrity of the simulated system. The broadening in the dynamic RDFs is because of the thermal vibrations of atoms about their mean positions and the successive peaks in the RDFs represent the consecutive position of neighboring atoms of that particular species with respect to the reference atom. It is observed that the Sr–O distance is slightly higher than the La–O distance across the composition, which is attributed to the larger ionic radius of Sr²⁺ and its relatively weaker Coulombic attraction with O²⁻.

Fig. 3 The radial distribution function (RDF) between select species in La_1−xSr_xFeO_3−x/2, where x = 0.1 to 0.6, from NPT-MD simulations at 1200 K. The black bars represent the corresponding X-ray locations. The plots are displaced along the y-axis for clarity.

3.2 Ion transport

As noted earlier, the substitution of Sr²⁺ at the La³⁺ sites results in the creation of oxygen vacancies in the matrix, which are quintessential for oxide ion transport in these systems. It shall be noted that the interstitial sites are not energetically favorable for the oxide ions due to their bonding characteristic with the cationic framework.⁵⁷ Induced by the vacancies, the thermally activated oxide ions hop off from one lattice site to another vacant site, leaving behind a vacancy that is now accessible to other oxygens in the neighborhood. As this process continues, the vacancies essentially migrate in the system in a complementary fashion, though the nature of the dynamics of the ion and its residual vacancies could be different.^33,58–60 We shall, however, note that while the vacancies are essential for oxide ion transport, other factors, such as the nature of the dopant ion (Sr), also start controlling the overall transport, as detailed later in the article.

The gross transport properties of the system, such as the self-diffusivity and ionic conductivity, can be obtained by calculating the mean squared displacement (MSD) of atoms with time, given by

MSD(t) = 〈|r̄(t + τ) − r̄(τ)|²〉 (3)

where r̄(t) is the position vector of the particle at time t and the angular bracket indicates averaging over all the ions and for numerous time origins τ. Shown in Fig. 4 are the MSDs of oxide ions in La_1−xSr_xFeO_3−x/2, for different compositions x = 0.1 to 0.6 at 1200 K. It is observed that MSDs of oxide ions decrease with increasing Sr concentration. The MSDs of the cations in the system (as demonstrated in Fig. S2 in the ESI†) remain flat for all compositions, indicating that their dynamics are limited to vibrations about the mean position, providing further evidence of the structural stability of the simulated system.

Fig. 4 Mean squared displacement of oxide ions in La_1−xSr_xFeO_3−x/2, where x = 0.1 to 0.6, from NPT-MD simulations at 1200 K.

The self-diffusion coefficient of oxide ions is estimated by fitting to the Einstein relation,

Furthermore, oxide ion conductivity is measured using the Nernst–Einstein equation given in ref. 58,where n is the number density of oxide ions, q is the formal charge of the oxide ions, k_B is the Boltzmann constant and T is the temperature of the system. The self-diffusivity (D) and the ionic conductivity σ of oxide ions at 1200 K for all the compositions are shown in Fig. 5.

Fig. 5 The self-diffusivity of oxide ions (black) and the corresponding conductivity (blue) for La_1−xSr_xFeO_3−x/2, where x = 0.1 to 0.6 from NPT-MD simulations at 1200 K.

As seen in Fig. 5, the calculated self-diffusivity and oxide ion conductivity decrease monotonically with increasing Sr concentration. This behavior can be understood in terms of the changes in the energetics of the local oxide ion environment with increasing Sr concentration, discussed in the later sections. It shall be noted that the oxide ion conductivities in earlier experimental reports are quite scattered, differing from one another by orders of magnitude, thus making it difficult for any direct comparison.^61–67 In fact, an earlier experimental study by Patrakeev et al. observed an increase in ionic conductivity with Sr doping, producing a maximum around x = 0.5 and dropping afterward, which is in contrast to the present observations.³⁴

3.3 Microscopic insights

In LaFeO₃, oxygen occupies the corners of the FeO₆ octahedra and each oxygen is shared by two such octahedra. Thus the structure can be imagined as a network of such octahedra linked through oxide ions. Oxygens are also four coordinated with the La sites that are partially replaced by Sr, as shown in Fig. 6. Thus, based on the number of La and Sr neighbors, five distinct oxygen environments are possible in the matrix. We choose to label these oxygen environments as Sm, where m = 0, 1, 2, 3, and 4, depending on the number of Sr present locally around an oxygen site. For example, S0 represents an oxygen environment where all the four neighbors of the oxygen are La, and S2 represents an environment that has two Sr and two La neighbors.

Fig. 6 The local cationic environment around the oxygens (red), shared by two FeO₆ octahedra (Fe is located at the center of the octahedra). The locally four coordinated cations (La³⁺/Sr²⁺) are shown as green balls.

3.3.1 Effect of Sr doping on the oxygen environment. The potential energy of the individual oxide ions shall be defined aswhere U_ij is the pairwise interaction between ith and jth atoms (given in eqn (1)). To investigate the effect of Sr dopant on oxygen environments, we calculate the MD averaged potential energy of oxide ions when occupying five different oxygen environments mentioned above (Fig. 6). We observe that these ‘site energies’ of oxide ions tend to decrease with the number of Sr in a given environment (see the left panel of Fig. 7). As seen, the difference in the energies of the five different oxygen sites (S0 to S4) starts widening with the composition. A matrix having larger undulations (such as the case with x = 0.6) in its energy landscape hampers the transport due to the Boltzmann factor.^59,68 The right panel of Fig. 7 shows the oxygen occupancies of the five different environments, S0 to S4, at 1200 K across the compositions. Here, the oxygen occupancy at a given environment is the average count of oxide ions per frame of the MD trajectory (in a spherical volume of radius 1 Å around the site), normalized by the number of such sites in the matrix for the specific composition. Consistent with the Boltzmann probability, the decrease in the energy of the oxygen sites having larger coordination numbers with Sr results in their higher oxygen occupancies. This suggests that the oxygen vacancies largely populate the La rich environments. Additionally, it shall be inferred that the oxygen vacancies are distributed randomly in the matrix, as Sr doping is carried out at random La locations.

Fig. 7 Left: Site energies of oxide ions at the five different local cationic environments, S0 to S4 (where the indices mark the Sr coordination of the sites), across the composition range, x = 0.1 to 0.6, from NPT-MD simulations at 1200 K. Right: The corresponding oxygen occupancies in these environments.

3.3.2 Migration pathways. In principle, the oxide ions can have two distinct local migration paths: (a) intra-octahedral (where oxide ions hop between two neighboring sites that are part of the same octahedra) and (b) inter-octahedral (where the path connects two neighboring oxygen sites that are part of two different octahedra). The preferred migration channel in the system is examined with the help of spatial density distribution, which maps the density of oxide ions on a fine three-dimensional mesh spanning the unit cell, employing the entire MD trajectory. Fig. 8 shows the corresponding iso-density surface of oxide ions in La_1−xSr_xFeO_3−x/2, at 1200 K for x = 0.3, at three different iso-values, blue (high), red (intermediate), and yellow (low). As seen, for larger iso-values (blue and red) the spatial density of oxygens is limited to local vibrational motion. At a lower iso-value (yellow), the oxygen density can be seen to develop along the edges of the FeO₆ octahedra, signifying intra-octahedral migration. However, the absence of density connecting possible inter-octahedral sites (highlighted in black dashed cross-wires in Fig. 8) suggests that this path is not favorable. This observation is in agreement with earlier computational studies, wherein prohibitively high barriers were observed for inter-octahedral hops in the orthorhombic phase of LaFeO₃.^47,57

Fig. 8 The diffusion pathways of oxide ions in the rhombohedral unit cell of La_1−xSr_xFeO_3−x/2, for x = 0.3, are demonstrated with the help of three different oxygen iso-densities (in the order, blue > red > yellow), from NPT-MD simulations at 1200 K. The ideal X-ray structure, with the FeO₆ octahedra (highlighted in light blue color) with the oxygen sites (red balls) and the location of the (La/Sr) cations (green balls), is shown for reference. The four inter-octahedral channels within the unit cell (visible from this perspective) are marked by dashed cross-wires in black. For the higher density iso-values, as in panels (a) and (b), the oxide ion densities are largely localized around their lattice sites. At a lower iso-value (c), the development of migration channels along the edges of the FeO₆ octahedra (intra-octahedral path) is seen. However, the density along the inter-octahedral paths remains disconnected, indicating that inter-octahedral oxygen diffusion is not favorable.

3.3.3 Oxide ion migration barriers. Having identified the oxide ion migration paths, we calculate the energy barrier of individual ion hops along the migration paths. These individual energies of oxide ions (as in eqn (6)) are shown as a function of distance from the nearest oxygen site (disregarding whether it is a departing or arrival site), up to the saddle point located midway between the two intra-octahedral oxygen sites. The individual oxide ion energies are averaged over bins along the intra-octahedral paths employing the entire NPT-MD trajectory at 1200 K. The energy profiles for the different compositions (Fig. 9) show an increase in the migration barrier with Sr concentration, justifying the observed decrease in oxygen diffusivity. As noted earlier, while the oxygen vacancies created due to Sr²⁺ doping at La³⁺ by themselves should enhance the oxide ion diffusion in the system, the larger ionic radius of the dopant Sr²⁺ leads to the observed increase in the microscopic migration barrier, which shall be demonstrated later in this section.

Fig. 9 Microscopic energy profiles of oxide ion diffusion in La_1−xSr_xFeO_3−x/2, across the compositions x = 0.1 to 0.6, at 1200 K. The distances are measured from the oxygen lattice site (at x = 0) and shown up to the transition state (saddle point), at x = 1.5.

The unanticipated dip in energy profiles, which emerges variably around 0.6–0.9 Å depending on the composition, observed in a recent DFT study as well,⁴⁹ is attributed to the local distortions of the oxygen sub-lattice, wherein oxygens temporarily relocate closer to the saddle point at instances when the adjacent oxygen site is vacant. A vacancy at the neighboring site results in a temporal imbalance of forces on an oxide ion due to the repulsion from the oxide ions in the other direction. However, a further shift is controlled due to the repulsion of bottleneck cations (detailed in the next sub-section); thus, an intermediate minimum is established. This mechanism is established by examining the potential energy profiles under two possible scenarios (presented in Fig. S3 in the ESI†): (a) when the neighboring oxygen site is occupied – the minimum at the intermediate distance vanishes; (b) but when the neighboring site is vacant – the intermediate minimum is prominently seen.

3.3.4 Cooperative nature of oxide ion hops. As demonstrated earlier in Fig. 8, oxide ion migration in LaFeO₃ occurs along the edges of the FeO₆ octahedra. During the hops, oxygen passes through a triangular bottleneck, located midway in the migration path (see illustration in Fig. 10), that is made up of two La or Sr cations and one Fe-cation, as noted in earlier studies.^47,48,69Fig. 11 illustrates the instantaneous dimension (circumradius) of the individual bottlenecks as a function of the distance of the approaching oxygen, averaged over the stored MD trajectory at 1200 K. As seen in Fig. 11, the radius of the bottlenecks (R_BN) formed by the framework cations is sensitive to the distance of the oxide ion (R_C–O) from its circumcentre (marked “C”). When an oxide ion approaches the bottleneck, the short-range repulsion between the oxygen and cations constituting the bottleneck (which dominates over the electrostatic attraction between them at shorter distances) pushes away the cations. This mechanism provides the necessary increase in the cross-sectional area of the bottleneck to facilitate the ion hops. The present study thus reveals an interesting dynamic correlation between the framework cations and the mobile oxygens. It shall be noted that in one of the recent studies of Li conducting NASICON-type solids, we observed a similar cooperative mechanism facilitating ion migration.⁷⁰ Hence, this behavior could plausibly be operative across a wide range of superionic solids in varying degrees.

Fig. 10 Schematic diagram of oxide ion diffusion from one site to another through a triangular bottleneck made up of two La/Sr (green) and one Fe (brown). The black arrow shows the oxygen diffusion path along the edges of the FeO₆ octahedra.

Fig. 11 The circumradius of the bottleneck, R_BN, as a function of the distance, r_C–O, of an approaching oxide ion for all the compositions, averaged over the NPT-MD trajectories at 1200 K. The nested image is the schematic representation of the triangular bottleneck (black line), and its circumcircle is centered at “C” and an approaching oxide ion is marked “O”. Here, the origin (x = 0) is chosen at the center of the circumcircle (“C”) and the nearest oxygen site is around x = 1.5 Å.

3.3.5 Localization of oxide ions. For a better understanding of the role of these bottlenecks in the oxide ion transport in La_1−xSr_xFeO_3−x/2 solids, we examined the energetics and dimensions of the three individual types of bottlenecks, namely (i) La–Fe–La, (ii) La–Fe–Sr and (iii) Sr–Fe–Sr. It is noted that the dimensions of the bottlenecks (see Fig. S4 in the ESI†) are, by and large, insensitive to the combination of cations forming them. However, as shown in Fig. 12, the barriers increase progressively as La is replaced by Sr. Consequently, the estimated percentage of ion hops across these three individual bottlenecks (shown in the brackets) reduces as predicted by the Boltzmann probability. The variation of energy barrier with the number of Sr ions involved in the triangular bottleneck formation was also observed by Mastrikov et al. in their DFT study.⁴⁸ The increase in the energy barrier for the Sr dominant bottlenecks is attributed to the larger size of the Sr ions compared to La. This is prompted by the fact that the dimensions of all three kinds of bottlenecks were very similar; thus, larger Sr ions reduce the window size available to oxygens for hops. Thus in La_1−xSr_xFeO_3−x/2 solids with the increase in Sr content, the number of bottlenecks with Sr corners increases, thus posing higher barriers for oxide ion migration. For instance, as the oxygen environment consists of all Sr (S4), the bottlenecks in all the eight migration channels emerging from the sites are of Sr–Fe–Sr-type causing large barriers, which results in the localization of oxide ions over longer periods. This is further compounded by the lowering of the energies of oxygen sites having a larger number of Sr neighbors (Fig. 7), noted earlier in Section 3.3.1. This results in a wider distribution of the energetics of the sites, causing larger undulations in the potential energy landscape in the matrix, hampering the oxide ion transport.^59,68 This explains the reduction in oxide ion conductivity in La_1−xSr_xFeO_3−x/2 with Sr doping.

Fig. 12 The nature of the microscopic potential energy profiles of the oxide ions while passing through the three kinds of bottlenecks in La_1−xSr_xFeO_3−x/2, namely La–Fe–La, La–Fe–Sr and Sr–Fe–Sr, is demonstrated for x = 0.6. These are averaged over the NPT-MD trajectory at 1200 K. The origin of distance is the nearest oxide ion site. The estimated percentage of ion hops (based on the counts at the saddle-point) through the individual bottlenecks is indicated in brackets.

The sluggish oxide ion mobility in Sr rich environments can be demonstrated from the timescales of oxygen hop from the five different oxygen environments (S0 to S4). For this, we exploit the self-part of the van-Hove correlation function G_s(r̄, t), given by

where r̄_i(t) is the position vector of the ith particle at time t, N is the number of oxide ions, and δ is the Dirac delta function. G_s(r̄, t) thus describes the probability density that an ion initially (at t = 0) at its original site is found at a distance r (that is, in a spherical shell of volume 4πr²dr around the original site) after time t. Fig. 13 shows the corresponding probability, 4πr²G_s(r̄, t), that an ion is found in a shell of radius r and thickness unity, after an elapsed time t, in the five different environments for two different time intervals. As a note in passing, in a recent study, we have utilized this function to estimate the residence time of oxygens in yttria doped ceria (Y_xCe_1−xO_2−δ).⁷¹

Fig. 13 The evolution of the self-part of the van-Hove correlation function of the oxide ions in La_1−xSr_xFeO_3−x/2 from the five distinct oxygen environments (S0 to S4) at two different time intervals is demonstrated for the composition x = 0.6. These profiles are averaged over the NPT-MD trajectories at 1200 K.

Fig. 13 suggests that at short time intervals of 10 ps, only the first peak is prominent for all five oxygen environments, while for a larger time interval of 1 ns, the oxide ions start jumping off from their original site to the nearest sites (at a distance of ∼2.7–3.0 Å) leading to the development of second peaks (at the expense of the first). It can be seen that during the 1 ns, a significant fraction of the oxide ions jump off from the S0 environments (having no Sr in the neighborhood), as evident from the height of the second peak. However, as the number of Sr in the environment increases, this fraction decreases systematically, demonstrating that ion hops are becoming sluggish progressively.

4 Conclusions

A comprehensive investigation of various microscopic aspects related to oxide ion mobility in Sr doped LaFeO₃ (La_1−xSr_xFeO_3−x/2, x = 0 to x = 0.6) is carried out employing molecular dynamics simulations at 1200 K. The structural features of the rhombohedral phase examined in the present simulation are in agreement with experimental observations across the composition range. The present simulation predicts that the oxide ion conductivity of the La_1−xSr_xFeO_3−x/2 system drops with Sr concentration. However, the available experimental results on the ionic conductivities of these solids from different groups are rather scattered, both quantitatively as well as qualitatively, making it hard to make comparisons.

It is noted that the oxide ion diffusion in the system is through intra-octahedral hops, as opposed to inter-octahedral, along the edges of the FeO₆ octahedra. The microscopic energy barrier for the oxide ion diffusion along this pathway increases steadily with Sr doping corroborating the nature of the ionic conductivity observed. Further insights on the mechanism of ion transport are derived through analysis of the local oxygen environments and the bottlenecks along the migration channels. The individual energy of oxide ions decreases in the Sr-rich environments favoring higher occupancies of such sites, resulting in oxygen vacancies populating the La-rich environments. Thus, as a result of the random doping of Sr at the La sites, the oxygen vacancies remain positionally disordered in the La_1−xSr_xFeO_3−x/2 matrix.

Interestingly, we observe a cooperative mechanism wherein the triangular bottleneck formed by the cations opens up as the oxygen approaches it for an impending jump. The three kinds of bottlenecks, due to the different La/Sr combinations forming the corners of the triangle, are analyzed individually for their dimensions and energetic barriers. While the dimensions of the triangular bottlenecks are insensitive to the cations (La/Sr) forming its corners, the bottlenecks start posing substantially higher energy barriers as La is replaced by Sr. The larger ionic radius of Sr, as compared to La, results in narrowing of the effective window size, causing higher barriers for oxide ion migration. The resulting slowdown of ion hops in Sr-rich environments is explicitly demonstrated using time correlation functions. Thus doping results in the localization of oxide ions due to the increase in the fraction of Sr-rich environments in the matrix, which accounts for relatively higher oxygen occupancies and poses higher migration barriers. This explains the slower oxide ion transport in La_1−xSr_xFeO_3−x/2 systems with Sr concentration.

The study highlights that an increase in the effective bottleneck radius, which involves an intricate interplay between the lattice parameters and dopant size, could enhance the oxide ion transport. Thus, smaller dopants at A and B sites that do not cause significant contraction of the lattice need to be explored. Alternatively, the oxide ion transport could possibly be enhanced if the fraction of Sr-rich environments can be controlled through tailored synthesis strategies producing uniform distribution of dopants (Sr) in the matrix. The need for further experimental studies is also evident from the scatter of the conductivity data in the literature.

Author contributions

Not applicable.

Data availability

The data supporting this article have been included as part of the ESI.†

Conflicts of interest

No conflict of interest or competing interests to declare.

Acknowledgements

We thank the Centre for Development of Advanced Computing (CDAC), Pune, and Param Kamrupa, IIT Guwahati, managed by the National Super Computing Mission (NSM), India, for the generous allowance of CPU hours. The use of the LAMMPS⁵⁶ package for MD simulations and VMD⁷² and VESTA⁷³ software for visualization is acknowledged.

References

1. A. B. Stambouli and E. Traversa, Renewable Sustainable Energy Rev., 2002, 6, 433–455.
2. N. Q. Minh, J. Am. Ceram. Soc., 1993, 76, 563–588.
3. S. C. Singhal, Solid State Ionics, 2000, 135, 305–313.
4. R. M. Ormerod, Chem. Soc. Rev., 2003, 32, 17–28.
5. N. Mahato, A. Banerjee, A. Gupta, S. Omar and K. Balani, Prog. Mater. Sci., 2015, 72, 141–337.
6. P. Singh and N. Q. Minh, Int. J. Appl. Ceram. Technol., 2004, 1, 5–15.
7. D. J. Brett, A. Atkinson, N. P. Brandon and S. J. Skinner, Chem. Soc. Rev., 2008, 37, 1568–1578.
8. L. Malavasi, C. A. Fisher and M. S. Islam, Chem. Soc. Rev., 2010, 39, 4370–4387.
9. S. M. Haile, Acta Mater., 2003, 51, 5981–6000.
10. T. Elmer, M. Worall, S. Wu and S. B. Riffat, Renewable Sustainable Energy Rev., 2015, 42, 913–931.
11. J. Savioli and G. W. Watson, Curr. Opin. Electrochem., 2020, 21, 14–21.
12. S. B. Adler, Chem. Rev., 2004, 104, 4791–4844.
13. C. Sun, R. Hui and J. Roller, J. Solid State Electrochem., 2010, 14, 1125–1144.
14. A. B. Muñoz-Garca, M. Pavone, A. M. Ritzmann and E. A. Carter, Phys. Chem. Chem. Phys., 2013, 15, 6250–6259.
15. D. Ding, M. Liu, Z. Liu, X. Li, K. Blinn, X. Zhu and M. Liu, Adv. Energy Mater., 2013, 3, 1149–1154.
16. A. J. Jacobson, Chem. Mater., 2010, 22, 660–674.
17. E. D. Wachsman and K. T. Lee, Science, 2011, 334, 935–939.
18. L. Yang, S. Wang, K. Blinn, M. Liu, Z. Liu, Z. Cheng and M. Liu, Science, 2009, 326, 126–129.
19. J. C. Ruiz-Morales, D. Marrero-López, M. Gálvez-Sánchez, J. Canales-Vázquez, C. Savaniu and S. N. Savvin, Energy Environ. Sci., 2010, 3, 1670–1681.
20. Y. Zhu, W. Zhou, R. Ran, Y. Chen, Z. Shao and M. Liu, Nano Lett., 2016, 16, 512–518.
21. T. Ishihara, Perovskite oxide for solid oxide fuel cells, Springer Science & Business Media, 2009.
22. O. Tripathy and P. P. Kumar, J. Mater. Sci., 2017, 52, 6542–6553.
23. P. R. Chandran and T. Arjunan, Int. J. ChemTech Res., 2015, 7, 488–497.
24. T.-L. Wen, H. Tu, Z. Xu and O. Yamamoto, Solid State Ionics, 1999, 121, 25–30.
25. J. Richter, P. Holtappels, T. Graule, T. Nakamura and L. J. Gauckler, Monatsh. Chem., 2009, 140, 985–999.
26. A. Orera and P. Slater, Chem. Mater., 2010, 22, 675–690.
27. Y. Sakaki, Y. Takeda, A. Kato, N. Imanishi, O. Yamamoto, M. Hattori and Y. Esaki, Funtai oyobi Funmatsu Yakin, 1999, 46, 293–299.
28. K. Huang, H. Y. Lee and J. B. Goodenough, J. Electrochem. Soc., 1998, 145, 3220.
29. H. Kamata, A. Hosaka, J. Mizusaki and H. Tagawa, Solid State Ionics, 1998, 106, 237–245.
30. A. Chroneos, B. Yildiz, A. Tarancón, D. Parfitt and J. A. Kilner, Energy Environ. Sci., 2011, 4, 2774–2789.
31. R. De Souza and J. Kilner, Solid State Ionics, 1998, 106, 175–187.
32. P. G. Bruce, Solid state electrochemistry, Cambridge University Press, 1997.
33. A. B. Muñoz-Garca, A. M. Ritzmann, M. Pavone, J. A. Keith and E. A. Carter, Acc. Chem. Res., 2014, 47, 3340–3348.
34. M. Patrakeev, J. Bahteeva, E. Mitberg, I. Leonidov, V. Kozhevnikov and K. Poeppelmeier, J. Solid State Chem., 2003, 172, 219–231.
35. J. Mizusaki, T. Sasamoto, W. R. Cannon and H. K. Bowen, J. Am. Ceram. Soc., 1983, 66, 247–252.
36. T. Kawada and H. Yokokawa, Key Eng. Mater., 1996, 125, 187–248.
37. A. Chroneos, D. Parfitt, J. A. Kilner and R. W. Grimes, J. Mater. Chem., 2010, 20, 266–270.
38. E. Robens, R. Rauschen, J. Kaub, J. P. Parras, D. Kemp, C. L. Freeman and R. A. De Souza, J. Mater. Chem. A, 2022, 10, 2388–2397.
39. F.-F. Lu and H.-K. Tian, Phys. Chem. Chem. Phys., 2023, 25, 18973–18982.
40. M. Jaipal and A. Chatterjee, J. Phys. Chem. C, 2017, 121, 14534–14543.
41. P. P. Dholabhai, J. B. Adams, P. Crozier and R. Sharma, Phys. Chem. Chem. Phys., 2010, 12, 7904–7910.
42. M. S. Islam, D. J. Driscoll, C. A. Fisher and P. R. Slater, Chem. Mater., 2005, 17, 5085–5092.
43. S. K. Moharana and P. P. Kumar, Solid State Ionics, 2023, 403, 116409.
44. S. Madhual, K. Pramanik and P. P. Kumar, Phys. Chem. Chem. Phys., 2022, 24, 18281–18290.
45. P. P. Dholabhai, S. Anwar, J. B. Adams, P. Crozier and R. Sharma, J. Solid State Chem., 2011, 184, 811–817.
46. D. A. Andersson, S. I. Simak, N. V. Skorodumova, I. A. Abrikosov and B. Johansson, Proc. Natl. Acad. Sci. U. S. A., 2006, 103, 3518–3521.
47. A. Jones and M. S. Islam, J. Phys. Chem. C, 2008, 112, 4455–4462.
48. Y. A. Mastrikov, R. Merkle, E. A. Kotomin, M. M. Kuklja and J. Maier, Phys. Chem. Chem. Phys., 2013, 15, 911–918.
49. T. Das, J. D. Nicholas and Y. Qi, Phys. Chem. Chem. Phys., 2020, 22, 9723–9733.
50. A. M. Ritzmann, A. B. Munñoz-Garca, M. Pavone, J. A. Keith and E. A. Carter, Chem. Mater., 2013, 25, 3011–3019.
51. Y.-S. Zheng, M. Zhang, Q. Li, Y.-A. Zhu, Z.-J. Sui, D. Chen and X.-G. Zhou, J. Phys. Chem. C, 2018, 123, 275–290.
52. A. Bonkowski, J. A. Kilner and R. A. De Souza, RSC Appl. Interfaces, 2024, 1, 699–710.
53. A. Fossdal, M. Menon, I. Wærnhus, K. Wiik, M.-A. Einarsrud and T. Grande, J. Am. Ceram. Soc., 2004, 87, 1952–1958.
54. F. Schwarz, S. Barthel and A. Mace, Chem. Mater., 2024, 36, 11359–11376.
55. S. Dann, D. Currie, M. Weller, M. Thomas and A. Al-Rawwas, J. Solid State Chem., 1994, 109, 134–144.
56. A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown, P. S. Crozier, P. J. in't Veld, A. Kohlmeyer, S. G. Moore and T. D. Nguyen, et al. , Comput. Phys. Commun., 2022, 271, 108171.
57. F. H. Taylor, J. Buckeridge and C. R. A. Catlow, Chem. Mater., 2016, 28, 8210–8220.
58. P. P. Kumar and S. Yashonath, J. Chem. Sci., 2006, 118, 135–154.
59. H. Mehrer, Diffusion in solids: fundamentals, methods, materials, diffusion-controlled processes, Springer Science & Business Media, 2007, vol. 155.
60. C. Kittel and P. McEuen, Introduction to solid state physics, John Wiley & Sons, 2018.
61. H.-D. Wiemhöfer, H.-G. Bremes, U. Nigge and W. Zipprich, Solid State Ionics, 2002, 150, 63–77.
62. F. Bidrawn, S. Lee, J. M. Vohs and R. J. Gorte, J. Electrochem. Soc., 2008, 155, B660.
63. J. E. ten Elshof, H. J. Bouwmeester and H. Verweij, Solid State Ionics, 1995, 81, 97–109.
64. V. Vashook, M. Al Daroukh and H. Ullmann, Ionics, 2001, 7, 59–66.
65. T. Ishigaki, S. Yamauchi, K. Kishio, J. Mizusaki and K. Fueki, J. Solid State Chem., 1988, 73, 179–187.
66. M. Søgaard, P. V. Hendriksen and M. Mogensen, J. Solid State Chem., 2007, 180, 1489–1503.
67. T. Ishigaki, S. Yamauchi, J. Mizusaki, K. Fueki, H. Naito and T. Adachi, J. Solid State Chem., 1984, 55, 50–53.
68. S. Anand, B. Ouyang, T. Chen and G. Ceder, Phys. Rev. Mater., 2023, 7, 095801.
69. J. Kilner and R. Brook, Solid State Ionics, 1982, 6, 237–252.
70. K. Pramanik, K. Sau and P. P. Kumar, J. Phys. Chem. C, 2020, 124, 4001–4009.
71. S. Madhual and P. P. Kumar, J. Mater. Sci., 2023, 58, 4499–4512.
72. W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33–38.
73. K. Momma and F. Izumi, J. Appl. Crystallogr., 2011, 44, 1272–1276.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d5cp00131e

---