The decisive role of electrostatic interactions in transport mode and phase segregation of lithium ions in LiFePO₄

Abstract

Understanding the mechanism of slow lithium ion (Li⁺) transport kinetics in LiFePO₄ is not only practically important for high power density batteries but also fundamentally significant as a prototypical ion-coupled electron transfer process. Substantial evidence has shown that the slow ion transport kinetics originates from the coupled transfer between electrons and ions and the phase segregation of Li⁺. Combining a model Hamiltonian analysis and DFT calculations, we reveal that electrostatic interactions play a decisive role in coupled charge transfer and Li⁺ segregation. The obtained potential energy surfaces prove that ion–electron coupled transfer is the optimal reaction pathway due to electrostatic attractions between Li⁺ and e⁻ (Fe²⁺), while prohibitively large energy barriers are required for separate electron tunneling or ion hopping to overcome the electrostatic energy between the Li⁺–e⁻ (Fe²⁺) pair. The model reveals that Li⁺–Li⁺ repulsive interaction in the [010] transport channels together with Li⁺–e⁻ (Fe²⁺)–Li⁺ attractive interaction along the [100] direction cause the phase segregation of Li⁺. It explains why the thermodynamically stable phase interface between Li-rich and Li-poor phases in LiFePO₄ is perpendicular to [010] channels.

---

1. Introduction

Olivine LiFePO₄ is one of the prevalent cathode materials in lithium-ion batteries for electric vehicle applications due to its distinctive advantages, such as low cost, high safety, and outstanding structural stability.^1–3 However, LiFePO₄ suffers from sluggish kinetics of charge transport, leading to poor intrinsic electronic conductivity (10⁻⁹–10⁻¹⁰ S cm⁻¹) and ionic diffusivity (10⁻¹²–10⁻¹⁷ cm² s⁻¹).^4–9 Substantial progress has been made over the past two decades to improve the rate performance of LiFePO₄ by surface coating,^10–12 ion doping,^13–16 and control of grain size.¹⁷ It is worth noting that the electronic conductivity and ion diffusion rate increase simultaneously, which signifies a synergy between ion transport and electron transfer (ET).

Both experimental and computational studies have revealed the phenomenon of Li⁺ diffusion coupled ET in LiFePO₄.^18–23 Ellis et al.²¹ provided the first experimental clue to coupled Li⁺ and electron mobility based on Mössbauer spectrum measurements. The small electron polaron hopping is precisely related to Li⁺ disorder and vibrational (phonon) modes in the lattice. Maxisch et al.¹⁹ and Sun et al.²⁰ clarified large binding energy and strong correlation between Li⁺ and e⁻ in migration paths from first-principles calculations. Besides, Fraggedakis and Bazant²⁴ developed a theory of ion–electron coupled transfer, revealing that ions and solvent molecules fluctuate cooperatively to facilitate non-adiabatic ET.

In our earlier research,^6,25 we developed a kinetic model for Li⁺ transport in LiFePO₄, which unravels that the ion mobility in LiFePO₄ is coupled with electron transfer between Fe²⁺ and Fe³⁺. The electronic coupling effect is incorporated into the pre-exponential factor of the ionic diffusion coefficient (D_Li⁺) by considering the electronic transmission coefficient (κ_el) in Landau–Zener theory,^26–28 given by

where . v_n is the effective frequency of nuclear motion and is the electron hopping frequency. V_AB is the electronic coupling matrix element, which describes the coupling strength between electron orbitals of the redox centers. This kinetic model quantitatively reveals that the very small value of V_AB in LiFePO₄ results in lower ionic diffusivity. Li⁺ hopping is in conjunction with a nonadiabatic electron transfer. In the meantime, we employed a neural network-based deep potential (DP) to calculate hopping barriers of Li⁺ with various Li–Li coordination environments.⁴ The results show that the gradient distributions of Li⁺ cause great asymmetry in the barriers of forward and backward ion hopping. It results in a very slow Li⁺ diffusion rate and a diverse variation of D_Li⁺ based on the Kinetic Monte-Carlo (KMC) simulation.

Despite the abovementioned progress, a better understanding of the influence of electronic coupling on ion transport is needed. Specifically, the way electron transfer affects ion transport and its impact on the ion hopping barrier (E_a) remain unclear. Herein, we aim to develop a theoretical framework for exploring the mechanism of ion transport in the LiFePO₄ lattice and further unveiling the effect of electronic coupling on the activation barrier. Properly considering the relaxation of the lattice structure and the electrostatic interaction between Li⁺ and e⁻ (Fe²⁺), a model Hamiltonian is constructed for this purpose. Potential energy surfaces show that ion–electron coupled transfer is the optimal path with the lowest energy barrier. Additionally, the model is further employed to interpret the phase segregation behavior of LiFePO₄ at equilibrium in terms of electrostatic interaction.

