Accelerating cathode material discovery through ab initio random structure searching

The availability of reliable, safe, and accessible energy storage solutions is crucial for a world transitioning toward renewable energy sources from fossil fuels. Li-ion batteries (LIBs), initially developed in the 1990s, have become the leading answer to these challenges. While the manufacturing costs of LIBs have been reducing over the years, the reliance on transition metals that are expensive and sensitive to supply chain issues must be overcome for large-scale applications, such as electric vehicles (EVs) and grid-scale storage. Developing the next generation cathode material to meet these requirements is a multifaceted challenge,^1,2 which involves both understanding and optimizing existing materials as well as discovering new ones. The three main classes of commercialized cathode materials, layered LiCoO₂,³ spinel LiMn₂O₄,⁴ and olivine LiFePO₄,⁵ were all proposed three decades ago, highlighting the slow nature of cathode material identification and research. Current high-performance lithium cathode materials under intensive research include Li-rich NMC,⁶ disordered rock salts,⁷ and ε-VOPO₄.⁸ 

The explosive increase in computational power in the last two decades has made computational materials research a valuable tool.⁹ First-principles calculations can give invaluable insights for various properties of cathode materials,¹⁰ such as voltages,^11,12 Li diffusion barriers,¹³ disordering,⁷ and anionic redox.¹⁴ However, the crystal structure of the material of the interest must be known in the first place. The combination of standardized density functional theory (DFT) calculation routines and online databases has made it possible to perform the so-called “high-throughput” screening studies that aim to locate new materials based on data mining of hundreds and thousands of known materials.^15–17 On the other hand, it is an open question whether such an approach can truly find “new” materials beyond those generated by simple elemental substitution. In the meantime, in other fields, such as high-pressure research where experimental data are not as plentiful, first-principles calculations have become an invaluable tool for predicting the crystal structure of unknown materials with little or no experimental data.¹⁸ To predict stable crystal structures, the global minimum in a high dimensional potential energy surface (PES) must be located. Various algorithms have been developed for tackling this problem by exploiting intrinsic features of the PES. These methods include basin hopping,¹⁹ minima hopping,²⁰ genetic algorithms,²¹ particle swarm optimization,²² and random structure searching.^23,24 Several new data-driven approaches involving machine learning techniques have also been proposed recently, attempting either to learn the configuration spaces for generating sensible candidate structures^25,26 or to use machine learning potentials as surrogates for otherwise expensive first-principles calculations.^27–29 

In this article, we use several old and new cathode materials: LiFePO₄, LiFeSO₄F, and Fe-containing oxalates as examples to show that ab initio random structure search (AIRSS)^23,24 is a simple yet efficient tool for exploring the configuration space of complex materials and predicting the structure of existing and novel complex cathodes. We focus on Fe in this work as it is non-toxic, very abundant, and low-cost, making Fe-based LIBs viable for EV and grid-scale energy storage applications. Details of the working principles of AIRSS can be found in the literature.²⁴ In short, this method generates random but physically sensible candidate structures, followed by geometry optimizations using density functional theory calculations. This process is repeated until the collection of low energy structures has been repeatedly found. Compared to other structure search approaches, the simplicity and the lack of iterative improvement processes ensure that AIRSS is highly parallel. While first-principles calculations are usually thought to be expensive, the actual computational cost used for search can be significantly reduced, thanks to the crystal symmetry being kept during geometry optimization and the use of less strict convergence settings. As a result, a short time-to-result can be achieved.

To illustrate the utility of AIRSS in the field of cathode simulations, we first show that the two known experimental structures of quaternary LiFePO₄ can be located easily with AIRSS. Second, for quinary LiFeSO₄F, our search has located several experimentally observed polymorphs, together with previously unknown polymorphs that are predicted to be lower in energy than the known tavorite and triplite phases. One of these new polymorphs is predicted to possess both the high voltage found in the triplite phase and the three-dimensional Li diffusion network featured by the tavorite phase. Finally, we demonstrate how AIRSS can be successfully applied to determine new cathode materials with an analysis of fluorinated novel Li-stuffed Fe oxalates that have improved specific capacity and stability compared to the existing material Li₂Fe(C₂O₄)₂.^30,31

