Ab initio random structure searching for battery cathode materials

Cathodes are critical components of rechargeable batteries. Conventionally, the search for cathode materials relies on experimental trial-and-error and a traversing of existing computational/experimental databases. While these methods have led to the discovery of several commercially viable cathode materials, the chemical space explored so far is limited and many phases will have been overlooked, in particular, those that are metastable. We describe a computational framework for battery cathode exploration based on ab initio random structure searching (AIRSS), an approach that samples local minima on the potential energy surface to identify new crystal structures. We show that by delimiting the search space using a number of constraints, including chemically aware minimum interatomic separations, cell volumes


I. INTRODUCTION
Lithium-ion batteries (LIBs) are electrochemical energy storage devices characterized by high energy density, good rate capability, and a long shelf-life. They have dominated the market for portable electronics since their commercialization in the 1990s. 1,2 This technology has been challenged by the increasing demands for higher energy density, due to the growing electric vehicle market, and grid-scale energy storage. One of the most important performance parameters, the energy density of a battery, is determined by its average discharge voltage, lithium storage capacity, and weight (or volume). Among these factors, both the voltage and the capacity depend heavily on the physiochemical properties of cathode materials used. Therefore, the search for cathode materials with high discharge voltages and large capacities is one of the central aims of battery research. 3,4 Beyond that, due to the demand for highpower output by electric vehicles and grid-scale energy storage, a high-rate capability of a battery is also essential. In this context, The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp cathodes that transport both lithium and electrons rapidly are critically needed. Conventionally, new cathode materials are discovered through ad hoc design based on the understanding of known chemistry and crystallography. For example, the intercalation reaction between guest ions with solid hosts has been known for over a century. 4 By studying the fundamental properties of the LixTiS 2 series, Whittingham used the layered TiS 2 as a cathode and demonstrated the first rechargeable lithium battery. 5 However, such cathodes suffer from low redox potential vs Li/Li + , so the energy density must be improved. To achieve this, Goodenough and co-workers replaced S 2− with O 2− and enabled the use of the Co 3+/4+ redox couple based on the knowledge that the S 2− 3p band lies at a higher energy as compared with that of O 2− 2p. 6 This led directly to the discovery of LiCoO 2 , one of the most successful cathode materials for LIBs to date. 7 Despite the commercial success, the search for cathode materials based on ad hoc design and experimental trial-and-error has proven difficult. In fact, only a very limited number of materials have been confirmed to be competitive, despite the tremendous efforts made in searching for them. A new paradigm for cathode material discovery is critically needed.
Computational materials science has become a powerful tool in materials discovery. [8][9][10] Thanks to developments in electronic structure methods, in particular, density functional theory (DFT) and related computational approaches, the first-principles computation of the properties of practical materials has become a routine in materials research. [11][12][13][14] At the same time, the dramatic increase in computational capacity of our research infrastructure has empowered researchers to carry out high-throughput screening of materials for those with desired properties. 8 Databases have been assembled to archive properties such as formation energies and bandgaps. 15,16 In the context of cathode materials, phase stability, discharge voltage, and cation diffusion properties can now be calculated and analyzed based on these precomputed data entries. 17,18 Pioneering efforts have also been made to predict new cathode materials. For example, elemental substitution 19,20 and motif templating [21][22][23][24] have been proposed as a means to effectively generate candidate structures, followed by property screening. While such approaches led to the discovery of a number of cathodes, 19,20,25 the data records typically are limited in number. The compositional and structural space explored is severely constrained, and many possible combinations of chemical species have been overlooked. More importantly, even for a single given composition, crucial information is omitted that would aid in the exploration of new materials, namely, polymorphism, or different low energy structures. These metastable materials can sometimes be synthesized and display desirable properties. In fact, materials with polymorphism have been widely studied in a number of areas including high pressure physics, organic chemistry, and electrochemistry. [26][27][28] In the context of battery cathode materials, such polymorphism may provide exciting opportunities for materials with extraordinary electron and ion transport characteristics.
Here, we address these issues by showing how ab initio random structure searching (AIRSS) can be used as an efficient tool for the exploration of novel cathode materials for batteries. [29][30][31] The method is based on a random sampling of the first principles potential energy surfaces (PESs) of atomic arrangements and has been successfully applied in a number of areas including high pressure physics, superconductors, semiconductors, and energy storage. 30,[32][33][34][35] Here, we demonstrate that by constraining the search space through the choice of parameters such as chemically aware minimum interatomic distances, cell volumes, and symmetries, the efficient identification of both stable and metastable cathode materials is possible. Specifically, we explore existing cathode compositions such as LiCoO 2 , LiFePO 4 , and LixCuyFz to showcase this method and show how to achieve an efficient sampling of the PES. Building on this, we apply the AIRSS method to lithium-stuffed transition-metal oxalates in the search for novel cathode materials. Screening the low energy outcomes of the search, we identify several oxalate polymorphs as good candidates with decent rate capability and energy density.