2. Theory development

2.1 Model specifications and Hamiltonians

Li⁺ transport in LiFePO₄ corresponds to an ion-vacancy exchange process accompanied by an electron–hole exchange procedure, which is described as a self-exchange reaction.^6,29–32

Li⁺_iFe²⁺PO₄(m) + V_i±1Fe³⁺PO₄(m) → V_iFe³⁺PO₄(n) + Li⁺_i±1Fe²⁺PO₄(n) (2)

Li⁺_i represents Li⁺ at the site i, and V_i±1 is a vacancy located at site i ± 1. Li⁺ can hop between the nearest-neighbor lattice sites (Fig. 1). m and n represent the lattice environment of Li⁺–e⁻ (Fe²⁺)/V–hole (Fe³⁺) pairs, which is determined using the concentration and site of Li⁺ in LiFePO₄. In the case of the same lattice environment before and after Li⁺ hopping, one has m = n, which is defined as the symmetric ion hopping process in LiFePO₄. Otherwise, the process is referred to as the asymmetric ion hopping process.

Fig. 1 Schematic illustration of ion hopping and electron transfer in LiFePO₄ and the model framework.

Computational and experimental studies^33–35 demonstrated that LiFePO₄ has an olivine structure belonging to the space group Pnma, in which Li⁺ migrates along the b-axis ([010]) direction. Phosphorus ions occupy tetrahedral sites, and Li ions are located in chains of edge-sharing octahedral sites, which are adjacent to FeO₆ octahedra. In the transition state, Li ions occupy tetrahedral sites. As shown in Fig. 1, the ion transport process includes Li⁺ hopping between different octahedral sites, electron transfer from Fe²⁺ to Fe³⁺, and the relaxation of the lattice framework, i.e. the vibration of P–O and Fe–O bonds. These processes induce changes in lattice relaxation energy and electrostatic interaction energy between ions and electrons, which are major origins of activation barriers for ion transport.

The ion transport process can be described by the coordinate of Li⁺, d, which is the b-axis lattice parameter, and the overall harmonic polarization reaction coordinate, x. Then, we construct model Hamiltonians of reactant and product systems as²⁴

H_sys,R(x,d) = H_sys,R(x_R,d_R) + H_rex,R(x) + H_els,R(x,d) (3a)

H_sys,P(x,d) = H_sys,P(x_P,d_P) + H_rex,P(x) + H_els,P(x,d) (3b)

where H_sys,R/P(x_R/P,d_R/P) is the equilibrium energy of the reactant/product system, corresponding to the state where both ions and electrons are located in the stable configuration of reactant R/product P. The change in volume and configuration entropy of LiFePO₄ is negligibly small throughout the ion transport process,^36,37 so the volumetric and entropic effects are not considered. The functions H_rex and H_els describe the change of the lattice relaxation energy and the electrostatic interaction energy, respectively.

2.1.1 Lattice relaxation energy term H_rex. The change in P–O and Fe–O bond lengths triggered by electron transfer can be approximated as a simple harmonic vibration, so the lattice relaxation energy of reactant and product systems is represented by parabolic curves in terms of the polarization coordinate x in the near equilibrium region,

H_rex,R(x) = k_R(x − x_R)² (4a)

H_rex,P(x) = k_P(x − x_P)² (4b)

where k_R and k_P are fixed force constants for reactants and products, respectively. The force constants relate directly to the configuration reorganization energies (λ) for the vibration of all bonds during ion transport and electron transfer.

2.1.2 Electrostatic interaction energy term H_els. Based on coulomb forces or charge–charge interactions,³⁸ the change of Li⁺–e⁻ (Fe²⁺) electrostatic interaction energy readswhere z_Li⁺ is the charge of Li⁺ and q_n represents the partial charges of an electron at different sites. ε₀ and ε_r are the permittivity of vacuum and the dielectric medium. The channels of Li⁺ migration in LiFePO₄ are similar to vacuum, thus ε_r = 1. r(d) is the distance between Li⁺ and e⁻ (Fe²⁺). The lattice structure of FePO₄ is regarded as an electrically neutral background. The first term on the right hand side (rhs) describes the electrostatic interaction energy of a Li⁺ at any position and an electron located in the stable configuration of reactant R/product P. The second term on the rhs refers to the electrostatic interaction energy of both a Li⁺ and an electron located in the stable configuration of reactant R/product P.

2.2 Model solution for the energy surface

Substituting eqn (4) and (5) into (3), the model Hamiltonians read,