Ab initio random structure searching (AIRSS)^23,24 was used for locating low energy structures of a given composition. This method involves generating random but physically sensible structures, followed by geometry optimizations using first-principles calculations. The random structures generated are constrained by species-wise minimum separations (Table S1) and include up to four randomly chosen symmetry operations. The polyanions are introduced as rigid units when generating the structures, but they are allowed to relax in subsequent geometry optimizations. The plane wave density functional theory code CASTEP³ was used for geometry optimizations. During the search, core-corrected ultrasoft pseudopotentials^32,33 are on-the-fly generated as defined in the built-in QC5 library. A plane wave cutoff energy of 300 eV is used, and the reciprocal space sampling is performed using Monkhorst–Pack grids with a maximum spacing of 0.07 2π Å⁻¹. VASP^34–37 was used for further relaxations and property calculations with a plane wave cutoff energy of 520 eV; gamma-centered Monkhorst–Pack grids were used with a maximum spacing of 0.04 2π Å⁻¹, with projector augmented wave³⁸ potential dataset Perdew–Burke–Ernzerhof (PBE) v.54 (potentials used are Li_sv, Fe_pv, O, S, Mn_pv, and C). The Phonopy package is used for the calculations of phonon dispersions through the finite-displacement method.³⁹ First-principles calculations are performed using the PBE exchange–correlation functional⁴⁰ unless otherwise stated in the text. We also performed calculations using the PBEsol functional⁴¹ to validate the energy rankings among different polymorphs since it gives lattice constants closer to the experimental values. Hubbard U corrections⁴² are used with U_eff = 4.0 eV for Fe d electrons.⁴³ For Li diffusion analysis, ab initio molecular dynamics (AIMD) simulations and climbing-image nudged-elastic band^44,45 calculations were performed. The MD simulations were carried out at an elevated temperature of 1200 K to enhance the sampling of the Li distribution. The temperature was controlled by coupling the system to a Nose–Hoover thermobath.^46,47 For nudged-elastic band (NEB) calculations, a single Li vacancy is created in a supercell of √2 × √2 × 2 for LiFeSO₄F, containing 128 atoms in total. This ensures the minimum distance from an atom to its periodic image to be larger than 10 Å. The force convergence criterion was set to be 0.02 eV Å⁻¹. Additional electrons have been added to avoid complications arising from charge localization in the vacancy containing cells. Since the barrier height is primarily affected by the local structure, this approximation is expected to give minor effects.⁴⁸ For oxalates, the bandgap values are calculated using the HSE06 hybrid functional based on the PBE+U relaxed structure. The AiiDA framework^49,50 was used for managing the calculations and preserving calculation provenance. The pymatgen,⁵¹ ase,⁵² and sumo⁵³ Python packages are used for manipulating structure and general analysis. The VESTA⁵⁴ software is used for visualizing crystal structures.

Cathode materials are known for their chemical complexity—they typically involve at least three elements. The increasing number of elements results in a combinatorial increase in the configuration space to be explored, making structure prediction a challenging task. Here, we show that with the application of a few physical constraints, as described in the section titled Methods, the ground state structure of the known cathode materials can be obtained efficiently using AIRSS, together with many other known or unknown polymorphs.

LiFePO₄ (LFP) is a well-studied cathode material that has been commercialized.⁵ It was originally identified in the late 1990s and is now found in LIBs for electric vehicles. It has an olivine structure with space group Pnma. The Fe ions are octahedrally coordinated, and FeO₆ octahedra are corner-sharing and arranged in a zig-zag fashion. The PO₄ groups share edges with the FeO₆ octahedra. The primitive unit cell contains four formula units, giving 28 atoms in total. Our search found both the experimentally known olivine phase⁵ and the high-pressure Cmcm phase,^55,56 shown in Figs. 1(a) and 1(b), respectively. Figure 2 shows the energy differences and volumes of the relaxed phases with energy differences less than 50 meV per atom as compared to the ground state olivine phase. The high-pressure Cmcm phase is 13 meV per atom higher in energy, with a smaller volume. Structures from the Materials Project¹⁵ database are depicted in the plot as diamonds (MP), and those from the AIRSS search are marked by dots. All structures are assumed to have ferromagnetic (FM) spin arrangements. The low-temperature ground state of LFP has anti-ferromagnetic (AFM) spin arrangements.⁵⁷ However, the energy differences between FM and AFM spin configurations are found to be only in the order of 1–2 meV for the two experimental phases. Such a small energy change would have a negligible effect on the ranking of the lower energy structures.

In addition to the two experimental phases, we also found other low energy polymorphs, and many of them are not present in the Materials Project database, as shown in Fig. 2. To the best of our knowledge, there are no details available as to how structures in the latter were obtained, other than the experimental phases, which are based on the data from the Inorganic Crystal Structure Database (ICSD). There is a cluster of apparent low energy structures with 20% higher volume compared to the Pnma phase, as indicated by the red circle in Fig. 2. These structures are characterized by networks of tetrahedrally coordinated Li and Fe atoms. Our searches also found a few polymorphs with similar features, although the apparent lowest energy phase inside this cluster is not reproduced. Care should be taken when comparing the DFT energies of candidate structures with very different densities, especially with the PBE functional that is known to systematically overpredict the lattice constants of solids.⁵⁸ After re-relaxing these phases with the PBEsol functional,⁴¹ which is designed for solids and reproduces the lattice constants better, the cluster of less dense structures becomes ∼30 meV per atom higher in energy than the olivine phase (Fig. S1). This makes them less interesting as candidates for further investigations. We chose to use PBE initially as it was more widely used in the literature, but our results here show that PBEsol is a better choice for future studies involving ranking different polymorphs. Zhu et al. have previously reported finding the olivine phase using a structure prediction approach that involved motif based random structure generation with symmetry,⁵⁹ which is similar to our approach here. Rather than performing DFT geometry optimization, they included an additional step of constructing embedded atom model (EAM) potentials and use them for pre-screening, followed by final DFT relaxation on selected phases. While the fitted EAM potentials can accelerate the geometry optimization, it is not clear whether they can faithfully reproduce features of the actual potential energy surface, which is essentially for biasing the search toward finding realistic low energy structures.

