Oxide Ion Dynamics in Hexagonal Perovskite Mixed Conductor Ba 7 Nb 4 MoO 20 : A Comprehensive Ab Initio Molecular Dynamics Study

Hexagonal perovskite Ba 7 Nb 4 MoO 20 -related materials are very promising solid electrolytes with high oxide ion conductivity and redox stability, making them potentially applicable in solid oxide fuel cells. Optimizing the properties of this family of materials necessitates atomic-level understanding of the oxide ion dynamics leading to high conductivity. Here we report extensive ab initio molecular dynamics simulations of Ba 7 Nb 4 MoO 20 investigating oxide ion motions, which allowed the observation of a continuous diffusion pathway for oxide ions in the ( ab ) plane, but also revealed significant contribution of the oxygen atoms from crystallographic sites located outside this plane, to the long-range dynamics. To probe the timescale of oxide ion diffusion, complementary quasielastic neutron scattering experiments were carried out, and showed that oxide ion dynamics in Ba 7 Nb 4 MoO 20 , even at 950 o C, are too slow to be observable on a nanosecond timescale. Based on the atomic-level understanding of structure-property relationships afforded by this detailed computational study, we propose new materials design strategies with potential to significantly increase oxide ion conductivity in Ba 7 Nb 4 MoO 20 -related hexagonal perovskites, which target to simultaneously increase the number of oxide ion charge carriers and rotational flexibility of the (Nb/Mo)O x polyhedra.


Introduction
Materials with high oxide ion conductivity are of considerable research interest owing to their applications in oxygen sensors, oxygen-permeable membranes and as electrolytes in solid oxide fuel cells (SOFC), devices for green electricity generation.][3] This limitation can be overcome with new materials with excellent oxide ion conductivity in the temperature range between 300-500 o C.
The structure of Ba7Nb4MoO20 can be described as an ordered intergrowth of palmierite Ba3M2O8 units and 12R perovskite Ba4M3O12 units (where M = Nb or Mo), forming a 7H hexagonal perovskite (Figure 1). 26,31At room temperature (Figure 1a), there are four crystallographically independent metal (M) sites modelled as jointly occupied by Mo and Nb (due to very similar scattering power for both X-rays and neutrons not allowing refinement of the Mo/Nb occupancies).The intrinsic cation vacancies are located on M1 and M2 sites, with fractional occupancies at 50 o C of 0.899(3) and 0.101(3), respectively. 26With increasing temperature, the occupancy of the M1 site increases while that of M2 decreases, resulting in occupancies of 0.931(3) and 0.069(3) at 530 o C. 26 M3 and M4 sites, located in the 12R perovskite layers, as well as the four unique Ba sites are fully occupied at all temperatures. 26Three of the five unique oxygen sites are fully occupied at all temperatures -the O3 site located in the palmierite layer, as well as O4 and O5 sites in the 12R perovskite Ba4M3O12 units. 26y contrast, O1 and O2 sites, located in the (ab) plane of the palmierite Ba3M2O8 units are both partially occupied below 530 o C; above this temperature the O2 site is essentially unoccupied, 26 resulting in tetrahedral coordination environment for the M1 Nb/Mo site (Figure 1b).Possible pathways of oxide ion diffusion in Ba7Nb3.9Mo1.1O20.05,Ba7Nb3.85W0.15MoO20.0758][29] Similarly, bond-valence energy landscape calculations for Ba7Nb4MoO20 suggest the presence of diffusion pathways via empty O2 sites. 32In the same study, a brief (20 ps) ab initio molecular dynamics (AIMD) simulation, performed only at one temperature (1000 K), was used to investigate the proton conductivity mechanism and confirm the O1-O2-O1 inplane oxide ion dynamics. 32While this short simulation gave robust results with regard to the proton diffusion pathways, no long-range oxide ion diffusion could be simulated, with only a single O1-O2-O1 oxide ion jump observed.
Here we report an extensive AIMD computational study of oxide ion dynamics in Ba7Nb4MoO20.Performing simulations at three temperatures for up to 200 ps allowed the observation of a continuous diffusion pathway for oxide ions in the (ab) plane, but also revealed significant contribution of O3 oxygen atoms, located outside this plane, to the long-range dynamics.To complement this comprehensive computational study, neutron scattering experiments were used to probe the timescale of oxide ion diffusion in Ba7Nb4MoO20.Based on the atomic-level understanding of structure-property relationships afforded by this study, we propose new doping strategies with potential to result in a significant increase in oxide ion conductivity of Ba7Nb4MoO20-related hexagonal perovskite materials.