The energy diagram can be represented by a two-dimensional surface, which is called the potential energy surface (PES).

2.2.1 Iso-energetic curve and activation energy. Combining eqn (6a) and (6b), we can obtain the iso-energetic curve of reactants and products. It can also be defined as an activation energy curve. The details are discussed in Section 3. For the symmetric ion hopping in LiFePO₄, the force constants k_R = k_P = k. The iso-energetic curve is expressed aswithwhere ΔG⁰(d) is written as

ΔG₀(d) = H_sys,P(x_P,d_P) − H_sys,R(x_R,d_R) + H_els,P(x_P,d) − H_els,R(x_R,d) (9)

In the transition state of the optimal reaction path, the ΔG^≠_sym,R/P(x_T,d) arrives at one of the extrema, which is the minimum of the iso-energetic curve, that is, dΔG^≠_sym,R/P/dd = 0. Here, the value of ion coordinate d is expressed as d_T. The minimum activation energy in the forward and backward directions reads

E^R/P_a,sym(x_T,d_T) = ΔG^≠_sym,R/P(x_T,d_T) (10)

which is the energy barrier of the optimal reaction path.

The structure of LiFePO₄ is the same before and after the symmetric ion hopping process, meaning ΔG⁰(d_T) = 0. Then, eqn (10) reduces to

with d_T=(d_R + d_P)/2.

Similarly, for the asymmetric ion hopping process (k_P ≠ k_R), the iso-energetic curve is written as

with

The activation energies of the optimal reaction path read

E^R_a,asym(x_T,d_T) = ΔG^≠_asym,R(x_T,d_T) (14a)

E^P_a,asym(x_T,d_T) = ΔG^≠_asym,P(x_T,d_T) (14b)

with d_T satisfying dΔG^≠_sym,R/P/dd = 0.

If ΔG⁰(d_T) = 0, eqn (12) reduces to

2.3 Comparison with classical charge transfer theories

The model Hamiltonian describes ion transport coupled non-adiabatic electron transfer. It returns to the classical charge transfer theories by ignoring ion effects, as shown in Fig. 2. For the symmetric electron transfer (k_R = k_P = k), eqn (7) is rewritten aswith ΔG₀ is a constant, ΔG₀ = H_sys,P(x_P,d_P) − H_sys,R(x_R,d_R). Eqn (16) is the activation energy of non-adiabatic electron transfer in the classical Marcus theory,^39–43 which was originally developed for homogeneous ET kinetics in solutions. It assumes a very weak electron coupling between active species.⁴⁴ By contrast, the ET process proceeds in the adiabatic mode when electron coupling is sufficiently strong. It is described by introducing an electron transmission coefficient into the standard rate constant.^6,41 Adiabatic ET had also been studied by Schmickler⁴⁵ for electrocatalysis by combining the Anderson–Newns model⁴⁶ and Marcus theory. Schmickler developed a theoretical model constituting an important straddle from weakly coupling systems, studied using the Levich–Dogonadze–Kuznetsov theory⁴⁷ and the Marcus theory, to strongly coupling systems which electrocatalytic reactions belong to.⁴⁸

Fig. 2 Potential energy curves of reactant and product systems for (a) the symmetric and (b) the asymmetric electron transfer.

In the same way, for the asymmetric electron transfer (k_R ≠ k_P), eqn (12) reads

The applicability of eqn (17) is not limited by the degree of difference between k_R and k_P. It is different from the activation energy of the asymmetric Marcus–Hush (AMH) model adopted by Compton et al. in their investigation of electrode processes.^44,49–51 The AMH model takes into account the difference between inner-shell force constants of oxidized and reduced species by introducing an asymmetry parameter into the symmetric Marcus–Hush model. The accuracy is found to decrease as the difference between the force constants increases due to the first-order truncation of the approximating series from which it was derived.^41,49,52

2.4 Model parameters

Model parameters are listed in Table S1.† To describe the potential energy surfaces, several parameters related to the lattice environment are detailed here. The force constants for reactants and products are calculated using reorganization energies, k_R/P = λ_1/2/(x_P − x_R)². For the self-exchange reaction, eqn (2), the reorganization term is the sum of reorganization energy of individual reactants,^31,32

λ_1/2 = λ_{Fe²⁺→Fe³⁺} + λ_{Fe³⁺→Fe²⁺} (18)