LiFeSO₄F is a high voltage iron-based polyanion cathode material giving similar capacity as compared to LiFePO₄ but having the advantage of being more ionically and electronically conductive,⁶⁰ which may remove the need for resorting to nanosizing and conductive coating of particles as in LiFePO₄ cathodes.^60–63 Recham et al.⁶⁰ successfully synthesized the tavorite phase using FeSO₄·H₂O and LiF as the precursors with an ionothermal method that operates at relatively low temperatures. They found the tavorite phase to be highly reversible with a voltage of 3.6 V against Li metal. Shortly afterward, a triplite phase of the same composition was reported to have a record-high voltage of 3.9 V for Fe-based polyanion cathodes.^64,65 Unlike the tavorite phase, which has well-defined Li and Fe sites,⁶⁶ the triplite phase has Li–Fe site occupancy disorder, giving an entropic stabilization. This is consistent with it being synthesized by annealing the tavorite phase,⁶⁴ and later, it was shown that spark-plasma synthesis and a high-temperature solid-state route can produce it directly from the precursors.^67,68 Other polymorphs also exist for this class of material taking the general formula Li/NaMSO₄F (M = Co, Ni, Mn).^69–71 The large difference in operating voltage between the tavorite and triplite phases has attracted the attention of theoretical studies.^12,72–74 Chung et al.⁷³ performed first-principles calculations and showed that the difference in voltage arises from the stabilities of the delithiated phases, while the lithiated structures are very similar energetically. The work of Yahia et al.⁷⁴ pointed out that the high voltage of the triplite phase can be related to the cis arrangement of F⁻ ions, giving larger repulsion than the trans arrangement found in the tavorite phase. As the lithiated triplite and tavorite structures have similar stability, the strong F–F repulsion in the former is compensated by attractive Li–F interactions, which result in a less stable delithiated phase. The secondary inductive effect caused by neighboring Li atoms close to the O anions has been identified as an effective indicator of open circuit voltages among different Fe-based polyanionic cathodes.¹² The increased voltage of the triplite phase is attributed to having an increased number of Li neighbors around the Fe-centered octahedra when compared to that of the tavorite phase.

The rich set of structure–property relationships in polymorphs of LiFeSO₄F leads one to wonder if there are yet more phases to be discovered, potentially with even better electrochemical properties. The involvement of F further increases the complexity of the search, making it a quinary system. The two known experimental phases, tavorite and triplite, are displayed in Figs. 3(a) and 3(b), respectively. The relative energies and volume per formula unit of the structures found in the search are shown in Fig. 3(c). To mitigate the issue of PBE favoring less dense structures, as it was found in the search of LiFePO₄ above, we searched for ambient pressure as well as with a 10 GPa external pressure. There are two experimental tavorite structures with slightly different Li atom arrangements. The initial report of the tavorite phase proposed that the Li atoms occupy two half-occupied 2i Wycoff sites,⁶⁰ but later neutron diffraction data show that they occupy a single 2i site.⁶⁶ This structure is labeled Ta-I, and it has lower energy than the former in an ordered unit cell (Ta-II). Our search found two more variants that have slightly different Li sites but similar energies (phase F/F′). The triplite phase (space group C2/m) has been reported to exhibit Li–Fe site-occupation disordering.⁶⁴ To obtain a meaningful theoretical model, unique Li–Fe configurations in the primitive cell are enumerated by distributing four Fe and four Li atoms. These structures are labeled green crosses in Fig. 3(c). The lowest energy structure is shown in Fig. 3(a) (Tr-I), which has a cation arrangement consistent with the work of Chung et al.⁷³ The same structure is also obtained in the search.

