Lithium superionic behavior and defect robustness in LiNbOCl₄: a first-principles molecular dynamics study

Abstract

Lithium-ion conductors that exhibit high ionic conductivity, thermal robustness, and mechanical compliance are essential for advancing all-solid-state battery technologies. In this study, we systematically investigate Li-ion transport in pristine and LiCl Schottky-defected LiNbOCl₄ using density functional theory molecular dynamics (AIMD). Pristine LiNbOCl₄ demonstrates robust superionic behavior with low activation energy (0.236 eV) and high room-temperature conductivity (9.57 × 10⁻³ S cm⁻¹), facilitated by a rigid Nb–O–Cl framework and a disordered Li sublattice. Introducing LiCl Schottky defects slightly increases the activation energy to 0.241 eV and slightly reduces conductivity to 8.20 × 10⁻³ S cm⁻¹. While defects preserve global percolation networks and mechanical softness, they introduce localized structural disruptions at vacancy-adjacent polyhedra. Notably, the dynamic gating mechanism where coherent anion rotations transiently expand diffusion bottlenecks is impaired. Unlike the classical paddle-wheel mechanism involving rotating polyanion clusters, this mechanism describes a distinct mode of transient bottleneck expansion driven by coordinated motion of individual halide and oxide anions within an oxyhalide lattice. The disruption of this mechanism is reflected in the emergence of rotational incoherence, rapid bond angle decorrelation, attenuation of high-frequency O-based phonon modes (∼80 meV), and a strong reduction and spatial localization of anion reorientation events that are associated with enhanced Li-ion motion, as identified by event-triggered ensemble. Together, these effects suppress bottleneck breathing, limiting the transient widening of diffusion pathways and effectively increasing the migration barrier. Li space–time correlation analysis further reveals diminished temporal coherence and transport cooperativity in the defected structure. These findings underscore the importance of cooperative lattice dynamics in enabling low-barrier Li-ion transport and demonstrate that LiNbOCl₄ retains high conductivity even under defect-induced perturbations, establishing it as a defect-tolerant candidate for next-generation solid electrolytes.

---

1. Introduction

All-solid-state batteries (ASSBs) are attracting considerable attention as a promising alternative to conventional lithium-ion batteries (LIBs), primarily due to their enhanced safety, higher energy density, and suitability for advanced applications such as electric vehicles and portable electronics.¹ In contrast to traditional LIBs, which rely on flammable liquid electrolytes, ASSBs employ solid electrolytes (SEs) to facilitate lithium-ion (Li-ion) transport between the anode and cathode. The use of SEs eliminates the risks associated with liquid electrolyte leakage and flammability, thereby significantly improving the thermal and mechanical stability of batteries for large-scale applications.^1,2

An ideal SE must possess high Li-ion conductivity, chemical stability, and interfacial compatibility with both the cathode and anode materials. Generally, SEs can be classified into three main categories: oxides, sulfides, and halides, each offering unique advantages and challenges. Oxide-based SEs are recognized for their excellent chemical and electrochemical stability but often exhibit relatively low Li-ion conductivity.^3,4 Sulfide-based SEs, conversely, demonstrate high ionic conductivity but suffer from poor chemical and electrochemical stability, rendering them vulnerable to degradation under operating conditions.^5–7 Halide-based SEs offer a balance between these extremes, combining moderate to high ionic conductivity with good chemical and oxidative stability.⁸ However, their relatively low reductive stability may limit their applicability in certain battery architectures.⁹

Lithium (Li) oxyhalides have emerged as a new class of SEs that potentially integrate the favorable properties of both oxides and halides. By introducing oxygen into the halide lattice, oxyhalides are expected to retain high ionic conductivity akin to halides while enhancing chemical and electrochemical stability similar to oxides.^10–13 In particular, oxyhalide-type LiNbOCl₄ has been the focus of recent studies; it exhibits an ionic conductivity of about 10 mS cm⁻¹ at room temperature.^10,13 Structurally, LiNbOCl₄ has been described as a framework of corner-O-sharing NbO₂Cl₄ octahedra forming parallel chains, sometimes written as (NbCl₄O₂/2)∞, where O₂/2 indicates two bridging O atoms each shared by two Nb centers (i.e., ½O per Nb).^10,13 These chains create loosely packed regions that host disordered Li-ions; the Li sublattice does not form distinct layers but instead occupies partially filled interstitial sites.¹⁰ This structural motif differs from classical layered compounds like LiCoO₂, where well-defined alternating CoO₆ and LiO₆ octahedral layers exist.¹⁴ Instead, in LiNbOCl₄, the layered nature arises from the anisotropic connectivity of the rigid Nb-based octahedral network, while Li mobility is facilitated within and between these slabs due to structural openness and ionic disorder.¹⁰ In addition to its promising ionic conductivity, LiNbOCl₄ exhibits a low elastic stiffness which is critical for maintaining good interfacial contact and thus, a good structural integrity during battery operation.^10,13 These combined attributes place LiNbOCl₄ in direct competition with established SEs such as garnet-type Li₇La₃Zr₂O₁₂ (LLZO)¹⁵ and sulfide-based materials.^16–18 To date, a comprehensive understanding of the Li-ion migration mechanism within the LiNbOCl₄ structure remains under active investigation. Computational approaches, particularly those based on density functional theory (DFT), have proven invaluable for probing ion dynamics, defect chemistry, and thermodynamic stability in SE materials.^10,13,19,20 By integrating theoretical insights with experimental observations, it becomes possible to fine-tune material properties and guide the discovery of new electrolyte compositions with enhanced performance.

In this study, we systematically explore the fundamental properties of LiNbOCl₄ SE using DFT calculations. Combined with molecular dynamics (MD) simulations, we investigate Li-ion transport in both pristine and LiCl Schottky-defected structures. The initial part of this work focuses on establishing a computational framework by benchmarking various van der Waals (vdW) correction schemes within DFT to accurately model LiNbOCl₄'s structural and Li-ion transport properties. Structural and thermodynamic stabilities are assessed via convex hull analysis. We further analyze Li-ion dynamics through temperature-dependent mean square displacement (MSD) analysis, extract diffusion coefficients, and quantify ionic conductivity using the Nernst–Einstein relation. Activation energies for Li migration are derived from Arrhenius behavior, allowing a direct comparison between pristine and defected systems, as well as available experimental data. Our results provide useful insights for the design of high-conductivity and practical SEs for ASSBs.

2. Computational methods

2.1 Initial structure and atomic model generation

The initial structural model for LiNbOCl₄ is based on experimental crystallographic data reported for the tetragonal phase (space group I4/m).¹⁰ The lattice parameters are set to a = b = 8.9109 Å and c = 3.9542 Å, with angles α = β = γ = 90°. The atomic site coordinates and occupancies are summarized in Table S1 (in the SI). Li atoms occupy three partially filled 8h Wyckoff positions (Li(1), Li(2), and Li(3)) with partial occupancies of 0.04, 0.17, and 0.04, respectively, reflecting experimentally observed Li disorder.

To capture the intrinsic disorder brought upon by the Li partial occupancy sites, a 2 × 2 × 3 supercell is constructed and 1000 distinct Li configurations are generated by random sampling of Li sites. This procedure preserves the experimentally observed site occupancy ratios. The atomic structure is illustrated in Fig. 1, where Cl, Nb, O, and Li atoms are depicted with distinct colors to highlight the spatial distribution of partially occupied Li sites.

Fig. 1 Atomic structure representation of LiNbOCl₄, highlighting atomic configurations and Li occupancy distributions across crystallographic sites. Chlorine (Cl) atoms are depicted in green, niobium (Nb) atoms in blue, oxygen (O) atoms in red, and lithium (Li) atoms with partial occupancy in magenta. This visualization illustrates the atomic arrangement within the LiNbOCl₄ crystal lattice, with Li disorder emphasized according to occupancy variations, guided by symmetry considerations.

Ewald energy calculations are performed on each configuration using the pymatgen library²¹ to evaluate electrostatic stability. Based on this Ewald energy screening step, the 10 configurations with the lowest electrostatic energies are selected for subsequent full DFT relaxation. This approach ensures a representative sampling of low-energy Li distributions consistent with experimental findings.

To introduce the Schottky-type LiCl defect, one Li and one Cl atom are removed from each 168-atom LiNbOCl₄ supercell (Fig. 2), corresponding to a vacancy concentration of ∼1.1 at%. To systematically assess the defect energetics, 1000 distinct Li–Cl vacancy arrangements are generated in the same relaxed supercell. For each configuration, (i) the Ewald electrostatic energy and (ii) the minimum periodic separation between the Li and Cl vacancies, d_VLi–VCl are evaluated. When the configurations are ordered based on ascending Ewald energy, short vacancy separations (d_VLi–VCl ≲ 3 Å) populate the low-energy region, whereas structures with d_VLi–VCl ≳ 5 Å occur predominantly at higher Ewald energies (Fig. S2). Because the Ewald sum isolates the long-range coulombic contribution at a fixed composition and volume, this trend indicates that compact Li–Cl vacancy pairs are electrostatically preferred. The lowest-energy, short-separation configuration is therefore selected as the representative LiCl-Schottky model and it is relaxed fully using the same convergence criteria as for the pristine supercell, enabling a direct comparison of the relaxed geometries. The resulting lattice response is modest: after full relaxation; the a lattice parameter increases from 18.8488 Å (pristine) to 18.9803 Å (+0.70%), and b from 18.8935 Å to 19.1622 Å (+1.42%), while c remains essentially unchanged at ≈12.04 Å (−0.002%). The off-diagonal elements of the cell matrix are ≤∼0.02 Å, indicating negligible shear distortions. These values, taken directly from the relaxed cell matrices, show that the LiCl Schottky defect induces only a small, nearly isotropic in-plane expansion and leaves the c axis essentially intact. This relaxed lowest-energy configuration is used in all subsequent AIMD simulations.