where λ_{Fe²⁺→Fe³⁺} = E_FePO4+Li − E_LiFePO4 is the energy required to alter Fe–O bond length in LiFePO₄ under equilibrium, x⁰₂, to the distance, x⁰₃, in the equilibrium state of FePO₄. λ_{Fe³⁺→Fe²⁺} = E_LiFePO4–Li − E_FePO4 is the energy required to change the Fe–O bond length in FePO₄ from x⁰₃ to x⁰₂. Here, E_FePO4+Li is the energy of LiFePO₄ in the FePO₄ configuration and E_LiFePO4 is the energy of stable LiFePO₄ from the DFT calculations. Similarly, E_LiFePO4–Li and E_FePO4 are the energies of FePO₄ in the LiFePO₄ configuration and stable FePO₄, respectively. More details are shown in Fig. S1† and can be found in the original papers.⁶

The distance of Li⁺ and e⁻ (Fe²⁺), r_n(d), is calculated using

where (a,d,c) and (a_n,d_n,c_n) are the three-dimensional coordinates of Li⁺ and e⁻ (Fe²⁺). The trajectory of Li⁺ is fitted with the Fourier series, as shown in Fig. S2.†

It is essential to describe the distribution of an electron. Previous theoretical work^19,53 and numerous experimental studies^21,54 proved that charge carriers form small polarons in LiFePO₄. The initial/final states of Li⁺ transport contain a Li⁺–e⁻ pair with the electron located at a Fe site adjacent to the Li⁺. There are six adjacent Fe sites centered on the Li⁺ (Fig. 3). The Li⁺–e⁻ pair corresponds to an assembly of a variety of microscopic states with the electron appearing at one of the six Fe sites.⁵⁵ The probabilities of each microscopic state are different owing to the variable Li⁺–e⁻ (Fe²⁺) distance during the Li⁺ transport. Steady-state DFT calculations only capture one of these microscopic states.^55–58 For instance, the electron is localized on the Fe₁ or Fe₄ site for Li_a and Fe₄ or Fe₇ site for Li_b due to the minimum Li⁺–e⁻ (Fe²⁺) distance, d_Fe1–Lia = d_Fe4–Lia = d_Fe4–Lib = d_Fe7–Lib = 3.317 Å (Fig. 3a). Nevertheless, when the electron is localized on the Fe₄ site, it may not transfer during Li⁺ hopping from site a to site b. This disobeys the coupled ion–electron transfer mechanism. Over the dynamic course, the localized electron can transfer from one Fe atom to another in a period of time that is much shorter than that of Li⁺ transport. Macroscopically, Li⁺ is subjected to an average electrostatic effect from different Fe sites.

Fig. 3 Partial charges of an electron for (a) the isolated state and (b) the aggregated state of Li⁺ in the LiFePO₄ structure. The dimension unit of Li⁺–e⁻ (Fe²⁺) distance is Å.

Therefore, we approximately assign an electron localized on six Fe sites during Li⁺ hopping (Fig. 3). The charge distribution of electrons is assumed to follow two principles based on the short-range electrostatic interaction between Li⁺ and e⁻ (Fe²⁺). Firstly, it can be inferred that the probability of an electron localized at a Fe site is inversely proportional to the distance of Li⁺ and e⁻ (Fe²⁺). Secondly, electrons are more likely to be in the middle site of the two adjacent Li⁺. Consequently, the partial charge function of an electron, q_n, is constructed as

where p_n(d) is the distribution probability of an electron being localized at one of six Fe sites closest to Li⁺, which depends on the value of r_n(d). α_n is the partition coefficient of electrons related to the aggregation state of Li⁺. If there is an isolated Li⁺, α_n = 1 for all Fe sites, as shown in Fig. 3a. Two Fe sites in LiFePO₄ display point symmetry around the middle Li⁺, thus partial charges at these two Fe sites are equal. When two sides of Li⁺ are flanked by a Li⁺ and a vacancy, the values of α_n take 1.2 and 0.8 for the Fe site located in the middle of two adjacent Li⁺ and the middle of a Li⁺ and a vacancy, respectively (Fig. 3b). There is the structure of Li⁺–e⁻ (Fe²⁺)–Li⁺.

3. Results and discussion

3.1 Symmetric ion hopping process

Ion transport in pure Li_xFePO₄ and Li_1−xFePO₄ (x ≈ 0) are two typical examples. Li⁺ is isolated in FePO₄, while there is the structure of Li⁺–e⁻ (Fe²⁺)–Li⁺ in LiFePO₄. We take the reactant system as the reference, namely, H_sys,R(x_R,d_R) = 0. DFT-calculated parameters are λ^FePO_1⁴ = λ^FePO_2⁴ = 1.167 eV and λ^LiFePO_1⁴ = λ^LiFePO_2⁴ = 1.380 eV. Two symmetric ion hopping processes are depicted in Fig. 4 and S3.†