New low energy polymorphs have been found by our computational search, and many of them have comparable or even lower energy compared to the experimental phases, as shown in Fig. 3(c). Structures A–D have FeO₄F₂ octahedra arranged in a configuration containing one-dimensional edge-sharing chains. Phase A [Fig. 4(a)] has the lowest energy. It has F ions arranged in a trans fashion in the FeO₄F₂ chains, which are separated by isolated SO₄ groups and tetrahedral coordinated Li atoms. The flexibility in the octahedra chains and separators results in different polymorphs with similar energies. The B phase [Fig. 4(b)] has similar features, but the FeO₄F₂ chains are arranged in a different pattern. Both have larger volumes compared to the known tavorite/triplite phases. On the other hand, the C phase [Figs. 5(a) and 5(b)] is more efficiently packed. In fact, it closely resembles the sillimanite Al₂SiO₅ structure.⁷⁵ A variant of this structure, phase D [Fig. 5(c)], has a different occupation of Li and S sites between the FeO₄F₂ chains. Previously, LiZn_xFe_1−xSO₄ has been found to crystalize in the sillimanite structure, where x can be as high as 0.15 using solid-state synthesis, or 0.1 using an ionothermal approach, before the tavorite/triplite phase starts to crystalize instead.⁶⁹ The reported structure of this phase has Li atoms in the octahedral site and Fe atoms in the tetrahedral sites, e.g., the opposite of phases C and D. To avoid confusion, from now on, we refer to the sillimanite structure found here as sillimanite-II and the previously reported LiZn_xFe_1−xSO₄ phase as sillimanite-I. The latter is also rediscovered in our search and labeled structure E [Fig. 5(e)].

A factor that was neglected in the initial search and the analysis so far is the spin configuration, as all sampled structures are assumed to be ferromagnetic. To assess the effect of spins, collinear magnetic configurations have been enumerated to find the low energy AFM states using the method implemented in the pymatgen package.⁷⁶ The energy differences (with respect to the FM state of phase A) for both AFM and FM state are shown in Fig. S2. The energy difference between the AFM and FM spin state for each phase is tabulated in Table I. Changing FM to AFM spin arrangement reduces the total energy, with the magnitude of the reduction typically about 1 meV, which is inconsistent with those found for LiFePO₄ above. Although in some cases, it can be up to 5 meV per atom. This can be understood as the potential energies are mainly contributed by the electrostatic energies and short-range repulsions in ionic crystals. As a result, assuming FM ordering in the initial search would have a relatively small effect on the overall ranking of the structures.

Since phases A–D all share similar features but differ in their volumes by up to 10%, it raises the question of whether the seemingly better stability of phase A, which has a larger volume, is an artifact of the PBE exchange–correlation functional. This was found to be the case for the polymorphs of LiFePO₄ as discussed above. Additional relaxations have been performed using the PBEsol functional, keeping the same FM spin configurations. The energy difference between phases A and C is now reduced by ∼5 meV per atom, and the latter becomes the lowest energy structure, as shown in Fig. S3. The change in energy here is much smaller than those for LiFePO₄ polymorphs. Unlike LiFePO₄, despite subtle changes in the rankings, the collection of low energy structures remains the same as those using the PBE functional. Further checks using the hybrid functional HSE06 show that phase C has lower energy compared to the known sillimanite-I structure (phase E), as well as the tavorite and triplite structures (Fig. S4).

Selected properties of the different polymorphs are tabulated in Table I. The average voltages of the tavorite phases are found to be 3.52 eV, which is in good agreement with other theoretical studies^72,73 and the reported value of 3.6 V measured from experiments.⁶⁰ We note that the voltages calculated are positively correlated with the choice of the U parameter in the generalized gradient approximation (GGA)+U. The FM spin configuration is used since the energy differences in AFM-FM spin configurations give an insignificant effect on the voltage. The initial atomic positions of the delithiated cells are shaken by an amplitude of 0.05 Å to break the symmetry. For the triplite phases, we found their voltages to range from 3.9 to 4.3 V depending on the Li–Fe orderings. The most stable Fe orderings are likely to be different in the lithiated and delithiated cells. These values are in consistent with the experimental reported 3.9 V.⁶⁴ For structure E (sillimanite-I), the average voltage is found to be 3.51 V, which is in good agreement with the experimental value of 3.6 V.⁶⁹ On the other hand, those in the sillimanite-II structure, phases C and D, have a higher voltage of 4.0 V. Their volumes are reduced by 5.7% and 4.9%, respectively, after delithiation, which are lower than that of the tavorite phase (Ta-I). On the other hand, the triplite structures have much smaller volume changes upon delithiation, which can be related to Li–Fe disordering and being more densely packed in the first place. The positive volume changes displayed in Table I are likely to be a result of only considering removing Li in specific Li–Fe orderings; hence, they are not truly representative of a truly disordered system.

We now turn our focus to the sillimanite-II structured C phase that has edge-sharing chains of FeO₄F₂ octahedra and separated by chains of corner-sharing SO₄ and LiO₃F tetrahedra. This phase has space group P2₁/c—a subgroup of Pnma in the original sillimanite Al₂SiO₅ phase. Its symmetry is further lowered to $P1̄$ when Fe atoms are in an anti-ferromagnetic spin arrangement. Finite-displacement phonon calculations show that it is dynamically stable since no imaginary frequencies are found across the first Brillouin zone, as shown in Fig. S5. Compared to the previously reported sillimanite-I structure (phase E), the sillimanite-II structured phase C has a higher voltage (4.0 V vs 3.5 V). The increased voltage means a higher energy density for cathode applications.