Fig. 2 Structural representation of the Schottky-defected LiNbOCl₄ supercell highlighting the introduction of Li (V_Li) and chlorine (V_Cl) vacancies. Left: three-dimensional perspective showing the global arrangement of atoms and the position of vacancy pairs. Right: top-down view along the c-axis, emphasizing the local structural distortions induced by vacancies. Nb atoms are shown in blue, Cl in green, O in red, and Li in magenta. Vacancies are indicated by red circles and arrows, illustrating how localized removal of Li and Cl atoms perturbs the surrounding lattice and potentially affects ion transport pathways.

2.2 DFT calculation conditions

To ensure an accurate and reliable structural and transport property analysis of LiNbOCl₄, a systematic benchmarking of vdW correction schemes is carried out within the framework of DFT. All calculations are performed using the Vienna Ab initio Simulation Package (VASP)²² with the projector augmented-wave method.²³ Exchange–correlation effects are treated using the generalized gradient approximation as parametrized by Perdew, Burke, and Ernzerhof (PBE).²⁴ Various vdW corrections are evaluated, including DFT-D2,²⁵ DFT-D3,²⁶ DFT-D3(BJ),²⁷ optB88,²⁸ optPBE,²⁸ optB86b,²⁸ DF,²⁹ and DF2,³⁰ to assess their performance in reproducing experimental structural parameters.

All structural optimizations are conducted on a 2 × 2 × 3 supercell (168 atoms), with a plane-wave cutoff energy of 520 eV and a k-point grid spacing of 0.3–0.5 Å⁻¹. Electronic energy convergence is set to 1 × 10⁻⁵ eV, and ionic force convergence to 0.01 eV Å⁻¹, ensuring precise and stable relaxation.

2.3 Thermodynamic stability evaluation

The thermodynamic stability of pristine and Schottky-defected LiNbOCl₄ systems is evaluated through convex hull analysis based on formation energy calculations.³¹ Formation energies are derived from fully geometry-optimized structures obtained via DFT calculations and post-processed with the Pymatgen library.²¹

The energy above the convex hull (E_d) at 0 K is used as a measure of thermodynamic (meta)stability; it is computed for each structure using eqn (1):

E_d = E_phase − E_hull (1)

where E_phase is the formation energy per atom of the phase of interest and E_hull is the convex hull energy per atom at the same composition. The convex hull is constructed by identifying the set of phases with the lowest energies at each composition. Phases with E_d = 0 lie on the convex hull and are thermodynamically stable at 0 K, whereas non-ground-state phases with E_d > 0 are located above the hull.I n this case, E_d represents the decomposition energy, i.e., the energy difference between the phase and the closest linear combination of stable phases at the same composition.

Physically, E_d quantifies the thermodynamic driving force for decomposition: a zero E_d means the phase is stable, while a small positive E_d indicates a metastable phase that is still potentially synthesizable if the kinetic barriers to decomposition are high. Full derivation and implementation details for the convex hull construction are provided in the SI (Section S3).

2.4 Mechanical property evaluation via pressure–volume fitting

The mechanical compressibility of LiNbOCl₄ is evaluated through direct pressure–volume (P–V) fitting based on DFT-calculated structures.³² External hydrostatic pressures of 0, 2, 4, 6, 8, and 10 kbar are applied to the fully optimized supercell. At each pressure step, full structural relaxation is performed, allowing simultaneous optimization of both atomic positions and lattice vectors. This procedure enables a direct simulation of the material's volumetric response under hydrostatic stress.

The bulk modulus (B) is determined from the linear relationship between pressure and volume according to eqn (2):

where V is the equilibrium volume at zero pressure and dP/dV is the slope obtained from linear fitting of the P–V data. This direct fitting method provides an alternative to the commonly used Birch–Murnaghan equation of state³³ and is particularly well-suited for mechanically soft ionic materials such as LiNbOCl₄.

2.5 Thermal expansion analysis

The thermal expansion behavior of LiNbOCl₄ is investigated using constant-pressure AIMD simulations performed with the VASP in the isothermal–isobaric (NPT) ensemble with a thermostat damping parameter (γ) of 10 ps⁻¹. Temperature control is achieved using the Langevin thermostat,³⁴ with simulations at fixed temperatures of 600 K, 800 K, 1000 K, 1200 K, and 1400 K. Each simulation is run for 10 ps with a time step of 1 fs. The final 5000 steps of each trajectory are used to extract the average cell volume corresponding to each temperature.

The volumetric thermal expansion coefficient (α_V) is determined by linearly fitting the temperature dependence of the equilibrium volumes. It is calculated according to eqn (3):

where V₀ is the average volume at 600 K and dV/dT is the slope of the linear fit. This approach provides a direct assessment of the lattice's thermal response under conditions relevant to battery operation.

2.6 Li-ion transport analysis

Following thermal expansion analysis, AIMD simulations are performed to investigate Li-ion transport properties. Initial equilibration is conducted in the NPT ensemble using the Langevin thermostat for 10 ps to determine equilibrium volumes at each temperature. Fixed-volume supercells derived from this equilibration step are then used for production simulations in the canonical (NVT) ensemble, employing the Nosé–Hoover thermostat^35,36 with an integration time step of 1 fs. Each production run spanned 100 ps (100 000 steps) at target temperatures of 600, 800, 1000, 1200, and 1400 K.

Li-ion diffusivity is determined from the mean squared displacement (MSD) curves. The diffusion coefficients (D_Li) are calculated from the linear regime of the MSD using the Einstein relation, expressed as eqn (4):

where 〈|r_i(t) − r_i(0)|²〉 denotes the ensemble-averaged MSD of Li ions.

The Li-ion conductivity (σ_Li) is estimated using the Nernst–Einstein relation (eqn (5)):

where n is the number density of Li ions, q is the elementary charge, k_B is the Boltzmann constant, and T is the simulation temperature.

Activation energies (E_a) for Li diffusion are extracted by performing Arrhenius fits to the temperature dependence of the diffusion coefficients and ionic conductivities according to eqn (6):

where D₀ is the pre-exponential factor. Linear regressions of D_Li and ln σ_Li versus 1/T were used to obtain E_a from the slope. Trajectory density maps and structural snapshots were further analyzed to visualize Li-ion percolation networks and local diffusion pathways across different temperatures.

3. Results and discussion

3.1 Validation of structural models via van der Waals benchmarking

The benchmarking of dispersion-corrected density functionals is performed, and the results are summarized in Table 1. All values in Table 1 are obtained from geometry optimizations carried out in this work using identical initial structures and computational parameters; only the dispersion (vdW) treatment is varied across functionals. The performance of each functional is evaluated based on the deviations of lattice parameters, unit cell volume, and Nb–Cl bond lengths relative to experimental values. In the absence of vdW correction, the structure exhibits significant overestimations in both lattice constants and bond lengths, with a total deviation exceeding 6.8%. Among the tested functionals, DFT-D3(BJ), DFT-D2, and optB88 demonstrate the best agreement with experimental data. Notably, DFT-D3(BJ) produces the lowest total deviation (1.35%), with the a and c lattice parameters deviating by only 0.92% and 0.24%, respectively, and the unit cell volume within 0.04% of the reference value. The average Nb–Cl bond length is also consistent with experimental observations for related compounds, such as TlNbOCl₄ (2.38–2.41 Å).^10,13

Table 1 Comparison of vdW functional performance in predicting lattice parameters and Nb–Cl bond length. This table summarizes the calculated lattice parameters (a and c), unit cell volume, and average Nb–Cl bond length for various vdW functionals applied to the LiNbOCl₄ system. Percentage differences relative to the reference values¹³ are provided for the bond length, lattice parameters, and volume to assess the accuracy of each functional. The average is the combined percentage difference for all metrics, providing an overall measure of deviation

vdW functional	Lattice parameter		Volume (Å^3)	Average Nb–Cl bond length (Å)	Percentage difference to ref. 13 (%)
	a (Å)	c (Å)			Nb–Cl bond length	a	c	Volume	Average
Experiment^7	8.911	3.954	313.98	2.32	0.0	0.0	0.0	0.0	0.0
Without vdW (this work)	9.471	4.012	360.994	2.426	4.573	6.291	1.467	14.974	6.826
Df (this work)	9.227	4.042	345.391	2.454	5.784	3.543	2.227	10.004	5.39
df2 (this work)	9.042	4.02	332.811	2.472	6.56	1.47	1.659	5.997	3.922
dftd2 (this work)	8.881	3.965	317.641	2.418	4.203	0.34	0.263	1.166	1.493
dftd3 (this work)	8.903	3.978	319.747	2.422	4.401	0.084	0.595	1.837	1.729
dftd3-bj (this work)	8.829	3.964	313.85	2.418	4.216	0.918	0.243	0.041	1.354
optb88 (this work)	8.831	3.952	312.708	2.428	4.634	0.056	0.405	0.893	1.497
Optpbe (this work)	9.006	3.987	325.767	2.436	5.009	1.068	0.836	3.754	2.666
optb86b (this work)	8.803	3.947	311.183	2.419	4.263	1.213	0.178	0.891	1.636
r^2SCAN-D4 (this work)	9.521	4.013	362.859	2.426	4.569	6.844	1.483	15.568	7.116

---

These results confirm that explicit inclusion of vdW interactions is critical for accurately reproducing the lattice parameters and Nb–Cl bond lengths vs. the experiment. Based on this benchmarking, the DFT-D3(BJ) functional was selected for all subsequent static and dynamic simulations, including structural relaxations, AIMD simulations, and defect calculations.

3.2 Thermodynamic stability