II. METHODOLOGY
A. General searching framework The general searching framework is illustrated in Fig. 1. It consists of the choice of composition, the crystal structure search using AIRSS, and the post-search screening. Composition-wise, we chose LiCoO 2 , LiFePO 4 , LixCuyFz, and Li 2 TM(C 2 O 4 ) 2 (TM = Fe, Co, Ni, V, and Mn) in this study. LiCoO 2 and LiFePO 4 were selected because they are the most representative materials of two cathode families: simple transition metal oxides and polyanionic compounds. It is worth noting that, to determine the stability of a phase, not only the composition of interest needs to be considered but also all the other possible ones within the given chemical space. While for most systems, the stable phases can be directly taken from databases, for systems that are less studied, an extensive search over the entire compositional space is necessary to construct the relevant phase diagram. In this case, we demonstrate the search of a phase diagram using LixCuyFz because Cu 2 F is a high-capacity conversion-type cathode, which undergoes interesting conversion reactions during charge and discharge. Finally, Li 2 TM(C 2 O 4 ) 2 was chosen to highlight the capability of the current searching scheme because Li 2 Fe(C 2 O 4 ) has been reported to be a potential cathode with interesting anion redox properties while all the other members are less well studied.
After choosing the composition, crystal structure predictions were carried out using AIRSS. The details of this approach are provided in Subsection II B, and the impact of selecting the few search parameters on its efficiency is discussed in Sec. III A. After the AIRSS search, candidate structures were selected for property screening based on their phase stability. For a cathode material, the properties of interest are the voltage, the capacity, the cation transport barrier, and the electronic bandgaps. The first two determine the energy density of a cathode, while the latter two provide a proxy to predict the rate capability. These property calculations are standard in battery research and are briefly discussed in Sec. II C.

B. Ab initio random structure searching
At the core of the searching scheme, AIRSS samples the local energy minima of a PES by first generating random "sensible" structures and then locally optimizing their geometries using firstprinciples calculations. The efficiency or even the success of the algorithm, as we will explain later, is highly dependent on choosing sensible parameters and constraints based on physical quantities, e.g., reasonable density, no close contacts, and proper neighbor relations. The general workflow of an AIRSS search is depicted in Fig. 1(a). Generating random structures with few constraints will provide the widest coverage of the search space and therefore reduce the possibility of completely missing out the global minimum. However, it will also result in low search efficiency, i.e., the rate of encountering lower energy structures is reduced. Considering that cathode materials operate at ambient pressures, it is natural to use chemical ideas to constrain the search space and steer the search toward chemically "sensible" regions. This is done by applying constraints during structure generation. To generate a structure that is chemically sensible, both the bond lengths and the neighboring relation between atoms need to be constrained. Two parameters are effective in this context. One is the minimum separation between atoms (which may be different for each pair of species and so chemically aware). To acquire structures that fulfill the distance constraints, they are generated randomly and are optimized using a hard sphere potential to reach the desired species-dependent minimum separations. It is important to note that such minimum separations are species-pairwise and an example is shown in Table I. Such distance constraints avoid close contacts between atoms, which may lead to numerical problems for the first principles code during the geometry relaxation. Furthermore, by fixing the cell volume (per formula unit/f.u.), such species-pairwise minimum separations can be used to avoid unwanted neighboring relations. For example, by setting the minimum separations between Li and Li and between Li and O to 2.8 and 2.0 Å, respectively, in a lithium containing oxide, OLix clusters will be favored instead of Lix neighboring. Should one have knowledge on the chemical system, the values of the minimum separations and cell volumes can be measured from known structures. For the case where the chemistry is new, one can perform a preliminary AIRSS search using random minimum separations and then measure the resulting values from the structure with the lowest energy. In practice, to avoid overly constraining the system, the structural rejection criteria based on the minimum separation matrix are slightly loosened. Symmetry is also exploited during the AIRSS search as most known crystals are characterized by a degree of symmetry. Moreover, adopting symmetry accelerates the DFT calculations through reducing the number of k-points required to sample the Brillouin zone and the number of degrees of freedom in the local geometry optimization. In an AIRSS search, one can generate initial structures with a specified number of randomly chosen symmetry operations. The effect of symmetry will be discussed further in Sec. III A 1.
As a random sampling approach, there are no rigorous convergence criteria to tell when to stop an AIRSS search. However, it is clear when the search should not be stopped. Searching should continue at least until the following requirements are fulfilled: (1) known marker structures, if any, are identified; (2) a number of