While the triplite phase also has a relatively high voltage, Li–Fe disordering limits the reversible extraction and insertion of Li, giving a lower reversible Li content as compared to the tavorite and the sillimanite-I phase.^60,64,69 In the sillimanite-I structure, the nearest neighbor Li atoms are along the c direction, and the second nearest neighbor Li atoms are much further apart as compared to the nearest neighbors (3.2 Å vs 5.2 Å). This implies that Li atoms are constrained to diffuse along 1D channels. In contrast, the first and second nearest neighbors of Li atoms are 4.3 and 4.5 Å away in sillimanite-II structured phase C. We have identified two potential Li-vacancy hopping pathways, as shown in Figs. 6(a) and 6(b) respectively. Path A involves the zig-zag movement of Li atoms along the [010] direction. To hop to the next site using a direct path, the Li atom must go through an edge and a face of the LiO₃F tetrahedra. The NEB calculations indicate that the transition state barrier height of this process is 0.48 eV, as shown in Fig. 6(c). The shape of the curve is a result of the Li atom taking an indirect route—it hops around the edge of the coordination tetrahedron rather than squeezing through it. The results of NEB calculations for path B are shown in Fig. 6(d). In this case, the overall mass transport takes place along the [001] direction, although individual Li atoms also take zig-zap paths. The barrier of a single vacancy hop is found to be 0.45 eV, and the symmetrical shape of the curve is a result of the initial and final images being related by the inversion symmetry. Because a single hop along path B would take the Li atom to the next chain of path A in the [101] direction, it is also possible to have net Li diffusion along the [100] direction. We also performed AIMD simulations and found that the distribution of Li forms a zig-zag 3D network (Fig. S6) in good agreement with the two NEB paths identified. The relatively low average barrier height, 0.45 eV, makes phase C a promising material, with good ionic conductivity. In comparison, the barrier in LiFePO₄ was found to be about 0.55 eV,⁷⁷ the tavorite LiFeSO₄F and triplite LiFeSO₄F were reported to have barriers of 0.57 and 0.6–0.8 eV, respectively,⁷⁸ and the barrier in Li₂FeSiO₄ was reported to be 0.91 eV.⁷⁹ We note that the exact values of the barrier height can be affected by the details of the methods used, i.e., the level of theory and treatment of vacancy defects.

The relatively small energy differences between sillimanite-I and sillimanite-II phases raise the question as to whether it is possible to have Li–Fe occupancy disordering. Enumeration of all symmetrically unique Li–Fe arrangement phases C–E shows deviations from the original ordering would result in the energy change of at least 0.15–0.2 eV per primitive cell (∼5 meV per atom), which is higher than that of the triplite phases (1.5 meV per atom). The results are shown in Fig. S7. While detailed studies of Li–Fe anti-sites in the sillimanite-II structure are beyond the scope of this work, such defects should not significantly impact the Li diffusivity due to the three-dimensional nature of the diffusion pathways.

Our search has uncovered a rich landscape of low energy polymorphs for LiFeSO₄F, and a new sillimanite-II phase is predicted to have both the high voltage of the triplite and fast Li diffusion kinetics of the tavorite. This leads to the question: Why has this phase not been synthesized yet? Solid-state synthesis is known to be limited by the slow reaction rate and limited mass transport. The most widely reported synthesis routes for LiFeSO₄F start from FeSO₄·H₂O and LiF, where the former is already in the tavorite structure, making the formation of tavorite LiFeSO₄F favored via topotactical replacement of O²⁻ in H₂O by F⁻ in conjunction with the intercalation of Li⁺ ions.⁶⁰ The triplite phase is generally synthesized from the tavorite LiFeSO₄F by annealing or directly from the anhydrous FeSO₄ and LiF using spark-plasma synthesis or ball milling.⁶⁷ This can be rationalized by the entropic contributions from the Li–Fe disordering, making the triplite phase favored at high temperature. In addition, the sillimanite-I phase also has structural similarities to the tavorite phase, as it also has the corner-sharing Fe(Zn)O₄F₂ octahedra.⁶⁹ The sillimanite-II phase, on the other hand, contains the distinct feature of edge-sharing FeO₄F₂ chains not found in the other polymorphs, and the small energy difference compared with the tavorite makes the thermodynamic driving force relatively low for its formation. Nevertheless, we note that the sillimanite-II phase is in fact structurally closely related to the anhydrous FeSO₄. Future synthesis works may attempt to exploit this similarity, although it will still be challenging to avoid forming the triplite phase if high-temperature synthesis is required.