The thermodynamic stability of pristine and Schottky-defected LiNbOCl₄ structures is evaluated through E_hull analysis, based on formation energies calculated from fully relaxed DFT geometries. The pristine LiNbOCl₄ compound exhibits a formation energy of −2.101 eV per atom and lies 0.0211 eV per atom above the convex hull, indicating metastability. The Schottky-defected structure, generated by removing one Li and one Cl atom per supercell (Li₂₃Nb₂₄Cl₉₅O₂₄), shows a marginally more-negative formation energy of −2.104 eV per atom and a reduced hull distance of 0.0185 eV per atom.

Both E_hull values fall within the commonly accepted metastability threshold (<50 meV per atom), suggesting that pristine and defected LiNbOCl₄ are potentially synthesizable under appropriate experimental conditions. Decomposition pathway analysis identified NbCl₃O and LiCl as the dominant competing phases, with relative phase fractions of ∼71% and ∼29% for the pristine structure, and 72% and 28% for the defected structure, respectively.

The slight reduction in E_hull with LiCl-Schottky defects indicates that the introduction of Li–Cl vacancy pairs does not significantly destabilize the host lattice and may even promote a degree of stabilization. This result supports the synthesizability of Schottky-defected LiNbOCl₄ which is also consistent with the dynamic stability observed during AIMD simulations (see the later section).

3.3 Mechanical property evaluation

Fig. S6 (see SI S6) presents the pressure–volume (P–V) relationships for pristine and LiCl Schottky-defected LiNbOCl₄ as obtained from DFT simulations. Both systems exhibit a linear decrease in volume under increasing hydrostatic pressure, allowing direct determination of the bulk modulus via linear regression of the P–V data.

For pristine LiNbOCl₄, the bulk modulus is determined to be 15.98 GPa, with a standard error of 0.78 GPa, corresponding to a 95% confidence interval of 14.45–17.51 GPa. In the Schottky-defected structure, the bulk modulus is calculated to be 15.91 GPa, with a standard error of 0.47 GPa and a 95% confidence interval of 14.98–16.84 GPa. The overlap in these intervals indicates no statistically significant difference in compressibility between the two structures.

This relatively low bulk modulus, significantly lower than the ∼70 GPa obtained via Birch–Murnaghan EOS³³ fitting, aligns with values reported for other mechanically soft solid electrolytes such as Li₃PS₄, Li₆PS₅Cl, and Na₃PS₄ of 15–25 GPa.^6,37,38 This mechanical softness, arising from the intrinsic lattice flexibility and structural anharmonicity of these materials, is known to facilitate dynamic bottleneck expansion and local lattice fluctuations that enable rapid Li-ion transport.^39,40 Consequently, the comparable bulk moduli observed in both pristine and defected LiNbOCl₄ structures emphasize the inherent mechanical compliance of the Nb–O–Cl framework, supporting its applicability as a robust solid electrolyte under practical conditions.

3.4 Thermal expansion analysis

The temperature dependence of the equilibrium cell volume, extracted from NPT-MD simulations in the range of 600–1400 K, is depicted in Fig. S7 (see SI S7) for both pristine and LiCl Schottky-defected LiNbOCl₄. The data reveal a linear increase in volume with temperature in both systems.

Both pristine and Schottky-defected LiNbOCl₄ exhibit volumetric thermal expansion coefficients (α_V) that fall within the range reported for known superionic conductors. Specifically, α_V values of 1.03 × 10⁻⁴ K⁻¹ for the pristine system and 8.44 × 10⁻⁵ K⁻¹ for the Schottky-defected system are obtained from linear fits to temperature-dependent volume data. For reference, the sulfide-based electrolyte Li₁₀GeP₂S₁₂ (LGPS) has been computationally reported to have a linear thermal expansion coefficient of α_L ≈ 3.2 × 10⁻⁵ K⁻¹ at 300 K, corresponding to α_V ≈ 9.6 × 10⁻⁵ K⁻¹ assuming isotropic expansion.⁴¹ Experimental XRD studies on LGPS also show near-linear lattice expansion trends below its decomposition point (∼700 K), with minor anisotropy between a and c lattice parameters.⁴¹

In contrast, oxide-based solid electrolytes such as cubic Al-stabilized Li₇La₃Zr₂O₁₂ (LLZO) exhibit lower expansion coefficients. Synchrotron XRD measurements report α_L ≈ 15.5 × 10⁻⁶ K⁻¹ for LLZO in the temperature range of 25 °C to 700 °C, corresponding to a volumetric coefficient α_V ≈ 4.65 × 10⁻⁵ K⁻¹.⁴² When compared directly, the α_V values for both pristine and defected LiNbOCl₄ are higher than that of LLZO and comparable to those of LGPS, placing LiNbOCl₄ within the upper range of thermal expansivity observed in SEs. The absolute difference in α_V between the pristine and Schottky-defected LiNbOCl₄ systems (1.03 × 10⁻⁴ vs. 8.44 × 10⁻⁵ K⁻¹) is smaller than the full range of α_V values reported across sulfide and oxide electrolyte families, indicating that the introduction of Schottky defects leads to a moderate change in thermal expansion behavior.

To elucidate the local structural behavior underlying these macroscopic thermal responses, temperature-dependent Voronoi volume distributions and nearest-neighbor (NN) bond lengths are quantitatively analyzed (see the SI). In pristine LiNbOCl₄, the average Li-site Voronoi volume increases from 18.8 ų at 600 K to 21.5 ų at 1400 K. In the Schottky-defected structure, the Li-site Voronoi volume expands more substantially, from 19.8 ų to 22.7 ų over the same temperature range, providing direct evidence of enhanced local structural flexibility around Li sites. However, analysis of the NN bond lengths reveals only marginal changes (∼0.03 Å at 1400 K) in the Nb–O and Nb–Cl backbone, consistent with the comparable bulk moduli observed for both systems.

Quantitative analysis of temperature-dependent Voronoi volumes and bulk moduli reveals a decoupling between local structural flexibility and macroscopic thermal expansion behavior in LiNbOCl₄. Specifically, although Schottky defects induce localized volumetric expansion around Li sites, as evidenced by increased Li-site Voronoi volumes, the overall thermal expansion coefficient is slightly reduced in the defected system. This reduction occurs despite the enhanced local flexibility, suggesting that the global thermal response is predominantly governed by the comparatively rigid Nb–O–Cl framework. Such hierarchical mechanical behavior, wherein local dynamic environments coexist with a robust structural backbone, is a defining characteristic of LiNbOCl₄ as a mechanically soft SE, highlighting its potential for maintaining mechanical integrity under thermal cycling.

3.5. Analysis of Li-ion dynamics and diffusion behavior

The Li-ion dynamics in pristine and Schottky-defected LiNbOCl₄ are systematically investigated through mean squared displacement (MSD) analysis using 100 ps NVT-MD trajectories. The temperature-dependent MSD profiles are presented in Fig. 3a (pristine) and Fig. 3b (Schottky-defected). At 600 K, Li-ion motion in both systems is characterized by inter-site Li jumps (from 1st–3rd NN Li sites), with contributions from intra-site vibrational and back-and-forth jump dynamics. The MSD curves exhibit an initial ballistic regime within ∼1–2 ps, marked by a steep, near-quadratic rise due to inertial motion. This is followed by a vibrational plateau lasting up to ∼25–30 ps, reflecting confinement within local coordination environments. After ∼30 ps, a slight upward trend emerges, indicating the onset of the linear diffusive regime. The presence of Schottky defects introduces local structural distortions but does not significantly improve diffusion, confirming that the thermal energy at 600 K is insufficient to overcome migration barriers in either system.

Fig. 3 Mean squared displacement (MSD) profiles of Li ions in LiNbOCl₄ from 100 ps NVT-MD simulations. Panels (a) and (b) show the temperature-dependent MSD curves for the pristine and LiCl Schottky-defected systems, respectively, at 600 K (blue), 800 K (green), 1000 K (orange), 1200 K (red), and 1400 K (purple). Panels (c) and (d) present atomic species-resolved MSDs at 1200 K for the pristine and defected systems, respectively. In both cases, Li ions (blue solid line) exhibit substantial diffusion indicative of high ionic mobility, forming the primary contributors to long-range transport. Nb atoms (yellow dashed line) and O atoms (green dash-dotted line) show minimal displacement, remaining close to their crystallographic positions and highlighting the robustness of the host lattice. Cl atoms (red dotted line) exhibit slightly higher displacement than Nb and O, especially in the defected system, reflecting moderate thermal vibrations that are enhanced by the presence of vacancies.

As the temperature increases to 800–1400 K, a clear progression from vibrationally confined motion to long-range diffusion is observed. The MSD curves increasingly display linear behavior, indicating thermally activated hopping and eventually robust superionic conduction. The pristine system maintains slightly higher MSD slopes across the temperature range, suggesting more efficient long-range transport pathways. Meanwhile, the defected system shows comparable diffusion behavior at high temperatures but may experience mild localization effects due to structural distortion near vacancies. Despite these variations, both systems exhibit sustained Li-ion mobility at elevated temperatures, confirming the dynamic stability of LiNbOCl₄ even under defect perturbation.

3.6 Estimation of Li-ion conductivity and activation energy

3.6.1 Diffusion coefficients and the Nernst–Einstein relation. Li-ion diffusion coefficients (D_Li) are calculated from the linear regime of the mean square displacement (MSD) profiles according to eqn (5). The corresponding ionic conductivities (σ_Li) are then derived using the Nernst–Einstein relation:where c is the Li-ion carrier density, z = 1 is the ionic charge, F is the Faraday constant, R is the ideal gas constant (8.314 J mol⁻¹ K⁻¹), and T is the absolute temperature. The carrier concentration c is computed based on the number of mobile Li ions in the simulation supercell and the equilibrium volume at each temperature, determined from NPT-MD simulations.