Fig. 4 (a) Structure diagram of Li⁺ hopping along the b-axis in FePO₄. (b) Potential energy surface based on the Hamiltonian model (left) and its top view (right) for ion hopping and electron transfer. (c) Potential energy curves of paths A, B, and C.

The potential energy surface as a function of x and d is informative, as shown in Fig. 4b. The potential energy surface consists of two 2D parabolas corresponding to the energy of the reactant system and product system, respectively. The two lowest points are the equilibrium state R(x_R,d_R) and P(x_P,d_P). The intersection of two parabolic surfaces, M − T^≠ − N, is the iso-energy curve described by eqn (7). In principle, electrons can tunnel from the reactant state to the product state for any value of ion coordinates d along the iso-energy curve. For example, in the case where d = d_R there is a probability for electron transfer without ion hopping (Path A). The final state Q(x_P,d_R) is the configuration of the electron tunneled to the product state and Li⁺ located at the reactant state (Fig. S4†). The energy barrier is 0.791 eV, which is prohibitively large due to the short-range electrostatic interactions of the Li⁺ and e⁻ (Fe²⁺) (Fig. 4c).

In contrast, there is also a chance for ion hopping without electron transfer (Path B), which must overcome the energy barrier of 0.757 eV. The final state O(x_R,d_P) represents that Li⁺ is located in the product state while the electron is still localized at the reactant state. The energies of states O(x_R,d_P) and Q(x_P,d_R) are both 0.757 eV higher than that of the initial state R(x_R,d_R). It is the electrostatic energy required to eliminate an ion–electron pair.

Reactant species fluctuate on the 2D energy surface during ion hopping and electron transfer. The fluctuating trajectory most likely follows the lowest energy path C. The transition state T^≠(x_T,d_T) is the lowest point of the iso-energy curve. Along this optimal reaction pathway, ion hopping is coupled with electron transfer. As shown in Fig. S5,† partial charges are located at six Fe sites around Li⁺ in the initial state. The values of partial charges vary with the distance between Li⁺ and e⁻ (Fe²⁺). Once Li⁺ arrives at the tetrahedral site and the lattice structure relaxes to the iso-energetic state of reactants and products, electron transfer will occur. Then Li⁺ hops to the final state with the redistribution of partial charges.

As shown in Fig. 4b and S3b,† the minimum activation energy of ion–electron coupled transfer is 0.253 eV and 0.336 eV for FePO₄ and LiFePO₄, respectively. The values are much lower than those of separate electron transfer (0.791 eV and 0.827 eV) or ion transport (0.757 eV and 0.816 eV). According to eqn (11), we calculate the value of H_els(x_T,d_T) in the transition state, which is the difference in the electrostatic interaction energy between the transition state and the initial state. The results are −0.0392 eV in FePO₄ and −0.009 eV in LiFePO₄. These negative values originate from a closer distance of Li⁺ and e⁻ (Fe²⁺), r_n(d), and the higher charges, q_n, located in the nearest Fe site in the transition state than that in the initial state. Compared with the higher energy barriers of paths A and B, these results clearly illustrate that ion–electron coupled transfer is driven by the electrostatic attraction of Li⁺ and e⁻ (Fe²⁺).

DFT calculated barriers correspond to the coupled ion–electron transfer process. The values serve as a reference for the model, while the model with an analytical nature can help understand the components of the DFT-calculated barriers. Different contributions, such as the configuration reorganization and the electrostatic interaction, can be separated and compared. The results in this model are in good agreement with those of NEB calculations (Fig. S6†), which are also close to literature values (0.24–0.27 eV for FePO₄ and 0.32–0.39 eV for LiFePO₄).^6,59–61

3.2 Asymmetric ion hopping process

The charge transport in Li_0.125FePO₄ of a 3a × 2b × 1c supercell is depicted in Fig. 5. In the initial state, Li⁺ accumulate in the b-axis ([010]) transport channel, while the hopping Li⁺ is adjacent to another Li⁺ along the a-axis ([100]) at the final state. The partial charges at the initial and final states are shown in Table S4.† The reorganization energy λ is different for the two structures with λ_initial = λ₁ = 1.571 eV and λ_final = λ₂ = 1.647 eV.

Fig. 5 (a) Structure diagram of the initial and final states of Li⁺ hopping. (b) The activation energy calculated using the NEB for asymmetric ion hopping in Li_0.125FePO₄, and (c) potential energy surface based on the Hamiltonian model (left) and its top view (right) for ion transport and electron transfer.