The commercial success of lithium–iron–phosphate cathodes has triggered extensive research interest in searching for polyanion compounds as cathode materials for batteries due to their low cost, long cyclic life, and high safety.⁸⁰ To date, a wide spectrum of polyanions have been intensively studied, including phosphate (PO₄)³⁻, sulfate (SO₄)²⁻, and silicate (SiO₄)²⁻.^81–84 Despite the significant effort, another family of polyanion compounds have been largely overlooked, i.e., the oxalates (C₂O₄)²⁻. In fact, to the best of our knowledge, only two iron-based oxalates, i.e., Fe₂(C₂O₄)₃·4H₂O⁸⁵ and Li₂Fe(C₂O4)₂,³⁰ have been reported as cathodes for LIBs, while the electrochemical performance of other transition-metal oxalates is largely missing. Interestingly, if one considers the polarizability of the oxalate group, it is comparable to that of (PO₄)³⁻. Therefore, these anions should exert a strong inductive effect on transition metals, thereby providing competitive redox potentials during battery discharge.

The oxalates are also characterized by a few other advantageous features when used as cathodes for batteries. One of the most obvious is the ease of synthesis. Due to the varied solubility of oxalic acid and oxalate metal salts, oxalate-based compounds can be synthesized under relatively mild conditions (usually below 200 °C) via solution-based methods such as hydrothermal synthesis. This exerts advantages over conventional inorganic polyanion compounds such as LiFePO₄ that requires solid-state sintering at elevated temperatures of up to 700 °C. In fact, it is well-known in synthetic and structural chemistry that the oxalates can crystallize into a wide variety of different polymorphs by serving as monodentate, bidentate, tridentate, or tetradentate ligands. By complexing with different metal centers, the oxalates could form one-dimensional chain-like, two-dimensional layered, and three-dimensional connected frameworks.³⁰ Such polymorphism provides yet another exciting opportunity to tune the structural stability and the electrochemical performance of the material. In particular, such structural features have a profound influence on cation mobility and can be used to enhance the rate capability of the cathode. For example, in a recent contribution, it is shown that by reducing the dimensionality of the host framework, superionicity can be achieved in lithium-rich anti-perovskites.^86–89 Beyond structural considerations, the polyanionic nature and the covalent characteristics of the C–O bonds in oxalate anions offer exciting opportunities to tune the electrochemistry as well. Recently, anionic redox (mostly oxygen) has been discovered as an important source of additional capacity beyond the transition-metal contribution during battery discharge.⁹⁰ Unfortunately, in conventional oxide materials such as Li(Li_xNi_{1−x−y−z}Mn_yCo_z)O₂ (Li-rich NMC), the stripping of electrons from the oxygen non-bonding states leads to cyclic instability and sometimes oxygen release. In this regard, polyanionic species that bear covalently bonded oxygen, in principle, should have better anti-O₂ release capability and better reversibility in anionic redox.

All these features make oxalates a valuable family of cathode materials to explore. Nevertheless, they are less studied as a cathode material and only started to gain attention very recently. Ahouari et al.⁸⁵ first looked at the redox and structural evolution of a commercially viable Fe (III) oxalate, e.g., Fe₂(C₂O₄)₃·4H₂O, during lithiation and showed that these compounds demonstrate an average voltage of 3.35 V vs Li/Li⁺. However, due to the high concentration of coordinated H₂O in such materials, the gravimetric capacity is relatively low and only reaches 98 mA h g⁻¹. Moreover, the study was restricted to relatively narrow voltage cutoffs, and the potential capacity contribution from oxygen was not looked at. Very recently, Yao et al. synthesized a family of Li-stuffed Fe (II) oxalates with a composition of Li₂Fe(C₂O₄)₂.³⁰ Interestingly, they showed that more than one Li can be taken out of the structure reversibly and anionic redox is incorporated during the process. We recently carried out a thorough search of polymorphs on the composition of Li₂Fe(C₂O₄)₂ and found several relatively stable phases that are potentially synthesizable.³¹ Despite these pioneering efforts, oxalate cathodes are still in their infancy. Many more compositions are yet to be looked at. The additional chemical space and configuration space that may yield possible polymorphs with excellent electrochemical performances have never been explored as well. In this contribution, we take advantage of the efficient AIRSS algorithm and conduct searches with two compositions that have not been reported before: LiFeC₂O₄F and Li₂FeC₂O₄F₂. The main motivation for incorporating with F⁻, other than expanding the chemical space, is to increase the theoretical capacity since F⁻ has a lower mass-to-charge ratio (19 e⁻/u) compared with C₂O₄ (44 e⁻/u). Fluoroxalate materials have been proposed for cathodes of potassium-ion batteries and found to have exceptional cyclability.⁹¹ By placing a special focus on the Fe (II) containing compounds, we look for cheap, high-rate, easy-to-synthesize, and reasonably energy-dense cathodes. Specifically, to keep the problem computationally tractable, we look at both compositions with the number of formula units no larger than four. The C₂O₄ groups are fixed as units for the sole purpose of generating the initial random structures, which are constrained to have two to four symmetry operations.