The extracted diffusion coefficients for both pristine and Schottky-defected LiNbOCl₄ are summarized in Fig. 4. We note that the diffusion coefficient extracted from the COM-corrected MSD (e.g. at 1000 K, D ≈ 1.5 × 10⁻⁴ cm² s⁻¹) lies within the range commonly reported by AIMD studies of Li superionic conductors at comparable temperatures (10⁻⁵–10⁻³ cm² s⁻¹).^19,20,43 The temperature dependence of D_Li (see Table S6 in the SI) clearly captures the enhancement of Li-ion mobility with increasing thermal excitation. Notably, the pristine LiNbOCl₄ consistently displays higher D_Li values across the studied temperature range, as compared to the Schottky-defected structure. Although LiCl Schottky defects introduce vacancy sites that could potentially serve as transient migration pathways, our simulations indicate a net hindrance to long-range diffusion. To clarify the origin of this reduction in conductivity, analyses are performed on local anion rotational behavior, dynamic bottleneck narrowing, and attenuation of high-frequency vibrational modes (see the later section).

Fig. 4 Arrhenius plots of temperature-dependent ionic conductivity (log(σT)) for pristine and Schottky-defected LiNbOCl₄, compared with previously reported FF-MD, AIMD and experimental results. Activation energies are annotated alongside linear fits.

3.6.2 Temperature dependence and Arrhenius behavior. The temperature dependence of σ_LiT for both pristine and Schottky-defected LiNbOCl₄ is analyzed using the Arrhenius formalism:where σ₀ is the conductivity prefactor, E_a is the activation energy for Li-ion migration, and R is the gas constant.

The Arrhenius plots from MD simulations are shown in Fig. 4. Linear regression of the simulated data yielded activation energies of 0.236 eV for the pristine structure and 0.241 eV for the Schottky-defected variant. Extrapolated to room temperature (300 K), the corresponding ionic conductivities are estimated to be 9.57 × 10⁻³ S cm⁻¹ and 8.20 × 10⁻³ S cm⁻¹, respectively.

These simulation-derived values are in close agreement with experimental measurements. Tanaka et al. reported an activation energy of 0.24 eV and room-temperature ionic conductivity of 10.7 mS cm⁻¹ for cold-pressed LiNbOCl₄ pellets.¹⁰ Similarly, Jeon et al. achieved an experimental conductivity of 8.4 mS cm⁻¹ at 25 °C using a hydrochloric acid-free synthesis route and measured an activation energy consistent with values between ∼0.20 and 0.25 eV, depending on synthesis conditions.¹¹

3.7 Trajectory density analysis

Trajectory density maps obtained from 100 ps NVT-MD simulations provide detailed insights into the dynamic behavior of each atomic species in pristine and Schottky-defected LiNbOCl₄ (Fig. 5). These visualizations elucidate both the vibrational stability of the lattice framework and the emergence of Li-ion transport pathways under thermal excitation.

Fig. 5 Trajectory density maps in the c-directions of all atoms and individual atomic species, Li (magenta), Nb (blue), O (red), and Cl (cyan) in LiNbOCl₄, obtained from 100 ps NVT-MD simulations. Panels (a) and (c) correspond to the pristine structure at 600 K and 1200 K, respectively, while (b) and (d) show the LiCl Schottky-defected structure at the same temperatures. Li ions exhibit pronounced delocalization, forming continuous three-dimensional diffusion networks. In contrast, Nb and O atoms remain near their equilibrium positions, reflecting a stable and rigid crystal framework. Cl atoms show moderate vibrational broadening that becomes more pronounced with increasing temperature and in the presence of vacancies.

For both pristine and defected systems, Nb atoms exhibit minimal spatial displacement across the temperature range with trajectory densities sharply localized near their crystallographic positions. Here, “average spatial displacement” refers to the time-averaged root-mean-square (RMS) displacement of species β defined as over the 0–100 ps trajectory at each temperature. Quantitative analysis shows that the average spatial displacements of Nb atoms remain below 1.7 Å up to 1400 K in the pristine structure, increasing modestly to ∼1.5 Å in the defected system at 1200 K and 1.49 Å at 1400 K (Table S2-1). These consistently low values highlight the mechanical rigidity of the Nb–O–Cl framework and confirm the structural stability of the host lattice.

O atoms display similarly restricted mobility in both systems, with time-averaged RMS displacements below 1.7 Å up to 1200 K, underscoring the robust backbone connectivity provided by the oxygen sublattice.

Cl atoms, in contrast, exhibit progressively enhanced thermal vibrations with increasing temperature, particularly in the Schottky-defected system. In pristine LiNbOCl₄, Cl time-averaged RMS displacements increase from ∼0.8 Å at 600 K to ∼3.4 Å at 1200 K. In the defected structure, Cl displacements increase significantly from ∼1.2 Å at 600 K to ∼3.5 Å at 1200 K. This increase is attributed to the additional accessible volume created by Cl vacancies, which enhances the vibrational freedom of neighboring Cl sites.

The rotational behavior of first-nearest-neighbor (1NN) Cl and O atoms around the Nb polyhedron of the (potential) vacancy site is analyzed using spherical coordinates (r, θ, ϕ) over the final 5 ps of MD simulations at 600 K (Fig. 6). Here, r denotes the instantaneous Nb–X distance (X = Cl, O); θ is the polar (tilt) angle measured from the crystallographic c axis; and ϕ is the azimuthal angle between the projection of the Nb–X vector onto the ab plane and the a axis. In the pristine system (Fig. 6a–c), the Cl neighbors, Cl-near (closest to the potential vacancy site) and Cl-far (second closest), exhibit moderate and coherent reorientation. The Cl-near atom has r fluctuating within 2.2–2.6 Å and θ and ϕ confined to ∼80–120° and ∼100°, respectively. The Cl-far atom displays similar r variations but with ϕ extending down to ∼−100°, suggesting slightly more flexible motion. The oxygen neighbor displays pronounced azimuthal reorientation (ϕ excursions approaching ±100°) at nearly constant, modest tilt (θ ∼ 5–25°) and limited Nb–O bond-length fluctuation (r ∼ 1.7–2.1 Å). Thus, O predominantly precesses about c within a small tilt cone, while the octahedral geometry remains intact.

Fig. 6 Spherical coordinate representation of anion rotational dynamics around vacancy-adjacent Nb polyhedra. Time evolution of spherical coordinate components r, θ, and ϕ for first-nearest-neighbor anions, Cl-near (blue), Cl-far (orange), and O (green), surrounding the Nb center in pristine (a–c) and Schottky-defected (d–f) LiNbOCl₄, tracked over the final 5 ps of the MD simulation at 600 K. The local coordinate frame and structural reference are shown in panel (g), where the spherical system is centered on Nb; θ is the polar (tilt) angle from the crystallographic c axis; ϕ is the azimuth about c; the Cl vacancy lies at θ = 90°, ϕ = 0°, oriented along the crystallographic a-axis. In the pristine structure, O shows large azimuthal excursions (ϕ up to ±100°) at modest, nearly constant tilt (θ ∼ 5–25°) and narrow Nb–O distance fluctuations (r ∼ 1.7–2.1 Å), i.e., precession about c within a small tilt cone. In contrast, the Schottky-defected system displays larger-amplitude, erratic, and spatially uncorrelated rotations (shown by red arrows), particularly in ϕ, with abrupt shifts and excursions approaching ±180° (shown by the red circle). These behaviors highlight the loss of rotational synchrony and emergence of anharmonic distortions in the vacancy-adjacent environment, consistent with disrupted local lattice dynamics.

In the Schottky-defected system, the local structural environment around Nb transitions from octahedral to a five-fold pyramidal coordination, typical of Nb centers lacking one anionic ligand. The Cl neighbors (Cl-near and Cl-far) in this environment display significantly larger and more erratic rotational dynamics, with r oscillating up to 2.8 Å, θ extending to ∼120°, and ϕ spanning nearly 180°. These pronounced and unsynchronized rotational motions highlight the disruption of collective coherence within the polyhedral framework due to the vacancy-provided volume. O neighbors in the defected structure also exhibit enhanced rotational flexibility compared to pristine, with r extending from 1.8–2.4 Å, θ reaching ∼50°, and ϕ displaying broad rotations of ±120°. While the amplitude of angular motion increases in the defect (broader |ϕ|, wider θ cone, and larger r fluctuations), the reorientation rate is lower (longer dwell between bursts; see Section S9 in the SI). At 1000 K, the overlapping (sliding-window) azimuthal switching frequency (1 ps window; 0.1 ps stride) averages 〈ν_ϕ〉_overlap ≈ 1.78 ps⁻¹ in the pristine cell and ≈1.47 ps⁻¹ in the Schottky cell. Together, these observations indicate that the enhanced local vibrational/rotational freedom in the defected system is less temporally synchronized, reflecting a loss of collective coherence within the polyhedral network.

To directly link these local anion dynamics to Li-ion transport, we next examine a representative long-range migration event. As illustrated in Fig. 7, the Li-ion ion exhibits long-range migration across the bc plane, accompanied by coordinated displacement of surrounding Cl atoms within the NbO₂Cl₄ framework. These results confirm that Li-ion diffusion occurs in concert with the dynamic framework response. The Cl atoms undergo significant angular distortions and cooperative motion during Li-ion jumps, indicating a soft vibrational environment.^44,45 This behavior is consistent with previous simulations that reported pronounced angular fluctuations in the Cl–Nb–Cl bonds and overlapping Li–Cl vibrational modes in the low-energy phonon region.^44,45 While Cl⁻ ions do not move translationally, their local dynamical motion contributes to transient bottleneck widening during Li-ion migration.