Fig. 5c displays the asymmetric potential energy surface of ion–electron coupled transfer. The activation energy of the optimal reaction path C is 0.246 eV, which is reasonably close to the 0.231 eV calculated using the NEB method (Fig. 5b). The value is much lower than those of path A (0.646 eV) and path B (0.817 eV) for separate electron transfer and ion hopping, as shown in Fig. S7.† The energy of the final state is about 0.192 eV less than that of the initial state, which indicates that the structure is more stable when Li⁺ amass in the [100] direction (a-axis) rather than the [010] direction (b-axis). This is mainly caused by the difference in electrostatic interactions in two directions, which is discussed in the next section.

3.3 The origin of ion aggregation in LiFePO₄

In the [010] direction, the distance between two adjacent Li⁺ is 3.045 Å, which is smaller than all three kinds of distances between the Li⁺ and e⁻ (Fe²⁺), 3.317 Å, 3.548 Å, and 3.702 Å (Fig. 6a). However, the distance between adjacent Li⁺ is 4.746 Å in the [100] direction, larger than the two kinds of distances between Li⁺ and e⁻ (Fe²⁺), 3.548 Å and 3.702 Å (Fig. 6b). The complex electrostatic interactions include the electrostatic repulsion of the Li⁺–Li⁺ couple and the electrostatic attraction of the Li⁺–e⁻ (Fe²⁺) pair. Based on the partial charge function, eqn (20), the electrostatic interaction energy between the two adjacent Li⁺ is calculated. In the [010] direction, the electrostatic interaction energy is 1.237 eV. The positive value indicates that the two Li⁺ are repulsive. In contrast, the two Li⁺ exhibit attractive interactions under the attraction of e⁻ (Fe²⁺) when they collect in the [100] direction with the electrostatic interaction energy being −0.029 eV.

Fig. 6 Structural diagrams of Li⁺ aggregation (a) in the [010] direction and (b) in the [100] direction; purple and green arrows represent electrostatic attraction and repulsion, respectively. The dimension unit is Å.

The structure of Li⁺–e⁻ (Fe²⁺)–Li⁺ is attractive in the [100] direction, which induces phase separation into Li-rich and Li-poor phases in Li_xFePO₄, as reported in previous studies.^62,63 Besides, the results also demonstrate that Li⁺ are more likely to amass in the [100] direction. It means that the phase interface of Li-rich and Li-poor phases is perpendicular to the [010] axis, which is consistent with observations of HRTEM images.⁶⁴ In the meantime, it should be pointed out that the distance of layer to layer in the c-axis ([001]) is greater than 5 Å, which results in negligibly small interlayer interactions between Li ions.⁴ Therefore, the structure along the c-axis is not considered.

4. Conclusions

In summary, a simple model Hamiltonian has been developed to explore the mechanisms of ion transport and non-adiabatic electron transfer in LiFePO₄. The model provides two analytical activation barriers as a function of lattice parameters for both symmetric and asymmetric ion hopping, respectively. The symmetric activation barrier is reduced to that of classical Marcus charge transfer theories by ignoring ion effects. Importantly, the model reveals that ion–electron coupled transfer is the optimal reaction path driven by the electrostatic attraction between Li⁺ and e⁻ (Fe²⁺). Separate ion transport or electron transfer rarely occurs due to the extremely large energy barrier, which originates from the electrostatic interaction energy required to eliminate an ion–electron pair. In addition, the model shows that Li⁺–Li⁺ repulsive interaction in the [010] direction and Li⁺–e⁻ (Fe²⁺)–Li⁺ attractive interaction in the [100] direction cause the phase segregation of LiFePO₄. The results verify that the interface of Li-rich and Li-poor phases is perpendicular to the [010] channel. This work provides a new perspective on designing high-rate performance electrode materials. One can adjust the electrostatic interaction to improve the ionic and electronic conductivity by doping, coating, and other techniques.

Data availability

The data that support the findings of this study are available on request from the corresponding author.

Author contributions

S. L. C. conceived the project. X. X. W., J. H. and S. L. C. derived the model equations. X. X. W. performed the DFT calculations. Y. W. L. discussed the derivation and results. X. X. W., J. H. and S. L. C. wrote the paper. All authors have given approval to the final version of the manuscript.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (22272122, 21673163, and 21832004). J. H. is supported by the Initiative and Networking Fund of the Helmholtz Association under the award number VH-NG-1709. We gratefully thank the Supercomputing Center of Wuhan University and Wuhan Supercomputing Center for generous computing resources to support our research.

Notes and references