Experimental
Synthesis.For neutron scattering experiments, two large (7 g) polycrystalline samples of Ba7Nb4MoO20 were prepared using stoichiometric amounts of BaCO3 (Merck 99 %), Nb2O5 (Fisher Sci 99.5 %) and MoO3 (Fisher Sci 99.5 %). 27The reagents were ball-milled under ethanol for 30 min at 300 rpm.The samples were left to dry and then heated at 900 o C for 12 h.The powder was ground thoroughly using an agate mortar and pestle, pressed into a rod using a hydraulic press and approximately 70 MPa, and heated at 1100 o C for 24 h.The progress of the solid state reaction and product purity were monitored using powder X-ray diffraction (PXRD), performed on a Bruker D8 ADVANCE diffractometer using Cu Kα radiation.The collected patterns were analyzed using the Rietveld method 33 implemented in Topas Academic, 34 using the structural data obtained from neutron diffraction data for Ba7Nb4MoO20 reported by Fop et al. 26 as a starting model.The refinement parameters included cell parameters, zero-point, one global temperature factor, and pseudo-Voigt peak shape parameters.The background was modelled using a 9 parameter Chebyshev polynomial.The final product was found to be 92(1) % pure, with Ba5Nb4O15 (2.3(3) %) and Ba3NbMoO8.5 (5(1) %) identified as the impurities present.The final Rietveld fit obtained is shown in Figure S1.