Fig. 7 Representative AIMD snapshot at 600 K (50–70 ps). The magenta trace marks a diffusing Li-ion ion in the bc-plane, while the surrounding framework atoms, Cl (teal), Nb (blue), and O (red), undergo cooperative displacements. This highlights a dynamic coupling between Li-ion migration and local octahedral distortions within the NbO₂Cl₄ framework. This visualization approach is consistent with prior studies of framework–ion interactions in similar systems.⁴⁵

Trajectory density maps and average spatial displacements demonstrate that Li-ion delocalization in both pristine and LiCl-Schottky-defected LiNbOCl₄ evolves significantly with temperature. Following the observed enhancement in Cl and O rotational amplitudes, particularly in the defected structure, where polyhedral coherence is locally diminished, the Li sublattice exhibits correspondingly increased dynamic behavior. In the pristine system, Li-ion motion increases steadily, with time-averaged RMS displacements increasing from 8.64 Å at 600 K to 22.54 Å at 1200 K, and further to 27.46 Å at 1400 K. This trend reflects the development of extended diffusion pathways, consistent with the emergence of three-dimensional network connectivity observed in trajectory maps along the c-direction (Fig. 5). Similarly, in the Schottky-defected system, Li-ion time-averaged RMS displacements increases from 8.80 Å at 600 K to 27.69 Å at 1200 K, suggesting comparable long-range mobility at elevated temperatures. Although the trajectory density distribution in the defected system appears less uniform at lower temperatures, it becomes increasingly connected with increasing temperature. Nonetheless, both systems show extensive thermal delocalization above 1000 K, confirming that defect introduction does not significantly hinder global Li-ion transport under high-temperature conditions.

3.8 Factors governing the conductivity reduction in Schottky-defected LiNbOCl₄

The ionic conductivity of SEs is typically expressed through the Arrhenius relationship, where the conductivity prefactor σ₀ encapsulates critical microscopic descriptors of ion transport. According to eqn (9),where n is the Li-ion concentration, z is related to the diffusion geometry and correlation, Ze is the ionic charge, a₀ is the jump distance, ν₀ is the vibrational attempt frequency, k_B is the Boltzmann constant and ΔS_m is the migration entropy. The prefactor thus integrates both geometric and entropic contributions to ion migration: the product (a₀)²ν₀ reflects the intrinsic hopping dynamics, while the exponential term captures the entropic gain from the availability and accessibility of diffusion pathways. Here, a detailed analysis of σ₀ is essential for understanding the fundamental origin of transport suppression in Schottky-defected LiNbOCl₄; it also enables us to disentangle geometric advantages from entropic and dynamical penalties.

Despite structural indicators suggesting improved diffusion accessibility in the Schottky-defected LiNbOCl₄, namely, an enlarged average bottleneck radius (1.9 Å vs. 1.5 Å in the pristine phase; see Section S8 in the SI) and a marginally increased Li-site Voronoi volume, it still shows a slightly higher activation energy (E_aSchottky = 0.241 eV vs. E_aPristine = 0.236 eV). Meanwhile, the extrapolated prefactor is noted to increase from 2.60 × 10⁴ to 2.72 × 10⁴ S cm⁻¹ K, but this minor gain is offset by the exponential penalty imposed by the higher barrier. The calculated Boltzmann factor ratio at 300 K, exp(−E_aSchottky/k_BT)/exp(−E_aPristine/k_BT) ≈ 0.818, confirms that the conductivity reduction is largely owed to the increased activation energy.

To probe the microscopic origin of the increased activation energy observed in the Schottky-defected system, we examine the rotational dynamics of the anionic sublattice and its influence on Li-ion mobility. The Bond Angle Correlation Function (BACF) further probes host lattice dynamics by evaluating angular stability of anion bonds around the Nb polyhedron of the potential vacancy site. Defined as

BACF(τ) = 〈v̂(t)·v̂(t + τ)〉 (10)

where v̂ is the normalized bond vector at time t, this function measures the temporal persistence of anion orientations. The BACFs of Nb–O and Nb–Cl coordination environments reveal a distinct divergence in angular stability between pristine and defected phases. As shown in Fig. 8, the BACF for Nb–Cl in the Schottky-defected structure decays rapidly, reaching zero and subsequently negative values within ∼9 ps. The appearance of negative correlations (minimum ≈−0.4) indicates bond vector inversion and loss of coherent rotational behavior, suggesting dynamic instability in the local anion environment. In contrast, the pristine phase maintains a high degree of rotational order, with BACF values remaining above 0.9 throughout the entire 80 ps simulation window. This persistent angular correlation reflects stable, cage-like oscillations of the anion framework around Nb polyhedra, consistent with quasi-harmonic behavior. The disrupted anion reorientation in the Schottky-defected structure likely contributes to increased Li-ion hopping barriers by narrowing diffusion bottlenecks, thereby providing a mechanistic basis for the elevated activation energy.

Fig. 8 Bond orientation correlation of Cl and O anions in pristine and Schottky-defected LiNbOCl₄. Bond angle correlation functions (BACFs) of first-nearest-neighbor anions, Cl-near (blue), Cl-far (orange), and O (green), coordinated to Nb in pristine (a) and Schottky-defected (b) LiNbOCl₄, evaluated over an 80 ps MD trajectory at 600 K. In the pristine system, BACF values remain near unity across all anions, indicating stable and long-lived angular orientation around the Nb polyhedron. In contrast, the defected system shows rapid decay in Nb–Cl BACF values, reaching minima near −0.4, reflecting inversion of bond vectors and breakdown of angular coherence. O atoms exhibit slightly damped but more persistent correlations. These trends confirm enhanced rotational disorder and lattice anharmonicity in the presence of Schottky defects.

These disruptions are further supported by qualitative insights from spherical trajectory analysis of Cl and O neighbors around Nb centers, as previously described (Fig. 6). In the pristine structure, anion rotations remain relatively constrained, with angular fluctuations exhibiting smooth and confined patterns indicative of coherent polyhedral breathing. In contrast, the Schottky-defected system displays irregular and spatially uncorrelated angular motions, particularly in the azimuthal component ϕ, with abrupt directional changes and excursions across a broader angular space. This loss of synchronicity reflects enhanced rotational disorder within the local coordination shells, consistent with the rapid BACF decay and inversion of bond vectors. These results signify a transition from harmonic to anharmonic rotation, where weakened local coordination leads to irregular anion dynamics. Importantly, such incoherent behavior does not translate into transport-effective dynamic gating: irregular, locally confined rotations are less efficient in generating the transient bottleneck openings required for Li migration.^46–48 Unlike the classical paddle-wheel mechanism involving rotating polyanion clusters, the dynamic gating mechanism introduced in this work describes a distinct mode of transient bottleneck expansion driven by coordinated motion of individual halide and oxide anions within an oxyhalide lattice (see SI S8), as quantified more explicitly by the event-triggered ensemble (ETE) analysis in Fig. 9.

Fig. 9 Event-triggered ensemble (ETE) Li-ion mean-squared displacement (MSD) in pristine and LiCl-Schottky-defected LiNbOCl₄ at 600 K. For each azimuthal reorientation (“φ-switch”) event at time t_i of the tracked NbX₆ polyhedron, the Li MSD is evaluated as 〈Δr²(τ)〉 = 〈|r(t_i + τ) − r(t_i)|²〉 and averaged over all events, with the delay τ = 0 marking the anion event (vertical dashed line). Li ions within 4.5 Å of the active polyhedron at t_i are classified as “near” (blue), and all remaining Li ions as “far” (orange). (a) In pristine LiNbOCl₄, near and far ETE–MSD curves are nearly indistinguishable for |τ| ≲ 0.5 ps, indicating that anion reorientation events are embedded in a coherently responding network that drives spatially delocalized Li motion. (b) In the Schottky system, Li-ions in the vacancy-adjacent environment remain more localized than the far-field population (near < far both before and after τ = 0), consistent with stronger caging around the defect and with the reduced long-time MSD slope and diffusivity of the Schottky-defected phase.

To directly test whether local anion reorientation acts as a dynamic gate for Li motion, we applied a time-aligned, event-triggered ensemble (ETE) analysis to the spherical-coordinate descriptors {r, θ, ϕ} introduced above. Around a tracked NbX₆ polyhedron, chosen adjacent to the LiCl vacancy in the Schottky case and at the same crystallographic site in the pristine cell prior to vacancy introduction, we construct several transport-relevant rotation/geometry descriptors (details in Section S10): an azimuthal switching rate ν_φ(t), which counts abrupt φ reorientation events per picosecond via unwrapping, median-centering, and a Schmitt trigger with hysteresis; a tilt-weighted azimuthal activity A₁(t) = |dφ/dt|sin θ, which emphasizes rotations that tilt the ligands toward the diffusion passage; a time-dependent bottleneck half-gap a(t), defined for a chosen anion pair that frames a local neck (larger a corresponds to a wider gate); and a normalized composite “gate” score G(t) that combines bond stretching and tilting through r and θ.

In the ETE framework, we first identify discrete anion-framework events around each NbX₆ octahedron, for which here we focus on azimuthal “ϕ-switch” events as the elementary anion-framework events. This procedure yields a reproducible list of timestamps {t_i} at which candidate transport-enabling openings occur, where the Schmitt state s(t) of the tracked NbX₆ octahedron flips sign, indicating an abrupt azimuthal reorientation (see SI S10). For each event time t_i we then compute the Li MSD relative to the configuration at the event,

〈Δr²(τ)〉 = 〈|r(t_i + τ) − r(t_i)|²〉,