1. J. B. Goodenough, Acc. Chem. Res., 2013, 46, 1053–1061.
2. J. B. Goodenough and Y. Kim, Chem. Mater., 2010, 22, 587–603.
3. B. Kang and G. Ceder, Nature, 2009, 458, 190–193.
4. Y. C. Hu, X. X. Wang, P. Li, J. X. Chen and S. L. Chen, Sci. China: Chem., 2023, 1–10.
5. Z. Ma, Z. J. Zuo, L. Li and Y. N. Li, Adv. Funct. Mater., 2022, 32, 2108692.
6. Y. Gao, J. Huang, Y. W. Liu and S. L. Chen, Chem. Sci., 2022, 13, 257–262.
7. K. Tang, X. Q. Yu, J. P. Sun, H. Li and X. J. Huang, Electrochim. Acta, 2011, 56, 4869–4875.
8. J. Xie, N. Imanishi, T. Zhang, A. Hirano, Y. Takeda and O. Yamamoto, J. Power Sources, 2009, 192, 689–692.
9. P. P. Prosini, M. Lisi, D. Zane and M. Pasquali, Solid State Ionics, 2002, 148, 45–51.
10. H. W. Zhang, J. Y. Li, L. Q. Luo, J. Zhao, J. Y. He, X. X. Zhao, H. Liu, Y. B. Qin, F. Y. Wang and J. J. Song, J. Alloys Compd., 2021, 876, 160210.
11. S. Heng, Q. Shi, X. Y. Zheng, Y. Wang, Q. T. Qu, G. Liu, V. S. Battaglia and H. H. Zheng, Electrochim. Acta, 2017, 258, 1244–1253.
12. Y. P. Xu and J. Mao, J. Mater. Sci., 2016, 51, 10026–10034.
13. H. Choi, J. Y. Seo and C. S. Kim, IEEE Trans. Magn., 2020, 57, 1–5.
14. O. M. Karakulina, N. R. Khasanova, O. A. Drozhzhin, A. A. Tsirlin, J. Hadermann, E. V. Antipov and A. M. Abakumov, Chem. Mater., 2016, 28, 7578–7581.
15. H. Q. Wang, A. J. Lai, D. Q. Huang, Y. Q. Chu, S. J. Hu, Q. C. Pan, Z. H. Liu, F. H. Zheng, Y. G. Huang and Q. Y. Li, New J. Chem., 2021, 45, 5695–5703.
16. H. C. Liu, Y. M. Wang and C. C. Hsieh, Ceram. Int., 2017, 43, 3196–3201.
17. M. A. Alsamet and E. Burgaz, Electrochim. Acta, 2021, 367, 137530.
18. L. Hong, L. S. Li, Y.-K. Chen-Wiegart, J. J. Wang, K. Xiang, L. Y. Gan, W. J. Li, F. Meng, F. Wang and J. Wang, Nat. Commun., 2017, 8, 1194.
19. T. Maxisch, F. Zhou and G. Ceder, Phys. Rev. B: Condens. Matter Mater. Phys., 2006, 73, 104301.
20. Y. Sun, X. Lu, R. J. Xiao, H. Li and X. J. Huang, Chem. Mater., 2012, 24, 4693–4703.
21. B. Ellis, L. K. Perry, D. H. Ryan and L. F. Nazar, J. Am. Chem. Soc., 2006, 128, 11416–11422.
22. R. Malik, A. Abdellahi and G. Ceder, J. Electrochem. Soc., 2013, 160, A3179.
23. G. H. Tao, J. Phys. Chem. C, 2016, 120, 6938–6952.
24. D. Fraggedakis, M. McEldrew, R. B. Smith, Y. Krishnan, Y. R. Zhang, P. Bai, W. C. Chueh, Y. Shao-Horn and M. Z. Bazant, Electrochim. Acta, 2021, 367, 137432.
25. X. X. Wang, J. Huang and S. L. Chen, Curr. Opin. Electrochem., 2023, 101358.
26. M. D. Newton and N. Sutin, Annu. Rev. Phys. Chem., 1984, 35, 437–480.
27. B. S. Brunschwig, J. Logan, M. D. Newton and N. Sutin, J. Am. Chem. Soc., 1980, 102, 5798–5809.
28. J. Franck and E. Dymond, Trans. Faraday Soc., 1926, 21, 536–542.
29. Z. M. Zhang, Y. Gao, S. L. Chen and J. Huang, J. Electrochem. Soc., 2019, 167, 013519.
30. Y. Gao, J. Huang, Y. W. Liu, J. W. Yan, B. W. Mao and S. L. Chen, Sci. China: Chem., 2019, 62, 515–520.
31. S. Ghosh, A. V. Soudackov and S. Hammes-Schiffer, J. Chem. Theory Comput., 2016, 12, 2917–2925.
32. B. S. Brunschwig and N. Sutin, Coord. Chem. Rev., 1999, 187, 233–254.
33. Q. F. Zhao, S. Q. Zhang, M. Y. Hu, C. Wang and G. H. Jiang, Int. J. Electrochem. Sci., 2021, 16, 1–10.
34. C. L. Gong, Z. G. Xue, S. Wen, Y. S. Ye and X. L. Xie, J. Power Sources, 2016, 318, 93–112.
35. S. E. Boulfelfel, G. Seifert and S. Leoni, J. Mater. Chem., 2011, 21, 16365–16372.
36. E. R. Logan and J. R. Dahn, J. Electrochem. Soc., 2021, 168, 120526.
37. K. Jalkanen, T. Aho and K. Vuorilehto, J. Power Sources, 2013, 243, 354–360.
38. J. N. Israelachvili, in Intermolecular and Surface Forces, ed. J. N. Israelachvili, Academic Press, San Diego, 3rd edn, 2011, pp. 53–70,  DOI:10.1016/B978-0-12-375182-9.10003-X.
39. R. A. Marcus, Rev. Mod. Phys., 1993, 65, 599–610.
40. R. A. Marcus, Annu. Rev. Phys. Chem., 1964, 15, 155–196.
41. R. A. Marcus, J. Chem. Phys., 1965, 43, 679–701.
42. R. A. Marcus, J. Chem. Phys., 1956, 24, 966–978.
43. S. L. Chen, Y. W. Liu and J. X. Chen, Chem. Soc. Rev., 2014, 43, 5372–5386.
44. E. Laborda, M. C. Henstridge, C. Batchelor-McAuley and R. G. Compton, Chem. Soc. Rev., 2013, 42, 4894–4905.
45. W. Schmickler, J. Electroanal. Chem. Interfacial Electrochem., 1986, 204, 31–43.
46. P. W. Anderson, Phys. Rev., 1961, 124, 41.
47. R. R. Dogonadze, A. M. Kuznetsov and T. A. Marsagishvili, Electrochim. Acta, 1980, 25, 1–28.
48. J. Huang, J. Chem. Phys., 2020, 153, 164707.
49. E. Laborda, M. C. Henstridge and R. G. Compton, J. Electroanal. Chem., 2012, 667, 48–53.
50. M. C. Henstridge, E. Laborda, Y. Wang, D. Suwatchara, N. Rees, Á. Molina, F. Martinez-Ortiz and R. G. Compton, J. Electroanal. Chem., 2012, 672, 45–52.
51. M. C. Henstridge, E. Laborda and R. G. Compton, J. Electroanal. Chem., 2012, 674, 90–96.
52. Y. Zeng, P. Bai, R. B. Smith and M. Z. Bazant, J. Electroanal. Chem., 2015, 748, 52–57.
53. S. P. Ong, V. L. Chevrier and G. Ceder, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 83, 075112.
54. K. Zaghib, A. Mauger, J. B. Goodenough, F. Gendron and C. M. Julien, Chem. Mater., 2007, 19, 3740–3747.
55. C. Cholsuk, S. Suwanna and K. Tivakornsasithorn, Mater. Today Commun., 2021, 26, 102043.
56. Y. Gu, M. Y. Weng, G. F. Teng, H. Zeng, J. S. Jie, W. J. Xiao, J. X. Zheng and F. Pan, Phys. Chem. Chem. Phys., 2019, 21, 4578–4583.
57. V. A. Dinh, J. Nara and T. Ohno, Appl. Phys. Express, 2012, 5, 045801.
58. T. W. Chen, Y. K. Ye, Y. Wang, C. Fang, W. C. Lin, Y. Jiang, B. Xu, C. Y. Ouyang and J. Zheng, Phys. Chem. Chem. Phys., 2023, 25, 8734–8742.
59. G. G. Xu, K. H. Zhong, J. M. Zhang and Z. G. Huang, J. Appl. Phys., 2014, 116, 063703.
60. H. Lin, Y. W. Wen, C. X. Zhang, L. L. Zhang, Y. H. Huang, B. Shan and R. Chen, Solid State Commun., 2012, 152, 999–1003.
61. Z. J. Liu and X. J. Huang, Solid State Ionics, 2010, 181, 907–913.
62. P. H. Xiao and G. Henkelman, ACS Nano, 2018, 12, 844–851.
63. F. Zhou, T. Maxisch and G. Ceder, Phys. Rev. Lett., 2006, 97, 155704.
64. Y. J. Zhu, J. W. Wang, Y. Liu, X. H. Liu, A. Kushima, Y. H. Liu, Y. H. Xu, S. X. Mao, J. Li, C. S. Wang and J. Y. Huang, Adv. Mater., 2013, 25, 5461–5466.

---

Footnote

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

---