Quasielastic Neutron Scattering.
Neutron scattering data were collected on the backscattering spectrometer IN16b at the Institut Laue-Langevin (ILL) in Grenoble, 35 using an incident neutron wavelength of 6.271 Å and instrument resolution of ∼ 0.8 µeV (FWHM).This gives access to dynamics on a nanosecond timescale.The Q-range accessible on IN16b is 0.1 -1.8 Å⁻ 1 , but due to the weak signal the Q-integrated data were used.Measurements were performed using a furnace, 36  protons, we performed measurements in vacuum, probing oxide ion diffusion, and in humid N2, probing proton diffusion.The energy transfer window used for all QENS measurements was ± 7 μeV.
In both atmospheres, elastic (E = 0 μeV) and inelastic (E = 2 μeV) scattered intensities were measured in steps of 25 K, so called elastic and inelastic fixed window scans (EFWS, IFWS). 37r the measurements performed under vacuum, 5.3 g of sample were placed in a Nb sample holder (diameter ∼ 8 mm, sample height ∼ 4 cm).During heating to 950 °C, EFWS (2 min) and IFWS (10 min) were collected, and at 950 °C a QENS measurement was performed with a total data collection time of 5 h.The sample was then cooled to 50 °C while again collecting EFWS (1 min) and IFWS (5 min).At 50 °C, QENS data were collected for 3 h.
For the measurement in humid atmosphere, approximately 8 g of sample were placed in a stainlesssteel sample holder (diameter ∼ 10 mm).N2 gas was bubbled through distilled water and blown on the sample via a thin metal capillary in the steel tube.A first QENS spectrum (5 h) was measured at room temperature (RT).The sample was then heated to 500 °C during which EFWS (2 min) and IFWS (10 min) were collected.During heating, QENS measurements were also performed at 350 °C, 425 °C and 500 °C for 6 h per temperature.The sample was then cooled to 225 °C, measuring EFWS (2 min) and IFWS (10 min), to collect QENS data at 225 °C for 6 h, and after further cooling at 175 °C for 4 h.Data were also collected at RT on the empty tube (QENS: 4 h; EFWS: 2 min; IFWS 10 min) and on a vanadium standard (QENS: 30 min) for data reduction and analysis, which was carried out using the Mantid 38 software.
Computational Methods.Density Functional Theory (DFT) simulations were performed with the Vienna ab initio simulation package (VASP). 39The projector augmented wave (PAW) method 40 with the GGA-PBE (Perdew-Burke-Ernzerhof generalized-gradient-approximation) exchange-correlation functional 41 were used for all simulations.Before any ab initio molecular dynamics (AIMD) simulations were performed, a geometry optimization of this fixed supercell was carried out, allowing ions to relax until the residual forces on them were less than 0.01 eV Å⁻ 1 .A structural model in space group P1 and containing 288 atoms, consisting of a 3×3×1 supercell of the high-temperature (530 °C) structure of Ba7Nb4MoO20 determined from neutron diffraction 26 was used as a starting point for the simulations.As the occupancy of the M2 site at this temperature is very low (0.069), it was assumed to be empty, resulting in full occupancy of the M1 site.Distribution of Mo atoms as well as the selection of one out of three partially occupied symmetry-equivalent O1 positions for each M1O4 tetrahedron was performed using a random number generator.The atomic content of the final simulation box was Ba63Nb36Mo9O180.AIMD simulations were carried out with a cut-off energy of 300 eV and the k-space was sampled at the gamma-point only.AIMD simulations of 180 ps (1000 o C), 100 ps (1500 o C), and 200 ps (2000 o C) using the NVT ensemble and a time step of 2 fs were performed.The data analysis was carried out using the MDANSE code, 42 and the space visited by oxide ions was visualized in cloud plots using LAMP. 43