C. Density functional theory calculations
As part of an AIRSS search, the geometry optimization and the energy evaluations are carried out based on first principles DFT calculations. In this work, the plane-wave DFT code CASTEP was used for each geometry optimization. 36,37 The Perdew-Burke-Ernzerho (PBE) exchange-correlation functionals and the on-the-fly generated ultrasoft pseudopotentials (QC5) were adopted. 37,38 A relatively low plane-wave cutoff energy of 340 eV was used. The reciprocal space was sampled on Monkhorst-Pack grids with a spacing of 0.07 Å −1 . The effect of spin-polarization and the Hubbard U correction on the sampling of the PES will be discussed in Sec. III A 1.
For further structural optimization and property calculations, spin-polarized DFT calculations were performed using PBE exchange-correlation functionals with Hubbard U corrections. For this, we used the VASP code to make use of the calibrated U parameters for phase stability and voltage prediction. [39][40][41][42][43] A plane-wave cutoff energy of 520 eV was used, alongside Monkhorst-Pack grids with a spacing of 0.04 Å −1 . Hubbard U corrections of 5.3, 6.2, 3.25, 3.9, and 3.32 eV were applied to the d-channel of Fe, Ni, V, Mn, and Co, respectively. 41,43 The bandgap values were estimated using the HSE06 hybrid functional based on the PBE + U relaxed structures. 44 Ferromagnetic spin arrangements are assumed in spin-polarized DFT calculations.

D. Property calculations
The average cathode voltage Va is calculated using the following equation: where ECat−Li x and ELi are the enthalpy per atom of the cathode with a Li concentration of x and the enthalpy per atom of the lithium metal, respectively. e is the charge of an electron. It is worth noting that the voltage of a cathode is relatively difficult to accurately determine from simple DFT calculations. One source of inaccuracy is the inexact exchange-correlation functionals employed in current DFT schemes. 45 While a number of studies have been carried out showing that high-level hybrid functionals give good approximations to the voltages, in this work, we chose to use the PBE functional with a Hubbard U correction, which is capable of predicting the voltages of cathodes with satisfactory accuracy and is sufficient for screening purposes. 46 Another source of inaccuracy in voltage prediction comes from an overestimated structural stability in DFT calculations. While in reality, the cathode may go through phase transition/structural degradation upon delithiation, DFT calculations are usually performed on small cells, and the relaxation scheme is not able to capture these transitions effectively. Therefore, the average voltage calculated by relaxing the delithiated structures is an estimate of the upper bound. Similar to the voltage, the capacity is difficult to be predicted accurately as well. In principle, the capacity is determined by the thresholds of Li concentrations below/above which the structure starts to collapse or degrade.
However, large scale structural degradation is difficult to capture using DFT-based relaxation schemes as mentioned earlier. In this work, we estimated the upper bound of capacity by looking at the stability of the delithiated structures by structural relaxation. When over-delithiated, the C 2 O 4 anions decompose into linear CO 2 molecules. Furthermore, the gravimetric energy density of the cathode is estimated by multiplying the average voltage and the capacity followed by normalization by its weight.
The rate capability of cathodes can be estimated by calculating the Li diffusion barrier and the bandgap, which give information about the ion and the electron transport, respectively. The Li diffusion barrier was estimated using climbing image nudged elastic band calculations, 47 and the bandgaps were calculated using the HSE06 functional as described previously. 44