as a function of delay τ, with τ = 0 marking the anion event. To resolve the spatial structure, Li ions are partitioned at each event into a “near” ensemble (within 4.5 Å of the tracked NbX₆ octahedron at t_i) and a complementary “far” ensemble. The far ions thus provide an internal reference for the background, system-wide motion away from the active octahedron. Averaging 〈Δr²(τ)〉 over all events yields near and far ETE–MSD curves, Δr²(τ)_near and Δr²(τ)_far (Fig. 9). By construction, both curves vanish at τ = 0, because displacements are measured relative to the event configuration, so all ETE curves exhibit a V-shaped minimum at τ = 0. For τ < 0, the curves describe how Li configurations approach the event from earlier times; for τ > 0, they quantify how rapidly Li positions decorrelate after the anion reorientation. Since events are selected based on a directional criterion on ϕ(t) (via the Schmitt trigger) rather than based on a symmetric oscillation about a reference angle, there is no requirement for the ETE curves to be mirror-symmetric about τ = 0. The physically relevant information lies in how the near and far branches differ on the pre-event (τ < 0) and post-event (τ > 0) sides.

This analysis reveals a clear contrast between pristine and LiCl-Schottky LiNbOCl₄ (Fig. 9). Over the 100 ps window at 600 K, the pristine cell exhibits several hundred ϕ-switch events around the reference octahedron, whereas the Schottky cell shows only a few tens of such events, indicating that dynamically active anion reorientation is strongly suppressed by the defect. In the pristine system, the near and far ETE–MSD curves are almost indistinguishable over |τ| ≲ 0.5 ps. Both ensembles display very similar growth of 〈Δr²(τ)〉 on the τ < 0 and τ > 0 sides, with only minor deviations. This similarity indicates that reorientation events of the reference NbX₆ octahedron are embedded in a collectively responding anion–Li network: when the octahedron reorients, the associated Li displacements are spatially delocalized, and the local environment of the tracked octahedron does not stand out as a special region of enhanced mobility.

By contrast, in the LiCl-Schottky phase, ϕ-switch events are rarer, and when they occur, far-field Li consistently exhibits larger MSD than vacancy-adjacent Li, both before and after the event (near < far for τ < 0 and τ > 0). This behavior can be naturally decomposed into a background contribution and an event-triggered contribution: the background describes the persistent difference in mobility between near and far Li that would exist even in the absence of a specific event, whereas the event-triggered part captures the additional change in MSD induced by the ϕ-switch itself. At the background level, the persistent inequality near < far on both sides of τ = 0 reflects the fact that Li adjacent to the vacancy reside in a structurally and dynamically constrained environment with intrinsically low mobility, whereas Li in less distorted regions of the framework sample a softer network with a higher underlying MSD. Superimposed on this background, the ϕ-switch produces a dynamic modulation of the local motion: for the near ensemble, the post-event branch 〈Δr²(τ > 0)〉_near lies systematically above the pre-event branch 〈Δr²(τ < 0)〉_near, showing that anion reorientation slightly loosens the local cage and promotes additional decorrelation of vacancy-adjacent Li after the event. However, the near curve never catches up to the far curve over the |τ| ≤ 0.5 ps window, indicating that even with this event-triggered “kick”, Li in the vacancy-adjacent region remain less mobile overall than Li in the far field.

In the present context, the ETE results thus provide a direct, time-resolved view of the dynamic gating mechanism. Here, “dynamic gating” refers to Li transport being controlled by transient openings of local bottlenecks formed by surrounding anions: rotations and vibrations of the NbX₆ octahedra temporarily widen or reorient these bottlenecks, allowing Li to move more easily before the gate recloses. In the ETE framework, each ϕ-switch is treated as a discrete anion reorientation event, and the near/far, τ-resolved MSD tracks how Li motion responds to that local gate dynamics in space and time. In pristine LiNbOCl₄, anion reorientations around the reference octahedron are frequent and elicit a delocalized Li response (near ≈ far), which together underpin the larger long-time MSD slope and higher diffusivity. In the LiCl-Schottky phase, comparable events around the vacancy-adjacent octahedron are both less frequent and less effective at mobilizing nearby Li (near < far), so their contribution to forming extended percolation pathways is limited, consistent with the reduced long-time MSD slope and lower computed diffusivity. These time-resolved observations show that the Schottky defect does not merely modify static bottleneck dimensions, but alters the frequency, coherence, and spatial reach of dynamically active anion reorientation, a change that, as we show below, is also reflected in the underlying vibrational spectrum (phonon density-of-states analysis).

The phonon density of states (DOS) at 600 K (Fig. 10) provides a complementary, frequency-domain view of these dynamical differences. In the defected system, all atomic species exhibit redshifts in their vibrational band centers, Li (32.66 → 31.14 meV), Nb (26.26 → 25.26 meV), Cl (20.48 → 20.14 meV), and O (57.74 → 56.17 meV), consistent with local lattice softening. Notably, a high-frequency peak at ∼80 meV, primarily arising from O vibrations and linked to fast vibrational modulation of Li-ion pathways, is significantly attenuated relative to the pristine structure. This suppression implies diminished access to thermally activated bottleneck widening, which in turn can increase the effective energy barrier.¹

Fig. 10 Phonon density of states (DOS) of (a) pristine and (b) LiCl Schottky-defected LiNbOCl₄ at 600 K derived from MD simulations. The total DOS is shown in black, with atomic contributions from Li (blue), Nb (orange), Cl (green), and O (red). The pristine system exhibits a sharper high-frequency peak (∼80 meV) primarily arising from O vibrations, which is significantly attenuated in the defected system. Additionally, all atomic species show redshifts in vibrational band centers in the defected phase, consistent with local lattice softening.

The partial phonon density of states (PDOS) in Fig. 10 reveals that Cl atoms contribute significantly to the low-energy vibrational modes below ∼18 meV, with a computed band center at around 20.14 meV. This value is in close agreement with previously reported results for LiNbOCl₄, where the Cl phonon band center was found near 19 meV.⁴⁵ Additionally, we observe a strong spectral overlap between the Cl and Li components in this low-frequency range, suggesting a substantial degree of vibrational coupling. These findings support the interpretation that Li-ion transport is dynamically facilitated by framework fluctuations, particularly involving the Cl sublattice. This framework-fluctuation picture is further supported, albeit indirectly, by experiments: Raman spectroscopy on LiNbOCl₄ shows distinct low-frequency bands assigned to Nb–Cl modes and reveals orientational disorder and dynamic rotational flexibility within the NbO₂Cl₄ units.⁴⁴ Taken together, these computational and experimental observations indicate that Li-ion migration proceeds in a soft, dynamically fluctuating NbO₂Cl₄ framework dominated by low-energy Cl motions.

To quantify how this framework-mediated dynamics is encoded in the macroscopic transport coefficients, we next examine the Arrhenius parameters, with particular emphasis on the conductivity prefactor σ₀. While the elevated activation energy clearly governs the conductivity suppression in the Schottky-defected phase, a more detailed decomposition of the conductivity prefactor σ₀ offers insight into how microscopic transport descriptors respond to defect-induced changes. Building on eqn (9), we evaluated the relative contributions of the effective jump distance a₀, attempt frequency ν₀, carrier density n, and correlation factor z, where the latter is taken to be approximately equivalent across both systems due to their nearly identical local environments. At 600 K, the Schottky-defected system exhibits a slightly shorter average jump distance (a₀ = 2.63 Å) than the pristine counterpart (a₀ = 2.78 Å), which would normally act to suppress the prefactor. However, this is compensated for by an ∼15.6% increase in the attempt frequency ν₀, from 3.98 × 10¹² s⁻¹ (pristine) to 4.52 × 10¹² s⁻¹ (Schottky). The carrier density n is marginally reduced in the Schottky-defected structure due to Li vacancies, but this effect remains under 5%. On the other hand, the ratio of correlation factors z between the two systems is expected to be close to unity and should also have a negligible contribution to the overall prefactor variation (see discussion below related to van Hove analysis).

When combining the different geometric and dynamical contributions, the overall increase in the conductivity prefactor σ₀ for the Schottky system remains under 5%. However, the extrapolated prefactor obtained from Arrhenius analysis shows a more noticeable enhancement. This discrepancy indicates an additional contribution from migration entropy. From the prefactor ratio and known parameters, we estimate that the entropic term exp(ΔS_m/k_B) increases by approximately 12% in the Schottky-defected system. While this reflects a slight broadening of the diffusion pathway landscape due to local structural disorder, the contribution remains modest and ultimately insufficient to overcome the conductivity suppression caused by the elevated migration barrier. These results highlight that, in this case, the noted entropic gain is insufficient to overcome the defect-related energetic penalty for Li ion migration in Schottky-defected LiNbOCl₄.

We further analyze the ion dynamics in pristine and Schottky-defected LiNbOCl₄ using the van Hove correlation function, G(r, t), derived from MD trajectories at 600 K (Fig. 11). This function quantifies the time-dependent probability of finding a particle at a position r at time t, given that another (or the same) particle was at the origin at time t = 0. It is mathematically defined as:

where the first term corresponds to the self-part, G_s(r, t), and the second to the distinct-part, G_d(r, t). The Dirac delta function δ(·) ensures positional matching, and 〈·〉 denotes the ensemble average over all ions.

Fig. 11 van Hove self and distinct correlation functions for Li-ion dynamics in pristine and Schottky-defected LiNbOCl₄. Two-dimensional van Hove correlation function heatmaps for Li-ion motion in pristine (left) and Schottky-defected (right) LiNbOCl₄ at 600 K over a 100 ps MD trajectory. The self-part G_s(r, t) (top) quantifies the probability of a Li-ion moving a distance r after time t, revealing long-range delocalization in the pristine phase and confined motion in the defected phase. The distinct-part G_d(r, t) (bottom) reflects pairwise spatial correlations between Li ions. In both systems, discontinuous and temporally fluctuating peaks near r ≈ 0–2 Å are observed, indicating uncorrelated, non-concerted hopping behavior. The defected system exhibits even more intermittent features, implying further disruption of spatial correlation and lower transport cooperativity.