In Fig. 7, we show the lowest energy structures obtained for Li₂FeC₂O₄F₂ and LiFeC₂O₄F. Li₂FeC₂O₄F₂ has a theoretical specific capacity of 274 mA h/g, assuming that all Li can be extracted, which is higher than that of the previous report Li₂Fe(C₂O₄)₂ (218 mA h/g). Both Li₂FeC₂O₄F₂ and LiFeC₂O₄F contain infinitely extending chains formed by six-fold coordinated Fe²⁺ where neighboring octahedra are connected by the oxalate group. The Li atoms are tetrahedrally coordinated by O and F. The corner-sharing Li tetrahedra form chains that are parallel to those made of Fe centered octahedra. The dynamic stabilities of the predicted structures are confirmed by finite-displacement phonon calculations, which have found no imaginary modes across the first Brillouin zone for both (Fig. S8). Our calculation shows that Li₂FeC₂O₄F₂ is metastable with a distance-to-hull of 55 meV per atom, which is slightly higher than that of the previously reported Li₂Fe(C₂O₄)₂ phase, which is 51 meV per atom above the hull. The metastability of these materials can be attributed to the (C₂O₄)²⁻ group. Even Li₂C₂O₄, which is known experimentally, is also predicted to have an energy above the convex hull (47 meV per atom). The existence of these metastable phases can be understood as breaking the C=C bond requires overcoming a large energy barrier. On the other hand, LiFeC₂O₄F is more unstable with a higher distance-to-hull of 83 meV per atom. Since the reaction LiFeC₂O₄F + LiF → Li₂FeC₂O₄F₂ has negative reaction energy, it is unlikely to be synthesized. Therefore, we focus only on Li₂FeC₂O₄F₂ from now on.

A cathode material must be both electronically and ionically conductive. Using the hybrid HSE06 functional, the bandgap of Li₂FeC₂O₄F₂ is computed to be about 2.8 eV, which is lower than other polyanion materials such as LiFePO₄ and is comparable to materials known to have decent electronic conductivities, for example, LiCoO₂.⁹² For ionic conduction, two Li diffusion pathways can be identified, one going along the Li-chains (paths 1 and 2) and the other across them (paths 3 and 4), as shown in Figs. 8(a) and 8(b), respectively. Both path 1 and path 2 have relatively low barriers of 0.39 eV with Li atoms moving through faces of the coordination tetrahedron. On the other hand, path 3 and path 4 have much higher barrier heights (0.82 and 1.31 eV), which can be attributed to long hopping distances and the need for the Li atom to go through tetrahedral edges. Hence, the movement of Li is likely to be constrained inside the 1D channels along the a direction shown in Fig. 8(a).

The stabilities of delithiated phases of Li_2−xFeC₂O₄F₂ (0 ≤ x < 2) are investigated by enumerating unique configurations of the lithiated phase with Li atoms removed in 2 × 1 × 1 and 1 × 2 × 1 supercells (four formula units). Each structure generated is subsequently relaxed following the specification outlined in the section titled Methods. Convex hull constructions show that the only stable composition in between the terminal compositions is LiFeC₂O₄F₂, resulting in two plateaus in the voltage profile, each originating from a two-phase reaction [Fig. 9(a)]. The voltage of the first step (0 < x < 1) is found to be 3.52 V, and the volume of the delithiated structure changes by ∼1%. This gives a theoretical energy density of 482 W h kg⁻¹, which is higher than the 425 W h kg⁻¹ of Li₂Fe(C₂O₄)₂ (TM redox only, 3.9 V). Further delithiation requires a much higher voltage of 5.3 V. In comparison, removing more than one Li from Li₂Fe(C₂O₄)₂ was reported to require only a minor increase in the voltage (up to ∼4.2 V),³⁰ and the extra capacity has been attributed to anionic redox. In contrast, no evidence of anionic redox is found for Li₂FeC₂O₄F₂ here based on GGA+U calculations. The Fe-projected magnetizations in delithiated structures are shown in Fig. 10(a). Within the range of 0 < x < 1, the increase in μ_b is consistent with Fe²⁺/Fe³⁺ redox, and the reversal beyond x > 1 indicates that the Fe³⁺/Fe⁴⁺ redox is active. The same trend is also found for the total magnetization of the entire cell. The inactivity of (C₂O₄)²⁻ is confirmed by tracking the C=C bond lengths, which stay almost unchanged as more Li atoms are taken away. In contrast, removing more than one Li per Fe in Li₂Fe(C₂O₄)₂ leads to spontaneous breaking of a (C₂O₄)²⁻ group into two CO₂-like parts, giving rise to increased average C–C distances, as shown in Fig. 10(b). The superior C=C bond stability in Li₂FeC₂O₄F₂ is further supported by the Fukui functions calculated for delithiated structures with x = 1, where the electron density of the C=C bond has a much smaller contribution compared to that in Li₂Fe(C₂O₄)₂, shown in Fig. S9. Similar results are obtained using the HSE06 hybrid functional, where the (C₂O₄)²⁻ groups only spontaneously break at very high delithiation levels (Fig. S10).