Results and Discussion
Timescale of Ionic Diffusion.To complement the extensive molecular dynamics study which is the main focus of this paper, quasielastic neutron scattering experiments were carried out.This section offers a brief discussion of the elastic (EFWS) and inelastic (IFWS) fixed window scans (Figure S2), https://doi.org/10.26434/chemrxiv-2023-2mf82ORCID: https://orcid.org/0000-0002-2151-169XContent not peer-reviewed by ChemRxiv.License: CC BY-NC-ND 4.0 collected on IN16b, aiming to provide an initial overview of the timescale of oxide ion and proton diffusion in Ba7Nb4MoO20.Generally, a gradual decrease of elastic intensity on heating is expected due to the Debye-Waller factor, which describes the loss of elastic intensity caused by an increase of thermal motions of atoms at higher temperatures.Any decrease of elastic intensity in addition to that is indicative of other effects, for example diffusion or dehydration.If the elastic intensity decreases due to diffusion processes, the IFWS is expected to increase at the same temperature.Dehydration, on the other hand, results in a loss of intensity in EFWS with no increase in inelastic intensity at the same temperature, due to the decrease of the total number of scattering atoms in the system; because hydrogen is a particularly strong neutron scatterer, this effect is very pronounced.
A large (∼ 50 %) decrease in elastic intensity occurs between 210 o C and 360 o C on heating under vacuum (Figure S2a).There is no increase of inelastic intensity in the same temperature range (Figure S2b), which suggests that this is caused by dehydration, i.e the loss of protons from the sample, and not an onset of dynamics.This observation is in excellent agreement with TGA data, 26 showing the dehydration of the sample between approximately 200 and 400 o C (Figure S3).However, the IFWS in both atmospheres do not show an onset of significant dynamics in Ba7Nb4MoO20, and no broadening was observed in the QENS data collected at 950 o C under vacuum (Figure S4), or 225 o C in humid N2 (Figure S5).This suggests that ionic diffusion in Ba7Nb4MoO20 is slower than the nanosecond timescale probed here.
Although QENS can be a powerful tool for the study of oxide ion diffusion, [44][45][46][47][48] the small and predominantly coherent neutron scattering cross section of oxygen makes both observation and analysis of QENS difficult in many cases.Ab initio molecular dynamic (AIMD) simulations are very compatible with QENS, 10,11,45,47 ⁻ 51 and well-suited for the study of atomistic diffusion.Consequently, computational results are the central emphasis of the investigation of oxide ion diffusion in Ba7Nb4MoO20 presented here.
Oxide Ion Dynamics Investigated by Ab Initio Molecular Dynamics.Mean square displacement (MSD) plots for all oxygen atoms at 1000, 1500, and 2000 °C are shown in Figure 2a.The flat MSD plot at 1000 °C is indicative of the absence of significant oxide ion dynamics, showing that higher temperatures are required to simulate diffusion of oxide ions.At 1500 °C, a slope in the MSD curve is apparent, suggesting the onset of long-range diffusion.The calculated diffusion coefficient of 0.9 × 10⁻ 7 cm 2 s⁻ 1 indicates that diffusion is slow.This is consistent with the absence of the QENS signal in the neutron scattering data collected on IN16b at 950 o C in Ba7Nb4MoO20.For comparison, the diffusion coefficient determined from AIMD for the excellent oxide ion conductor Bi0.931V0.087O1.587, on which dynamics were clearly observed on IN16b, 46 was approximately 5 × 10⁻ 7 cm 2 s⁻ 1 at 600 °C (i.e.five times higher than that for Ba7Nb4MoO20 at 1500 °C). 46Further increase of the simulated temperature to 2000 o C leads to a more pronounced slope in the MSD curve caused by significant dynamics, with the diffusion coefficient increasing to 4 × 10⁻ 7 cm 2 s⁻ 1 .While this simulation temperature is very high, AIMD calculations often require higher temperatures than experimental to observe long-range diffusion.This is due to the intrinsic constraints of this computationally expensive method, which significantly limit both the size of the simulation cell as well as the achievable simulation time.By accelerating the diffusion process, the increased temperature reduces computational efforts while ensuring that an adequate number of jumps are observed.Hence, it allows a more efficient analysis of the long-range dynamics, which are often slow compared to the timescale accessible by AIMD.53][54][55][56][57][58][59][60][61][62][63][64] This is also evident from the previous AIMD simulation on Ba7Nb4MoO20 at 1000 K, 32 where only one oxide ion jump was observed, despite the high conductivity measured experimentally at this temperature.In the present work the simulation at 2000 o C was analyzed in detail to gain atomic-level insights into the oxide ion dynamics mechanisms.The individual MSD curves for the oxide ions from different crystallographic sites occupied in the high-temperature structure (O1, O3, O4, O5) indicate that only those originally at O1 and O3 sites contribute to the dynamics in Ba7Nb4MoO20, while O4 and O5 oxygen atoms do not (Figure 2b).Cloud plots can be used to visualize the simulated trajectories by highlighting the space visited by oxygen atoms over the course of the simulation.Figure 3a shows the trajectories of all oxide ions originating on the O1 site in the simulation at 2000 o C, demonstrating a continuous pathway for oxide ion diffusion in the (ab) plane.8][29] However, our AIMD simulations reveal for the first time a significant contribution to the overall diffusion by the O3 oxygen atoms, located outside the plane of O1-O2-O1 hopping, as shown in the MSD curve and the cloud plots (Figures 2b and 3b).The O3 oxygen atoms jump in the direction essentially perpendicular to the (ab) plane, to the nearest O1 site and then diffuse via the O1-O2-O1 2D pathway (Figure 3b).6][67][68][69][70] Comparing the number of jumps shows that O3 dynamics become more https://doi.org/10.26434/chemrxiv-2023-2mf82ORCID: https://orcid.org/0000-0002-2151-169XContent not peer-reviewed by ChemRxiv.License: CC BY-NC-ND 4.0 significant with increasing temperatures (Table S1), but are prevalent (accounting for at least twothirds of all jumps) at both 1500 o C and 2000 o C, which can be explained by the activation energy for localized motions being generally lower than for long-range diffusion. 46,71