The self-part, G_s(r, t), captures the probability that a Li-ion has moved a distance r over time t, reflecting the extent of individual ion displacement. In the pristine system (as shown in Fig. 11), the G_s(r, t), distribution exhibits long-tailed profiles with nonzero probabilities persisting beyond r > 6 Å, indicative of long-range, percolative Li-ion diffusion. By contrast, the defected system shows more localized G_s(r, t) distributions, with peak intensities centered near r = 2 Å, persistent bands around r = 1 Å and rapidly decaying tails, consistent with confined motion, reduced jump distances, and increased back-jumping probability, hallmarks of restricted ionic mobility.

The distinct-part, G_d(r, t), evaluates the likelihood of finding a second Li-ion at distance r from the initial position of a reference Li-ion after time t, thereby characterizing pairwise correlations and cooperative motion. Our analysis reveals that both pristine and Schottky-defected LiNbOCl₄ exhibit discontinuous and temporally fluctuating peak intensities at 600 K near r ≈ 0–2 Å. These peaks do not persistently build up over time but rather appear and vanish intermittently, signaling weak ion–ion correlations. This behavior is slightly different from the persistent G_d(r, t) band in the same r range that is seen in other fast superionic conductors such as Na₃SbS₄.⁴⁹ Such behavior in the latter is indicative of relatively sustained jump correlations among Na⁺ ions, where consecutive ions rapidly occupy neighboring vacant sites, facilitating collective motion. In the LiNbOCl₄ system, however, these peaks emerge intermittently and decay without persistent growth, suggesting a reduced degree of correlation. Nevertheless, the non-zero intensities around r = 0 across multiple time frames indicate that the dynamics are not entirely random. Rather than purely Brownian, the Li-ion transport exhibits a mixture of correlated and uncorrelated jumps, where temporal coherence in motion may still exist but is weakened. This intermediate regime reflects partially cooperative behavior, which may still support efficient diffusion pathways without manifesting as fully concerted motion.

However, a closer comparison between pristine and defected LiNbOCl₄ reveals subtle but significant differences in their G_d(r, t) profiles. The pristine system shows relatively denser and earlier-emerging intensity near r = 0 Å across time, indicating a modest degree of temporally correlated motion. In contrast, the defected system exhibits fainter, more sporadic features with delayed onset and weaker central intensities, pointing to diminished temporal coherence in Li-ion movement. This contrast implies that while both systems operate within a partially cooperative transport regime, the presence of Schottky defects disrupts ion–ion correlations, weakening the transiently collective behavior. Such disruptions may arise from local site distortions or the loss of rotational coherence (as previously mentioned) in surrounding anion polyhedra, ultimately contributing to the observed reduction in conductivity and increased activation barrier.

These findings related to G_s(r, t) and G_d(r, t) indicate that the conductivity suppression in the Schottky-defected system stems not from a loss of collective motion per se, but from the disruption of local structural dynamics that facilitate such behavior. Although the static bottleneck geometry appears slightly more open in the defected system, this potential advantage is negated by incoherent anion rotations and a breakdown of synchrony within the surrounding polyhedral framework. This loss of rotational coherence weakens the dynamic gating mechanism, an essential contributor to transient bottleneck expansion, which in turn reduces the effective hopping volume accessible to Li-ions. Consequently, even within a regime of cooperative transport, the defect-induced disturbance in local lattice dynamics imposes a higher energy cost for migration, ultimately manifesting as an increased activation barrier and diminished ionic conductivity.

In summary, although the Schottky-defected LiNbOCl₄ system retains similar global elastic properties and diffusion topology to its pristine counterpart, its conductivity is reduced due to a localized increase in activation energy driven by disrupted anion rotational dynamics, and suppression of key phonon modes. These results highlight the critical role of cooperative anion motion in modulating the dynamic energy landscape of fast-ion conductors and provide insight into the defect tolerance limits of oxyhalide-based solid electrolytes.

4. Conclusion

In this work, we used AIMD to elucidate Li-ion transport in pristine and LiCl Schottky-defected LiNbOCl₄. Pristine LiNbOCl₄ exhibits robust Li mobility over 600–1400 K, with an activation energy of 0.236 eV and an extrapolated room-temperature ionic conductivity of 9.57 × 10⁻³ S cm⁻¹, in close agreement with the experiment. The rigid Nb–O–Cl framework provides mechanical stability, while the disordered Li sublattice supports percolating diffusion pathways, consistently captured by MSD, trajectory-density, and van Hove analyses.

Introducing LiCl Schottky defects preserves the global diffusion topology, mechanical softness, and metastability, but slightly increases the activation energy to 0.241 eV and lowers the room-temperature conductivity to 8.20 × 10⁻³ S cm⁻¹. Although static descriptors (e.g., a larger bottleneck radius and increased local free volume) suggest more accessible migration channels, dynamic analyses show the opposite trend: the defects induce pronounced anion rotational disorder and weaken coherent “gate-opening” fluctuations. Event-triggered ensemble analysis reveals that transport-enabling anion events become rarer and more spatially localized, limiting their ability to promote delocalized Li motion. Consistently, phonon DOS shows overall redshifts and a strong attenuation of the O-dominated mode near ∼80 meV associated with dynamic bottleneck modulation. Changes in the conductivity prefactor are modest (≤5% net increase, plus ∼12% from migration entropy) and cannot compensate for the exponential penalty from the slightly higher barrier; van Hove functions further indicate reduced temporal coherence and weaker Li–Li correlations in the defected phase.

Overall, our results demonstrate that LiNbOCl₄ possesses a highly defect-tolerant, mechanically compliant framework that maintains superionic Li conductivity even under LiCl Schottky perturbations. The slight conductivity reduction in the defected phase arises not from a loss of global percolation or catastrophic structural instability, but from a subtler reshaping of the dynamic energy landscape: rotational incoherence, attenuation of key phonon modes, and spatially localized gate-opening events collectively increase the effective migration barrier. These insights highlight the central role of cooperative anion motion and dynamic bottleneck modulation in governing ion transport, and they provide concrete design guidelines for engineering oxyhalide-based solid electrolytes with tailored defect chemistry and lattice dynamics.

Author contributions

Halimah Harfah: methodology, software, investigation, formal analysis, visualization, writing – original draft. Yoshitaka Tateyama: supervision, writing – review & editing. Kazunori Takada: supervision, writing – review & editing. Randy Jalem: conceptualization, supervision, project administration, funding acquisition, writing – review & editing, resources.

Conflicts of interest

There are no conflicts to declare.

Data availability

The data supporting this study are available from the corresponding author upon reasonable request.

Supplementary information (SI): Ewald and total energy screening of Li configurations and LiCl Schottky-defect structures, defect formation energies and convex hull stability analysis, DFT and AIMD computational details, average spatial displacements, Voronoi volumes, and nearest-neighbor bond-length distributions, mechanical properties, bulk modulus, and thermal expansion, Li-ion diffusion coefficients, static and dynamic bottleneck analyses, angle analysis, azimuthal switching frequencies, and transport-relevant gate descriptors. See DOI: https://doi.org/10.1039/d5ta05478h.

Acknowledgements

This work was financially supported in part by JST through Green Technologies of Excellence (GteX) grant number JPMJGX23S2, by JSPS KAKENHI grant number JP21K14729, MEXT as Materials Processing Science project (“Materealize”) grant number JPMXP0219207397 and the “Program for Promoting Research on the Supercomputer Fugaku” grant number JPMXP1020230325. The calculations were carried out both in the NIMS supercomputer (Numerical Materials Simulator) and the Fugaku supercomputer at the RIKEN through the HPCI System Research Project (project ID: hp240118, hp210105).

References