E. PES visualization
To present the results of our searches, we generated structure maps depicting the relative positions of each structure on the underlying PES. To do this, we first generated the Smooth Overlap of Atomic Position (SOAP) 48,49 description of each structure in a given dataset using the ASAP code, 50 characterizing each structure as a vector in a high-dimensional space. We used the universal SOAP parameters in ASAP, which are constructed based on the elements present in the structures, following the heuristics of Cheng et al. 50 Note that we included crossover terms between elements in these descriptors, and the global description of each structure is yielded by averaging over the local environments on an element-wise basis. The descriptor vectors were normalized so that each variable had a variance of unity across a given dataset.
We then applied the dimensionality reduction method Stochastic Hyperspace Embedding And Projection (SHEAP) to the highdimensional data to produce two-dimensional representations. 51,52 SHEAP constructs a weighted graph describing the similarity relationships in the source data. Following a test to combine equivalent structures, the data are cast into a low-dimensional space (2D in this case) using random projection. Another graph of weights is then constructed for the nodes defined by these mapped points. The layout of the map is then optimized by minimizing a cost-function that penalizes any mismatch between the weights.
In these maps, individual structures are represented by circles colored according to enthalpy, with areas proportional to the number of occurrences in the search. Circles are prevented from overlapping each other by a soft-core repulsion that is turned on only once the nearly optimal layout has been reached. The axes do not have any predetermined physical interpretation, but the relative positioning and clustering across the two-dimensional space reveals the structural similarity.