Conclusions and Implications on Improvement of Oxide Ion Conductivity of Ba7Nb4MoO20
We have investigated oxide ion dynamics in Ba7Nb4MoO20 using ab initio molecular dynamics and complementary neutron scattering.Our quasielastic neutron scattering experiments showed that oxide ion dynamics are not on the nanosecond timescale probed, indicating that long-range diffusion is slower and higher temperatures would be necessary to directly observe it.This finding is supported by the AIMD results, which do not show any oxide ion diffusion at 1000 o C on a timescale of hundred ps.
Modelling of possible ionic migration pathways reported so far [27][28][29] has revealed ionic hopping between the partially occupied O1 and O2 sites, suggesting that oxide ion conductivity could be improved by aliovalent doping of Ba7Nb4MoO20 to introduce extra oxygen in the structure and increase the number of charge carriers.The approach has been demonstrated in Ba7Nb3.9Mo1.1O20.05,Ba7Nb3.85W0.15MoO20.075,8][29][30] Our extensive AIMD simulations reveal for the first time a significant contribution of the out-of-plane O3 dynamics to the conductivity mechanism, through a process which essentially amounts to the rotation of the (Nb/Mo)Ox polyhedra followed by injection of O3 atoms into the (ab) plane.The contribution of this process is very significant, as evidenced by the fact that these jumps outnumber those in the (ab) plane.Based on this, we propose strategies for new compositions related to Ba7Nb4MoO20, to facilitate and enhance (Nb/Mo)Ox polyhedral rotations, as well as introduce extra charge carriers.For example, a systematic exploratory synthesis strategy, targeting improvements of both in-plane and out-of-plane oxide ion dynamics, could be to substitute some Mo 6+ with V 5+ (thereby introducing rotationally flexible VOx) 46,51,71,72 and chargecompensate by replacing Ba 2+ with La 3+ to create the Ba7-xLaxNb4Mo1-yVyO20 series; the best performers from this series could then be modified further to introduce extra oxide ion carriers, by varying the Nb 5+ :Mo 6+ ratio.This two-pronged approach is likely to result in more significant improvements of oxide ion conductivity relative to Ba7Nb4MoO20.

Figure 1 .
Figure 1.Crystal structure of Ba7Nb4MoO20 obtained from neutron diffraction data: (a) at 50 o C and (b) at 530 o C. 26

Figure 2 .
Figure 2. (a): MSD curves obtained from AIMD of all oxygen atoms at three simulation temperatures; (b): contribution of the individual oxygen sites to the MSD in the simulation performed at 2000 o C.

Figure 3 .
Figure 3. Visualization of the 200 ps trajectory obtained from AIMD simulations at 2000 o C. (a): The red arrows highlight the continuous pathway observed in the (ab) plane caused by O1 jumps via O2 sites and grey clouds show the space visited by all O1 atoms; (b): Cloud plots showing significant delocalization and dynamics of O3 oxygen atoms.The black arrows highlight the jump of an O3 oxygen atom to an adjacent Nb/Mox via two O1 sites.Ba atoms have been removed for clarity.
Diffraction patterns of the synthesized compounds from PXRD.EFWS and IFWS in vacuum and in humid N2, Comparison of EFWS with TGA.QENS data of the dry sample.QENS data of the humid sample.MSD of individual oxygen sites at 1000 o C, MSD of individual oxygen sites at 1500 o C.
Table summarizing jumps during simulation.