1. J. C. Bachman, S. Muy, A. Grimaud, H.-H. Chang, N. Pour, S. F. Lux, O. Paschos, F. Maglia, S. Lupart, P. Lamp, L. Giordano and Y. Shao-Horn, Inorganic Solid-State Electrolytes for Lithium Batteries: Mechanisms and Properties Governing Ion Conduction, Chem. Rev., 2016, 116, 140–162.
2. Z. Zhang, Y. Shao, B. Lotsch, Y.-S. Hu, H. Li, J. Janek, L. F. Nazar, C.-W. Nan, J. Maier, M. Armand and L. Chen, New horizons for inorganic solid state ion conductors, Energy Environ. Sci., 2018, 11, 1945–1976.
3. M. Umair, S. Zhou, W. Li, H. T. H. Rana, J. Yang, L. Cheng, M. Li, S. Yu and J. Wei, Oxide Solid Electrolytes in Solid-State Batteries, Batteries Supercaps, 2024, e202400667.
4. A. Banerjee, X. Wang, C. Fang, E. A. Wu and Y. S. Meng, Interfaces and Interphases in All-Solid-State Batteries with Inorganic Solid Electrolytes, Chem. Rev., 2020, 120, 6878–6933.
5. C. Xu, L. Chen and F. Wu, Unveiling the power of sulfide solid electrolytes for next-generation all-solid-state lithium batteries, Next Mater., 2025, 6, 100428.
6. Y. Kato, S. Hori, T. Saito, K. Suzuki, M. Hirayama, A. Mitsui, M. Yonemura, H. Iba and R. Kanno, High-power all-solid-state batteries using sulfide superionic conductors, Nat. Energy, 2016, 1, 16030.
7. T. Asano, A. Sakai, S. Ouchi, M. Sakaida, A. Miyazaki and S. Hasegawa, Solid Halide Electrolytes with High Lithium-Ion Conductivity for Application in 4 V Class Bulk-Type All-Solid-State Batteries, Adv. Mater., 2018, 30, 1803075.
8. B. Tao, D. Zhong, H. Li, G. Wang and H. Chang, Halide solid-state electrolytes for all-solid-state batteries: structural design, synthesis, environmental stability, interface optimization and challenges, Chem. Sci., 2023, 14, 8693–8722.
9. Q. Wang, Y. Zhou, X. Wang, H. Guo, S. Gong, Z. Yao, F. Wu, J. Wang, S. Ganapathy, X. Bai, B. Li, C. Zhao, J. Janek and M. Wagemaker, Designing lithium halide solid electrolytes, Nat. Commun., 2024, 15, 1050.
10. Y. Tanaka, K. Ueno, K. Mizuno, K. Takeuchi, T. Asano and A. Sakai, New Oxyhalide Solid Electrolytes with High Lithium Ionic Conductivity >10 mS cm⁻¹ for All-Solid-State Batteries, Angew. Chem., Int. Ed., 2023, 62, e202217581.
11. S. Jeon, K.-H. Park, W. Cho, G. Jeong, J. Yu, Y. J. Park and K. Kim, Hydrochloric acid-free synthesis of LiNbOCl₄ superionic conductor for all-solid-state Li batteries, Solid State Ionics, 2025, 421, 116791.
12. Q. Wang, Z. Shen, P. Du, Y. Zhou, P. Zhang and Y. Liu, New advances in solid-state electrolytes: from halides to oxyhalides, Inorg. Chem. Front., 2024, 11, 5810–5832.
13. S. Adams, Origin of Fast Li+-Ion Conductivity in the Compressible Oxyhalide LiNbOCl4, Energy Storage Mater., 2024, 68, 103359.
14. B. Gao, R. Jalem and Y. Tateyama, First-Principles Study of Microscopic Electrochemistry at the LiCoO2Cathode/LiNbO3Coating/β-Li3PS4Solid Electrolyte Interfaces in an All-Solid-State Battery, ACS Appl. Mater. Interfaces, 2021, 13, 11765–11773.
15. K. V. Kravchyk, H. Zhang and M. V. Kovalenko, On the interfacial phenomena at the Li7La3Zr2O12 (LLZO)/Li interface, Commun. Chem., 2024, 7, 257.
16. X. Yang, J. Luo and X. Sun, Towards high-performance solid-state Li–S batteries: from fundamental understanding to engineering design, Chem. Soc. Rev., 2020, 49, 2140–2195.
17. J. Lau, R. H. DeBlock, D. M. Butts, D. S. Ashby, C. S. Choi and B. S. Dunn, Sulfide Solid Electrolytes for Lithium Battery Applications, Adv. Energy Mater., 2018, 8, 1800933.
18. Y. Seino, T. Ota, K. Takada, A. Hayashi and M. Tatsumisago, A sulphide lithium super ion conductor is superior to liquid ion conductors for use in rechargeable batteries, Energy Environ. Sci., 2014, 7, 627–631.
19. Y. Mo, S. P. Ong and G. Ceder, First Principles Study of the Li ₁₀ GeP₂S₁₂ Lithium Super Ionic Conductor Material, Chem. Mater., 2012, 24, 15–17.
20. S. P. Ong, Y. Mo, W. D. Richards, L. Miara, H. S. Lee and G. Ceder, Phase stability, electrochemical stability and ionic conductivity of the Li_10±1MP₂X₁₂ (M = Ge, Si, Sn, Al or P, and X = O, S or Se) family of superionic conductors, Energy Environ. Sci., 2013, 6, 148–156.
21. S. P. Ong, W. D. Richards, A. Jain, G. Hautier, M. Kocher, S. Cholia, D. Gunter, V. L. Chevrier, K. A. Persson and G. Ceder, Python Materials Genomics (pymatgen): A robust, open-source python library for materials analysis, Comput. Mater. Sci., 2013, 68, 314–319.
22. G. Kresse and J. Furthmüller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169–11186.
23. P. E. Blöchl, Projector augmented-wave method, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953–17979.
24. J. P. Perdew, K. Burke and M. Ernzerhof, Generalized Gradient Approximation Made Simple, Phys. Rev. Lett., 1996, 77, 3865–3868.
25. S. Grimme, Semiempirical GGA-type density functional constructed with a long-range dispersion correction, J. Comput. Chem., 2006, 27, 1787–1799.
26. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu, J. Chem. Phys., 2010, 132, 154104.
27. S. Grimme, S. Ehrlich and L. Goerigk, Effect of the damping function in dispersion corrected density functional theory, J. Comput. Chem., 2011, 32, 1456–1465.
28. J. Klimeš, D. R. Bowler and A. Michaelides, Van der Waals density functionals applied to solids, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 83, 195131.
29. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth and B. I. Lundqvist, Van der Waals Density Functional for General Geometries, Phys. Rev. Lett., 2004, 92, 246401.
30. K. Lee, É. D. Murray, L. Kong, B. I. Lundqvist and D. C. Langreth, Higher-accuracy van der Waals density functional, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 82, 081101.
31. R. Jalem, Y. Tateyama, K. Takada and M. Nakayama, First-Principles DFT Study on Inverse Ruddlesden-Popper Tetragonal Compounds as Solid Electrolytes for All-Solid-State Li+-Ion Batteries, Chem. Mater., 2021, 33(15), 5859–5871.
32. E. Ziambaras and E. Schröder, Theory for structure and bulk modulus determination, Phys. Rev. B: Condens. Matter Mater. Phys., 2003, 68, 064112.
33. F. Birch, Finite Elastic Strain of Cubic Crystals, Phys. Rev., 1947, 71, 809–824.
34. G. J. Martyna, D. J. Tobias and M. L. Klein, Constant pressure molecular dynamics algorithms, J. Chem. Phys., 1994, 101, 4177–4189.
35. S. Nosé, A unified formulation of the constant temperature molecular dynamics methods, J. Chem. Phys., 1984, 81, 511–519.
36. W. G. Hoover, Canonical dynamics: Equilibrium phase-space distributions, Phys. Rev. A, 1985, 31, 1695–1697.
37. A. Hayashi, K. Noi, N. Tanibata, M. Nagao and M. Tatsumisago, High sodium ion conductivity of glass–ceramic electrolytes with cubic Na3PS4, J. Power Sources, 2014, 258, 420–423.
38. G. Bucci, T. Swamy, Y.-M. Chiang and W. C. Carter, Modeling of internal mechanical failure of all-solid-state batteries during electrochemical cycling, and implications for battery design, J. Mater. Chem. A, 2017, 5, 19422–19430.
39. A. K. Sagotra, D. Chu and C. Cazorla, Influence of lattice dynamics on lithium-ion conductivity: A first-principles study, Phys. Rev. Mater., 2019, 3, 035405.
40. M. T. Agne, T. Böger, T. Bernges and W. G. Zeier, Importance of Thermal Transport for the Design of Solid-State Battery Materials, PRX Energy, 2022, 1, 031002.
41. D. A. Weber, A. Senyshyn, K. S. Weldert, S. Wenzel, W. Zhang, R. Kaiser, S. Berendts, J. Janek and W. G. Zeier, Structural Insights and 3D Diffusion Pathways within the Lithium Superionic Conductor Li10GeP2S12, Chem. Mater., 2016, 28, 5905–5915.
42. A. A. Hubaud, D. J. Schroeder, B. J. Ingram, J. S. Okasinski and J. T. Vaughey, Thermal expansion in the garnet-type solid electrolyte (Li7-xAlx/3)La3Zr2O12 as a function of Al content, J. Alloys Compd., 2015, 644, 804–807.
43. X. He, Y. Zhu and Y. Mo, Origin of fast ion diffusion in super-ionic conductors, Nat. Commun., 2017, 8, 15893.
44. J. A. Newnham, J. Kondek, J. Hartel, C. Rosenbach, C. Li, V. Faka, L. Gronych, D. Glikman, F. Schreiner, D. D. Wind, B. Braunschweig, M. R. Hansen and W. G. Zeier, Correlation between the Coherence Length and Ionic Conductivity in LiNbOCl₄ via the Anion Stoichiometry, Chem. Mater., 2025, 37, 4130–4144.
45. B. Singh, Y. Wang, J. Liu, J. D. Bazak, A. Shyamsunder and L. F. Nazar, Critical Role of Framework Flexibility and Disorder in Driving High Ionic Conductivity in LiNbOCl ₄, J. Am. Chem. Soc., 2024, 146, 17158–17169.
46. Z. Zhang, H. Li, K. Kaup, L. Zhou, P. N. Roy and L. F. Nazar, Targeting Superionic Conductivity by Turning on Anion Rotation at Room Temperature in Fast Ion Conductors, Matter, 2020, 2, 1667–1684.
47. Y. Sun, Y. Wang, X. Liang, Y. Xia, L. Peng, H. Jia, H. Li, L. Bai, J. Feng, H. Jiang and J. Xie, Rotational Cluster Anion Enabling Superionic Conductivity in Sodium-Rich Antiperovskite Na3OBH4, J. Am. Chem. Soc., 2019, 141, 5640–5644.
48. K. Sau, S. Takagi, T. Ikeshoji, K. Kisu, R. Sato, E. C. dos Santos, H. Li, R. Mohtadi and S. I. Orimo, Unlocking the Secrets of Ideal Fast Ion Conductors for All-Solid-State Batteries, Springer Nature, 2024, preprint,  DOI:10.1038/s43246-024-00550-z.
49. R. Jalem, A. Hayashi, F. Tsuji, A. Sakuda and Y. Tateyama, First-Principles Calculation Study of Na+Superionic Conduction Mechanism in W- And Mo-Doped Na3SbS4Solid Electrolytes, Chem. Mater., 2020, 32, 8373–8381.

---