III. RESULTS AND DISCUSSION
A. Re-exploring known systems and the effect of searching parameters 1. LiCoO 2 -Impact of the species-pairwise minimum separations, cell volumes, and DFT parameters on search efficiency We start by performing an AIRSS search on LiCoO 2 , the most commonly used cathode material in LIBs. For this system, we The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp conducted the search using the minimum separation matrix and volume per formula measured from the known R3m high-temperature (HT) phase. 53,54 The values are listed in Table I. The number of formula units was chosen to take values of 1-5. The results of this search are summarized in Fig. 2(a). In total, 1200 random structures were generated and optimized. Of the relaxed structures, the trigonallayered HT phase had the highest encounter rate at 47%. The cubiclithiated-spinel low-temperature (LT) phase, which has a comparable energy, was also found multiple times (∼2%). 55 Other low energy structures found include the P2 and O2 phases of LiCoO 2 , which result from a slight variation of the CoO 6 layer stacking compared to the O3-type stacking of HT-LiCoO 2 . 56,57 All of these low energy structures are characterized by similar local coordinationcobalt ions are octahedral, while lithium ions are either octahedral or prismatic. This search also turned up many polymorphs with relatively high energies. These include several approximants of disordered rocksalts, in which cobalt and lithium ions are distributed in a non-layered manner, within a distorted MO 6 (M = Li and Co) framework. In these structures, the large lattice distortion gives rise to the higher energies. We expect these disordered rocksalts to be difficult to synthesize.
The above-mentioned search was carried out based on prior knowledge of LiCoO 2 . However, in exploratory searches, such prior knowledge may not exist. Now, we demonstrate that constraining the search space using the minimum separation matrix, cell volumes, and symmetry is critical to achieving efficient sampling. This is followed by a discussion on how one can use a preliminary AIRSS search to determine a suitable list of minimum separation values and cell volumes. Figures 2(b) and 3(c) summarize the search results obtained from the use of a matrix of minimum separations drawn uniformly from 1-3 Å and of a single fixed minimum separation for all species pairs of 1.5 Å, respectively. In comparison to Fig. 2(a), these searches resulted in many more high energy configurations with a significantly lower encounter rate for the lowest energy structures. In the case of a single fixed minimum separation, we can easily understand this by noting that the chemical environments of the cation and anion are not distinguished by this constraint. As a result, Li-Li and O-O close contacts are generated resulting in unphysical and high energy structures and chemistries.
As illustrated by the results of Fig. 2, an effective way to select suitable minimum separation values is to measure them from known structures with similar chemical compositions. However, in some cases, the detailed coordination environment is not known as a prior, and the minimum separations need to be selected using a first principles approach. A way around this is through a preliminary AIRSS search. Figure 3 shows the distribution of minimum separations arising in relaxed structures obtained from a search on LiCoO 2 using a random minimum separation matrix whose elements are drawn uniformly from 1-3 Å. No symmetry constraints were applied. We observed, for each pair of species, rapid convergence in the lower bound of the interatomic distances with respect to increasing the number of formula units, with the converged values corresponding to the minimum bond lengths of the system. Thus, for a "safe" search, one can use a relatively small cell to acquire these minimum separations, ruling out certain close contacts between pairs of species. However, these lower bounds do not consider the detailed local coordination of the system. For example, as shown in the inset of Fig. 3, the lower bound of Co-Co distance is ∼2.1 Å, which corresponds to Co pairs in an O-connected CoO 4 -CoO 4 configuration.
Compared with the CoO 6 octahedra, such local coordination is energetically less favorable. Therefore, for an "efficient" search, it is better to choose the minimum separations from the structure with the lowest energy, as indicated by the dashed lines in Fig. 3. This will bias the search toward having local coordinations that are more likely to occur in low energy structures. A search can be further speed up through the application of symmetry constraints. As illustrated in Figs. 2(a) and 2(d), for a well-chosen minimum separation matrix, the encounter rates for low energy structures are comparably high with and without the use of symmetry. However, the overall searching time is significantly longer (∼10 times) in the absence of symmetry due to slower DFT electronic steps during the geometry optimizations.
It is worth noting that the majority of the time taken for an AIRSS search is spent on the first principles geometry optimizations. Thus, speeding up the DFT calculations is highly beneficial. One way we achieved this was by adopting ultrasoft pseudopotentials (QC5) specially made for high-throughput calculations, which allow the use of low energy cutoff. Apart from that, the inclusion of spin polarization and Hubbard U parameters also affects the speed of calculation. Therefore, when necessary and tested, another way to accelerate the process is to perform non-spin polarized calculations during the search followed by a subsequent round of high accuracy calculations to refine the low energy structures. For LiCoO 2 , we compare the structural density of states obtained without spin, with spin, and with the Hubbard U correction. Here, the structural density of states is described as the proportion of structures at each energy. It evaluates the frequency of finding a local energy minimum with a specific energy on the potential energy surface. We found that the structural density of states obtained using different computational parameters is similar, and the low energy structures are all found with comparable encounter rates, as shown in Fig. 4(a). Figure 4(b) further shows the correspondence of energies obtained with non-spin polarized PBE and with spin-polarized PBE + U for the structures of LiCoO 2 . We observe a shift of stability order with an energy variation of ∼150 meV/atom. 32,33 For an AIRSS search carried out using non-spin polarized PBE, any structure with energy below this value would need to be recalculated to obtain the correct ordering. The magnitude of this threshold is of course system dependent. In extreme cases, the choice of DFT parameters may significantly alter the PES, possibly resulting in missing certain local minima. Therefore, we must note that The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp computational parameters (e.g., exchange-correlation functionals, spin-polarization, and correction schemes) must be chosen with caution. Preliminary and post-mortem searches need to be carried out to check whether the results are parameter dependent. If they are, the best practice is to report a confidence interval instead of a single prediction. In summary, our workflow for efficient AIRSS is as follows: First, obtain the species-pairwise minimum separations and the cell volume parameters, which constrain the search space, either by considering prior knowledge of the chemical system or by carrying out a preliminary search using random minimum separations; then, generate and relax random sensible structures with different numbers of formula units. Particularly in the case of many formula units, symmetry may need to be taken into consideration to accelerate the search. Additional re-optimization of low energy structures may need to be carried out if the search used relatively strong assumptions for the DFT calculations.