Our prediction of Li₂FeC₂O₄F₂ opens a new avenue for developing transition metal oxalate-based cathode materials and demonstrates the effectiveness of AIRSS in exploring cathode materials with new compositions. The incorporation of the F⁻ ion is shown to improve the specific capacities and stabilize the oxalate group from dissociation at high levels of delithiation. Despite no anionic redox being observed, our predictions provide a model system for future work to unravel the nature of the anionic redox in polyanions.

We have shown that ab initio random structure searching (AIRSS) can be used to predict the low energy polymorphs of complex cathode materials efficiently. For LiFePO₄, the two experimental phases are rediscovered in the search, along with many other unknown low energy polymorphs. Our searches for LiFeSO₄F have not only rediscovered several experimental phases but also found new sillimanite structured polymorphs with a sillimanite structure with edge-sharing FeO₄F₂ chains. One of these new phases is predicted to have both a higher voltage (∼4.0 V) and a 3D Li diffusion network with relatively low barrier heights, combining the advantage of both the existing tavorite and triplite phases. When applied to fluorinated Li–Fe oxalates, our search has found a Li₂FeC₂O₄F₂ phase with 1D Li diffusion pathways and an average voltage of 3.5 V utilizing the Fe²⁺/³⁺ redox couple. The incorporation of F⁻ increases the specific capacity and mitigates the decomposition of (C₂O₄)²⁻ groups at high delithiation levels, resulting in improved structural stability during cycling.

First-principles structure predictions are relatively computationally demanding. Given the eventual cubic scaling of plane wave DFT, the computational cost can increase significantly with the number of atoms in the unit cell, not to mention the fact that more structures will need to be sampled. Nevertheless, we found that systems up to 30–40 atoms can be addressed at moderate costs, and many inorganic crystalline materials have unit cells that fall within this range. While making predictions from existing databases though substitution is relatively “cost-free,” such an approach is not applicable to underexplored chemical spaces. Existing structural motifs and design rules^93,94 can be utilized to choose chemical systems with the potential to be high-performance cathode materials, followed by explorative searches to find synthesizable phases in targeted chemical spaces.

Cathode materials often exhibit site-occupation disordering because of similarities in ionic radii. It is difficult to addresses disordered phases directly by DFT, and therefore, finding them directly via an AIRSS search is unlikely; however, they can still be identified through the existence of an ensemble of relaxed structures with similar energies sharing common frameworks. Subsequently, methods as a special quasi-random structure⁹⁵ and cluster expansion^96–99 can be applied to test if the system is truly disordered.

Our approach is not limited to Li-ion intercalation cathodes, and crystal structures of Na/K containing cathodes can be searched similarly. Admittedly, finding new phases theoretically, i.e., demonstrating the existence of a low energy local minimum, does not guarantee that the phase can be realized, as kinetic barriers play a very important role in synthesis. Nevertheless, knowing the existence of possible new phases will undoubtedly help to inspire and guide future experimental works. We hope that by combining predictive computational approaches and experimental efforts, the discovery of new novel cathode materials will be greatly accelerated.

See the supplementary material for the details of the LiFePO4 search, tabulated species-wise separations, and additional calculation results.

This work was supported by the Faraday Institution (Grant No. FIRG017) and used the MICHAEL computing facilities. C.J.P. acknowledges support from the EPSRC through the UKCP consortium (Grant No. EP/P022596/1). D.O.S. acknowledges support from the European Research Council (Grant No. 758345). Through our membership of the UK’s HEC Materials Chemistry Consortium, which is funded by EPSRC (Grant Nos. EP/L000202, EP/R029431, and EP/T022213), this work used the ARCHER UK National Supercomputing Service (http://www.archer.ac.uk), the UK Materials and Molecular Modeling Hub that is partially funded by EPSRC (Grant Nos. EP/P020194 and EP/T022213). The authors also acknowledge the use of the UCL Myriad and Kathleen High Performance Computing Facility (Myriad@UCL, Kathleen@UCL), and associated support services, in the completion of this work.

The authors have no conflicts to disclose.

B.Z. and Z.L. contributed equally to this work.

The data that support the findings of this study are openly available in Zenodo at http://doi.org/10.5281/zenodo.5572347.