LiFePO 4 -Effect of adopting structural units on search efficiency
LiFePO 4 is a representative polyanion cathode material that is widely used in electric vehicles due to its low raw material cost. 58 LiFePO 4 is structurally more complicated than LiCoO 2 and therefore is a more challenging task for crystal structure prediction. In our application of AIRSS to this system, we generated random structures fulfilling the following requirements: (1) the structure has 2 to 4 symmetry operations; (2) the P and O species are generated as proper PO 4 units; (3) the interatomic distances are larger than the predefined species-pairwise minimum separation values, which are obtained from the known olivine phase (the P-O distances are specially treated so that only those measured from inter-PO 4 units are constrained to fulfill this requirement); (4) the volume per formula unit is within 15% as compared with that of the known olivine phase. The search results obtained from relaxing these random sensible structures are summarized in Fig. 5. Out of the 684 unique structures found, the experimentally observed olivine phase (Pnma) turned up four times, giving an encounter rate of ∼0.58%. Also found was a Cmcm structure, slightly higher in energy than the olivine phase, which has previously been reported under high pressure conditions. 59 Interestingly, many other phases with similar energies to this were located. Some of these are likely to be synthesizable under fine-tuned conditions and thus may be worth further computational or experimental study.
In our study of this system, we observed that the use of predefined chemical subunits (PO 4 ), as opposed to individual atoms (P and O), was effective in steering the search toward lower energy structures. This is illustrated in the inset of Fig. 5 through the structural density of states obtained with and without this constraint. Caution is needed with such constraints, however, since they restrict the sampling space. They should be used only when the chemistry is well understood. 60 In the case of LiFePO 4 , structures containing other polyphosphate anions are not energetically competitive and so can be disregarded, allowing us to accelerate the search by imposing P and O to be present as tetrahedral units of PO 4 .

LixCuyFz-A multi-compositional search for phase stability
The cases of LiCoO 2 and LiFePO 4 demonstrated the power of AIRSS for compounds of a given single composition. However, in exploratory searches, one does not know in advance which

ARTICLE
scitation.org/journal/jcp composition will give the best performance or the desired properties. Therefore, multi-compositional searches followed by post-search screening are necessary. More importantly, the determination of the stability of an unknown phase not only involves the calculation of its formation energy but also whether it is susceptible to decomposition into other phases. Therefore, it is necessary to predict the entire compositional phase diagram. For cathode materials, this is especially important when the chemistry of the system is less understood. The construction of the phase diagram is critical in determining the reaction pathways during lithiation since not all cathodes are of the intercalation type and topotactic phase transformation may happen. For example, transition metal fluorides, for example, FeF 2 , FeF 3 , and CuF 2, are typical cathode materials that go through a conversion reaction during discharge. 61 Here, we consider the system of Li-Cu-F, demonstrating how a ternary phase diagram can be generated by conducting an extensive multi-composition search in the chemical space. 62 In particular, we searched for structures with compositions with the formula LixCuyFz, where x, y, and z fulfill the following requirements: x + y + z = n, n < 10 and x, y, z ∈ Z+, where Z+ is the set of positive integers. This led to 163 compositions. For each composition, we generated at least 50 random structures for further first principles optimization. In total, 15 031 structures were generated and optimized. This number was found just enough to cover the stable phases in this study. In exploratory searches where the chemistries are unknown, many more structures may need to be considered. The number of formula units was restricted to being 1 and 2 in the current search. The minimum separation matrix and the cell volumes were determined from a preliminary search that used random species-pairwise minimum separations of 1-3 Å with ∼50 structures being sampled on each composition. The details are discussed in Sec. III A 1. All the structures generated were constrained to have 2 to 4 symmetry operations.
As shown in Fig. 6(a), multiple stable phases were obtained from the search. For the binary compounds, LiF, CuF 2 and Cu 3 Li were found to be stable. Such a result helps determine the terminal redox couple of CuF 2 as a cathode, i.e., the final thermodynamically stable reaction product should be Cu and LiF. The average voltage is calculated to be 3.64 V. For ternary compounds, multiple compounds were found to shape the convex hull, including Li 2 CuF 6 , LiCu 2 F 6 , and LiCuF 4, in close agreement with existing databases. 16 Since none of these stable phases are located within the CuF 2 -LiF-Cu triangle, CuF 2 is expected to adopt a single-step conversion reaction to form LiF and Cu upon lithiation. Interestingly, we also found a number of phases that are less stable (above the convex hull) but are potentially interesting based on kinetics. For example, the phase stability along the Li-CuF 2 tie line is shown in Fig. 6(b). The hypothetical Li-intercalated compound LiCuF 2 is high in energy (220 meV/atom above hull). Such an instability explains why CuF 2 goes through a conversion reaction, rather than intercalation, when lithiated. We also found phases containing F − anion lattices that are topotactical with those in LiF, as shown in Fig. 6(c). These have recently been proposed as possible buffer phases between the The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp reactant CuF 2 and the reaction product LiF. 63 Cu (I) containing fluorides have also been proposed as potential intermediate phases, especially during the charging process. Xiao and co-workers proposed Cu (I) containing phases with tetrahedrally coordinated Cu following electron diffraction analysis. 62 Our search reveals possible phases that are consistent with this. As shown in Fig. 6(e), the Li 3 CuF 4 phase has an energy above the hull of 60 meV/atom. Despite its relatively high energy, the formation of this at interfaces and under electrochemical conditions, is plausible.
To conclude, a multi-composition search can efficiently generate the phase diagram for prediction of phase stabilities and for analysis of the reaction mechanisms of a cathode.

B. Exploring new cathodes
Transition metal oxalates are a family of polyanionic compounds that possess a number of advantageous characteristics as cathodes for batteries but have been overlooked as compared to other polyanionic materials such as the phosphates PO 4 3− , sulfates SO 4 2− , and silicates SiO 4 4− . 64,65 Interestingly, the C 2 O 4 2− anion has a similar polarizability to PO 4 3− . As a result, it should be able to deliver reasonable redox potentials during lithiation, such as LiFePO 4 . Furthermore, oxalates can be synthesized under mild conditions via solution-based methods. In fact, they are known to nucleate into a chain-like, layered, and fully connected frameworks, thanks to complex ligand coordination chemistry. This offers exciting opportunities for oxalates to serve as cathode materials by tuning the structural features through metastability. Despite these advantages, oxalates have been overlooked as cathode materials, in part, due to the large molecular weight of the anion, which leads to relatively low energy density. 65 Recently, with the discovery of the reversible anion redox as an additional source of capacity, the oxalates are gaining attention from battery researchers. Jiao et al. synthesized a lithium-containing iron oxalate with composition Li 2 Fe(C 2 O 4 ) 2 . 66 In this material, Fe The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp has an oxidation state of +2. Interestingly, they found that more than 1 Li/f.u. can be taken out of the structure reversibly. Using spectroscopic evidence, they found that this over-delithiation is enabled by the reversible electron stripping from the C 2 O 4 2− anion. As compared to the oxides, especially Li-rich ones, the strong covalent O-C bonds in C 2 O 4 2− prohibit the O 2 2− -dimer formation and the O 2 evolution during anion oxidation.
Given the possible polymorphism and interesting chemistry in Li-stuffed transition metal oxalates, we carried out exploratory AIRSS searches on a range of systems of this nature. In particular, we considered Li 2 TM(C 2 O 4 ) 2 , where TM = Fe, Co, Ni, V, and Mn. The aim was to find low energy phases with the potential to be synthesized and display reasonable energy density and rate capability as cathodes for LIBs. Figure 7(a) summarizes our AIRSS searches for Li 2 Fe(C 2 O 4 ) 2 in the form of a SHEAP map. During the search, C and O were treated as C 2 O 4 units when generating the initial random structures, and the number of formula units was restricted to being 1-4. As illustrated in the SHEAP map, among the large number of local energy minima found, only a fraction are found to be low in energy, emphasizing the complexity of this PES. Representative low energy structures are labeled in Fig. 7(a). The structures with the highest encounter rates are layered, in which the Li + cations are distributed within the layers. These structures resemble the high-temperature phase of LiCoO 2 . A number of "checkerboard" structures, composed of one-dimensional structural units of octahedral-coordinated Fe and planar square-coordinated Li, were also encountered frequently. Table II lists the energetics of the low energy phases of Li 2 TM(C 2 O 4 ) 2, whose structures are given in Figs. 7(d)-7(i) using Fe compounds as an example. Interestingly, for Fe, Co, and Ni, we found that at least one of the layered polymorphs has an energy lower than the experimentally reported structure at the PBE level. For Fe, this polymorph is characterized by its well-connected Linetwork and is expected to conduct lithium more efficiently than the experimental phase. 67 This is confirmed through calculation of the diffusion barrier of Li in the two structures. In the layered structure, the barrier for Li diffusion is ∼320 meV, close to some superionic conductors, whereas for the experimentally reported structure, the value is ∼620 meV. 66,68 This indicates that the layered phase has the potential to display high-rate capability as a cathode.
In our search, as well as this layered ground-state structure, a number of metastable structures were found that may be synthesizable through complex nucleation kinetics in solution. Indeed, we calculated that the experimentally reported Li 2 Fe(C 2 O 4 ) 2 phase is also intrinsically unstable, with an energy above the hull of 58 meV/atom. Here, we used such a value as a reference and increased it slightly (by 20 meV/atom) as the threshold for selecting Li 2 TM(C 2 O 4 ) 2 structures on which to conduct further property calculations. These are the structures before the dashed line in Fig. 7(c). The average discharge voltages of structures below this threshold are listed in Table II. For Fe, all of them are over 4V vs Li/Li + , as expected, since the anion redox is involved. We estimated upper bounds for capacities of these materials by checking the structural integrity of the delithiated phases at different states of charge. We found that when the structure is over-delithiated, the C 2 O 4 2− anion decomposes into two CO 2 molecules. Using whether such decomposition happens during structural relaxation as a criterion, the estimated capacity upper bound is listed in Table II. Interestingly,   TABLE II. Properties of Li 2 TM(C 2 O 4 ) 2 polymorphs within an energy threshold ΔE of 20 meV/atom with respect to the most stable phase found. Here, ΔE is the relative energy compared with the most stable polymorph, Eg is the electronic bandgap, and Va is the average voltage vs Li/Li + . "Exp" and "CB" are abbreviations for "experimental" and "checkerboard," respectively. among the TMs we considered, V-oxalates seemed to give the highest capacity without structural collapse, whereas Ni-oxalates easily decomposed when delithiated. The stability of V-oxalates may be due to the multi-valence nature of V. Furthermore, by multiplying the capacity and the voltage, we computed an estimate for the energy density of oxalates. A number of polymorphs (Fe-Layered #1, Mn-Layered #2, Mn-Mixed, Co-Layered #2, and all Vpolymorphs) appear to have the potential to reach extraordinary gravimetric energy densities (>900 Wh kg −1 ), close to or even higher than that of LiNi 0.8 Co 0.1 Mn 0.1 O 2 . It is worth noting that, in these predictions, temperature effects are not considered. The entropy term arising from partially occupied Li sites/different Li arrangements at finite temperatures is not incorporated. While this will not significantly alter the energy densities predicted here and is sufficient for screening purposes, identifying disorder and computing the corresponding free energy is critical in analyzing the phase separation behavior of cathodes during discharge. Further methodological development is required in this context, such as the incorporation of lattice-based sampling methods, for example, cluster expansions.

IV. CONCLUSIONS
We have explored ab initio random structure searching as a tool for battery cathode discovery. We showed, through case studies of known cathode materials of LiCoO 2 , LiFePO 4 , and LixCuyFz, that efficient prediction of the low energy structures can be achieved. We demonstrated that using species dependent minimum interatomic separations and cell volumes to steer the search is critical to its efficiency. Further acceleration was achieved by imposing symmetry constraints and by using pre-defined chemical subunits for structure generation. Based on these principles, we have carried out explorative searches for Li 2 TM(C 2 O 4 ) 2 oxalates. A number of low energy polymorphs with layered structures were found. Post-search screening and property calculations identified promising cathode materials characterized by high voltages (>4V) and low Li diffusion barriers (∼300 meV), with the potential to deliver not only high energy densities (up to ∼900Wh kg −1 ) but also high-rate capability. This work serves as a basis for AIRSS-based cathode searches, providing a detailed framework for an efficient workflow.

AUTHORS' CONTRIBUTIONS
Z.L. and B.Z. contributed equally to this work.