Particle-and crack-size dependency of lithium-ion battery materials LiFePO

Lithium-ion batteries have become a widely-used commodity for satisfying the world’s mobile power needs. However, the mechanical degradation of lithium-ion batteries initiated by micro cracks is considered to be a bottleneck for advancing the current technology. This study utilizes a finite element method-based virtual crack closure technique to obtain particleand crack-sizedependent estimates of mixed-mode energy release rates and stress intensity factors. Interfacial cracks in orthotropic bi-materials are considered in the current study, whereas the crack extension along the interface is assumed. The results show that energy release rate, stress intensity factor, and the propensity of crack extension are particleand crack-sizedependent. In particular, our results show that for smaller plate-like LiFePO4 particles (100 nm  45 nm), a crack has lesser tendency to extend if crack-to-particle size is less than 0.2, and for 200 nm  90 nm particles, similar results are obtained for crack-to-particle sizes of less than 0.15. However, for larger particles (500 nm  225 nm), it requires an almost flawless particle to have no crack extension. Therefore, the current study provides insight into the fracture mechanics of LiFePO4 and the associated crack-to-particle size dependency to prevent crack extensions.


Introduction
Batteries have long been a very important commodity for satisfying the world's mobile energy storage needs.But more than ever before, advanced lithium-ion batteries are receiving an increased level of interest due to expanded power capabilities and the near limitless assortment of application opportunities, from vehicles and laptop computers to orbiting satellites and electric tools.As transportation markets are coming into play, the environmental impact of lithium-ion batteries and durability characteristics against climate variation are becoming very important in both the component material selection and design [1][2][3][4][5].In 1997, the concept of LiFePO 4 as a possible cathode for today's Li-ion batteries came to life [6].This material, which has an olivine structure, has been a focal point for much experimentation and discussion, as it boasts several attractive qualities including a relatively high theoretical capacity of 170 mAh/g, great structural stability, long cycle life, environmentally benign qualities, and most importantly, safety due to great thermal stability [6,7].However, LiFePO 4 has a very low electronic conductivity, a property essential for battery design, but it has been greatly improved through the successful application of conductive coatings or cation doping, making it comparable with other cathode materials [8][9][10][11].In addition to these attributes, the constituent materials are abundant and widely available at a relatively low cost [12], making LiFePO 4 extremely attractive for high-power vehicles, military and space flight operations [13,14], and everyday consumer electronics.
Phase transformation from lithium-poor (FePO 4 ) to lithium-rich (LiFePO 4 ) during intercalation induces different strains in each lattice direction: a, b, and c [15][16][17][18][19]. Stresses induced by the phase transformation create flaws or cracks that may become stress concentration areas when further cycling occurs [20][21][22][23][24][25][26].As a result, cracks are often observed along the phase boundary on the acplane [15,27]; structural failure due to fracturing of active material is a primary factor in battery degradation [7,22,[27][28][29][30].A combination of tensile and compressive stresses will induce microscopic crack propagation if critical values are reached.Over time, crack propagation will continue until the crack reaches a critical length, at which time complete fracture initiates.When this active electrode material is fractured away, the battery will effectively lose capacity and maximum power output diminishes.Thus, understanding the relationship between stresses and imperfections during lithium ion intercalation can be beneficial to better advance the current battery technology.
Olivine LiFePO 4 is a brittle material and linear elastic fracture mechanics could be applied to better understand mixed-mode fractures in lithium-ion battery materials.Extensive theoretical and computational studies have focused on the growth of interfacial cracks of bi-materials, and on quantifying the extension of a crack and the elastic field near a crack tip [31,32]; the energy release rates (G) and stress intensity factors (K 1 and K 2 ) of engineering materials are generally discussed in this regard.Hutchinson and Suo have provided analytical energy release rates and stress intensity factors for interfacial cracks in isotropic bi-materials under mixed-mode loading [33].Qian and Sun have extended the mixed-mode study to monoclinic and orthogonal bi-materials [34,35].Sequentially, the development and application of the finite element-based virtual crack closure technique (VCCT) facilitated the analyses on fracture mechanics for orthogonal (orthorhombic) bimaterials [36][37][38][39][40]. Raju has calculated energy release rates via the VCCT with higher order and singular finite elements for isotropic materials [41].Agrawal and Karlsson utilized the VCCT to obtain mixed-mode energy release rates and stress intensity factors for isotropic materials [40].For details about the overview of the VCCT, please refer to studies by Krueger et al., Rybicki, Dattaguru, Xie and Agrawal [36][37][38][39][40].
The current study aims to identify an optimal combination of the particle-and crack-size of LiFePO 4 for better battery performance.Because reducing the LiFePO 4 electrode particle size could effectively decrease stress inside the particles, in addition to reducing the diffusion path during electrochemical cycling [42,43].The current study is designed to provide insight into the sub-micro scale fracture mechanisms within the cathode material by accomplishing the following objectives to help aid in the understanding and development of a longer lasting and more powerful battery: (1) Perform a finite element simulation to evaluate fracture mechanisms inside a single LiFePO 4 particle due to the phase transformation; and (2) establish an energy-based approach to estimate the propensity of crack extension for lithium-ion battery materials.
We achieve these objectives by utilizing the VCCT to study fracture mechanics in LiFePO 4 battery materials under mixed-mode loading.Pre-existing cracks at the interface of orthotropic bimaterials are incorporated.Various particle and crack sizes are considered to illustrate the effects on the cracking information in LiFePO 4 due to phase transformation.Our aim is to present an approach towards improving battery structural design by providing a better understanding of the micromechanics, as well as a basis for future fatigue analyses, which may incorporate a variety of other battery chemistries and additional fatigue life parameters.

Model development
To study fracture mechanics in LiFePO 4 particles, VCCT with the ANSYS finite element software (ANSYS, Inc.Canonsburg, Pennsylvania) is utilized to obtain the energy release rates and stress intensity factors.Two-phase plate-like LiFePO 4 finite element models are generated based on the separate experimental observations [15,27,44] and different orthotropic material properties are adopted for both phases [28].An interfacial crack is considered to run parallel to the bc-plane (along the c-axis) where phase boundaries are present (Figure 1).Although crack extension is possible in other planes, literature and experimentation have reported that cracks are observed in the bc-plane between phases while advancing in the c-axis [15].In the current study, an assumption was made that the cracks will extend along their original direction.Moreover, the only experimental data that was available was for the surface energy and the prediction of crack turning or kinking was beyond the scope of the investigation.However, there exists a slight possibility that the shear mode fracture toughness may be very small in these materials and the crack path may deviate from a straight line.Therefore, we assume that a crack extends along the interface and remains parallel to the bc-plane [15,27,29] (Figure 1).Based on the experimental observation [15], previous simulations [45][46][47][48], and LiFePO 4 lattice constants of the unit cell-a = 10.3Å, b = 6.0 Å, and c = 4.7 Å [49]-three differently sized plate-like LiFePO 4 finite element models on the ac-plane are generated as follows: (a) 500 nm × 225 nm, (b) 200 nm × 90 nm, and (c) 100 nm × 45 nm, respectively.In general, a continuum mechanics-based finite element method is valid when the model size is above 25 nm, as several nano-scaled studies presented by Zhao et al. [22,50,51].Furthermore, lithium-ions diffuse into LiFePO 4 particles along the b-axis during the phase transformation; the elastic model is reduced to 2D by assuming plane stress and neglecting stresses in the b-direction, as described by Hu et al. [29].A total of 21 finite element simulations are conducted, wherein the initial crack is considered to start from the top face with variable lengths of crack sizes 0.05-0.8L/d, in which L represents the crack length and d represents the particle size along the caxis (Figure 1).The initial crack mouth opening is set at 0.5 nm in the a-axis in the range of experimentally observed crack dimensions [27].Strains were applied to the LiFePO 4 phase according to the misfit strains previously measured and reported in the literature: expansion occurs along the adirection (ε a = 5.03%) and contraction occurs along the c-direction (ε c = −1.9%)(Figure 1a) [52].Of note, ε b = 4.5%).Two-dimensional quadrilateral elements are used in the finite element analysis.Such eight-node parabolic elements allow for more flexibility and improved accuracy in contrast to four-node elements That is, the 8 noded parabolic quadrilateral element uses quadratic functions for the displacements, and the simulation shows the exact analytic solution for pure bending dominated problems even with a coarse mesh with only one element in depth [53,54].The finite element models are densely meshed around the crack tip to a size of 50 nm (Figure 1b).

Finite element method-based virtual crack closure technique (VCCT)
The VCCT is used along with the finite element modeling package ANSYS to determine the energy release rates and stress intensity factors under mixed-mode loading.The theory behind the VCCT is that the energy needed to separate or slide a crack surface is the same as the energy needed to close the surface back on itself [36].The tendency of a crack to extend, the energy release rates, and the stress intensity factors are calculated through command options in ANSYS [55] (Figure 2a), and they are currently not available on the graphical user interface of ANSYS.Stress intensity factors have been derived analytically for orthotropic bi-materials, as detailed in published studies [34,56,57].Energy release rates can be calculated with relative displacements (u, w), reaction forces (X, Z), and a unit thickness t = 1, as shown in Equation 1 and Figure 2 [39].From the VCCT and a few commands in ANSYS, the Mode 1 (tension mode) and Mode 2 (sliding shear mode) strain energy release rates for a flaw of given size in a two-phase LiFePO 4 particle could be determined (Figure 2b-c):

Crack extension
A crack will advance if the total energy release rate for Mode 1 and Mode 2, G T = G 1 + G 2 , is larger than approximately twice the surface energy of the particle, i.e., G T > 2 [29,58].A firstprinciple analysis by Wang et al. has reported a  value of 0.66 N/m for LiFePO 4 in the (100) crack face orientation [29,58].In the current study, we collect finite element results satisfying G T > 2, and the difference between these two energies (G T − 2) is used to predict crack extension.That is, crack extension is eminent and a crack will advance until the surface energy from the newly-formed crack faces bring the particle back to be equal to or larger than the strain energy release rate (Figure 3).In the current work, we follow steps in "Linear Elastic Fracture Mechanics" detailed by Fischer-Cripps [59]: for a crack to extend, the rate of strain energy release per unit of crack extension must be at least equal to the rate of surface energy requirement [60].Therefore, a 0.01-nm crack extension is set as a unit growth and the resulting crack growth (da) would be determined as follows: First, we calculate the difference between the G T and 2.Second, the additional energy release rate for a unit of 0.01 nm crack extension is determined based on the crack surface areas.Third, the updated additional surface energy for a 0.01 nm crack extension is determined.Fourth, the crack extension is calculated since the relieved internal energy due to crack propagation should be balanced by the additional strain energy release rate, as shown in Figure 3.

Results and Discussion
The results of the finite element method-based VCCT indicate that normal stresses dominate as compared to shear stress throughout the particle (stress fields not shown).In particular, the Mode 1 fractures predominated for all models, as the energy release rate G 1 is much larger than G 2 for all particle sizes (Figure 4).The fracture mode is a result of applied boundary conditions (misfit strains due to phase transformation in the current study: ε a = 5.03% and ε c = −1.9%),anisotropic elastic moduli for LiFePO 4 and FePO 4 [28], and as a result, the energy release rates and the stress intensity factors are indicators of the state of energy or stress around the crack tip.Since the stiffness value of LiFePO 4 from [28] along the a-axis is much smaller than in the other orthorhombic directions, the material is weaker in this direction, giving way to the splitting or opening Mode 1 type of fracture.Therefore, in the current study we focus on the discussion of Mode 1 fractures in LiFePO 4 battery materials from the mixed-mode boundary condition.
Several studies have shown that the energy release rate is proportional to the crack size since the energy release rate is the energy dissipated during fracture per unit of newly created fracture surface area [33].The variation of energy release rates with respect to the crack-to-particle ratios for two fracture modes are shown in Figure 4.It is revealed that G 1 is highly dependent on the crack size for all three particle sizes.Comparing to the 2 = 1.32 N/m, smaller particles (100 nm  45 nm) are able to better accommodate initial flaws without the crack propagation, i.e., L/d ≤ 0.2.In contrast, particles with a larger size (500 nm  225 nm) could only accommodate initial flaws as small as L/d ≤ 0.03 by a linear extrapolation.Therefore, it is the L/d point that needs to be considered to determine if crack propagation occurs during the lithiation/delithiation process.Below L/d = 0.5 and 0.6, the energy release rates increase with the increased crack lengths and then begin to level off.It is concluded that crack propagation acts to relieve internal stresses due to misfit strain during the phase transformation of LiFePO 4 materials, and it is particle-size and initial flaw-size dependent.Crack propagation due to mixed-mode loading is possible in other planes; however, experimentation has reported that cracks are observed between phases while propagating in the c-axis [15].The olivine crystal structure of LiFePO 4 and the Pauling's third rule support this observation: (100) planes are linked through stable, shared corner bonds between FeO 6 octahedra and PO 4 tetrahedra, while (010) planes are connected by weaker, edge-sharing bonds between FeO 6 octahedra [27,52,61].Lower energy release rates are observed in the Mode 2 fracture and exhibit less dependency on particle sizes.However, G 1 and G 2 both contribute to the G T for the crack extension.In contrast to the study by Hu et al. [29], a maximum value of the normalized energy release rate was identified; however, it is not clear which fracture mode is considered in their study.The approach in their study was to inversely look for a critical particle size d by solving boundary value problems, and the normalized energy release rate is a function of one component of the anisotropic stiffness matrix (c 11 ) and the volume misfit along the c-axis ( c ).The decreasing normalized energy release rate after L/d = 0.5 from Hu et al. could be attributed by several simplifications [29], such as the choice of the normalization where the other anisotropic effects of volume expansions were not incorporated.
Stress intensity factors (K) are also obtained from the command-coding interface of ANSYS via the VCCT (Figure 2a), and the particle-size-dependent relationship of stress intensity factors and L/d is shown in Figure 5a.For larger particles (500 nm  225 nm), it is observed that K 1 reaches a relative maximum value of 1.4 MPa-m 1/2 .For a particle size of 200 nm × 90 nm, K 1 reaches a relative maximum value of 1.0 MPa-m 1/2 .For smaller particles (100 nm × 45 nm), it is observed that K 1 reaches a relative maximum value of 0.6 MPa-m 1/2 (Figure 5a).Our results indicate that K 1 is particle-size dependent, which is consistent with the observations made by Krstic and Khaund [62].Note that these simulations are not based on Paris' law, but rather from systematic analyses of independent VCCT models (i.e., 16 models are solved to generate the 16 data points plotted above).It is observed that smaller particles exhibit faster crack propagation (i.e., a higher value for the exponent of ΔK 1 in the best-fit equations appearing in the plot above).
Further, to reveal the relationship between the crack advancement (da) and the change in the Mode 1 stress intensity factor for LiFePO 4 , the estimation is depicted in Figure 5b.The da is determined via the method described in Figure 3, and ΔK 1 is determined from ranges of the stress intensity factor of Mode 1.It is observed that da vs. ΔK 1 curves are particle-size dependent, wherein smaller particles exhibit faster crack extension (i.e., a higher value for the exponent of ΔK 1 in the best-fit equations appearing in Figure 5b): the crack growth is proportional to the K 1 to the power of 4.6229 for particle size of 100 nm × 45 nm, and for particle size of 500 nm  225 nm, the crack growth is related to the K 1 to the power of 1.3849 (Figure 5b).Detailed calculated data from the FEM-based VCCT are listed in Table 1.
Table 1.Results from the FEM-based VCCT of lithium-ion battery materials for each particle size with various crack-to-particle ratios: Mode 1 stress intensity factor (K 1 ), total energy release rate (G T ), twice of the surface energies (2γ), the energy difference between them (ΔE = G T -2γ), the required additional strain energy release rate per 0.01 nm crack growth (E R ), the required additional strain energy release rate per 0.01-nm crack growth (E s ), and the required crack advancement (da).Though the particle sizes chosen for the study retain the same aspect ratio, it is believed the surface energy has greater effects for smaller particles [63], where the same surface energy for three particle-size models are adopted, i.e.,  = 0.66 N/m [58].That is, for a particle with a smaller size, the newly created surface area has a stronger effect since crack propagation in a smaller particle acts to relieve more internal stresses.As a result, less strain energy release rates and stress-intensity factors are found in smaller particles (Figures 4 and 5a).In addition, a convergence test was performed in our finite element analysis.It has been observed that models beyond L/d = 0.6 provide similar strain energy release rates and stress intensity factors, suggesting that L/d = 0.6 is the critical value for both computational simulations and battery materials design.
In general, an empirical Paris' Law is used to describe crack propagation in engineering materials.However, fatigue experiments are not currently available for nanosized LiFePO 4 particles.The current study specifically aims to provide a method that estimates crack advancement based on the effects of particle and crack sizes under mixed-mode loading.The current study does not aim to estimate material parameters from the Paris' Law for LiFePO 4 particles.Our approach is in contrast to that of Deshpande et al. [64], who utilized Paris' Law and incorporated material parameters for LiFePO 4 to predict battery life (i.e., numbers of cycles).In their study, an isotropic core-shell model incorporating chemical and mechanical degradation was implemented, considering a surface crack that grows with each charge-discharge cycle [64].
The rate of discharge will determine the resultant stresses within the material, which will in turn alter the energy release rates and stress intensity factors at a crack in the particle.Therefore a currentrate-dependent crack propagation analysis for LiFePO 4 is currently being developed in our research group.In general, larger stresses over a shorter period of time may induce cracking more readily than for slower discharging rates; this is due to the limited volumetric expansion rate being unable to satisfying the rate of lithium diffusivity during fast discharging, similar to rush hour traffic in limited highway space.Zhu et al. [20] have reported current density dependent fracture mechanic characteristics for another battery chemistry, LiMn 2 O 4 , wherein isotropic spherical and ellipsoidal particles are considered.They have concluded that crack propagation occurs at the center of the particle even during relatively low current density discharging [20].
Molecular dynamic (MD) simulation has been a powerful tool to numerically solve the classical equations of motion for a system at a smaller-scale.It has been widely used to understand the behavior of liquid electrolytes for lithium-ion batteries to improve the performance of portable electronics devices, electric vehicles, and hybrid electric vehicles.However, such considerations represent limitations of our study, as we did not have the computational tools in the present study.We have chosen to conduct a continuum mechanics-based finite element analyses based available experimental and computational results from the literatures in which particles and models are on the same size order of ours.In this case, we could use the available surface energy values from the literature to calculate the crack extension, as shown in Table 1.Moreover, since our current study mainly focuses on interface-diffusion limited lithium-ion battery materials, an improved model by using MD could be developed to better predict interface diffusion and failure.Another limitation in the current study lies on the small element size at the crack tip in which a MD simulation is regarded as a better choice.To this matter, the logical approach in the desire to create truly multiple scale simulations that exist at disparate length and timescales has been to couple MD and FE.Unfortunately, the coupling of these methods is not straightforward since the major problem in multiscale simulations is that of pathological wave reflection, which occurs at the interface between the MD and FE regions.We did not have the coupled MD and FE computational tools in the present study.Therefore, the idea of meshing the finite element region down to a small-scale was one of our first attempts to eliminate spurious wave reflections, which has been generally reported at the MD-FEM interface.To this end, we are in the process by adopting MD simulations to improve our computational technique.

Conclusions
To the best of the authors' knowledge, this is the first study utilizing the finite element methodbased VCCT to obtain cracking information in LiFePO 4 battery materials.Various lengths of interfacial cracks in orthotropic bi-materials are incorporated.A method is established to calculate the amount of crack propagation for lithium-ion battery materials (Figure 3 and Figure 5b).It is observed that energy release rates, stress intensity factors, and crack propagation are all particle-and crack-size-dependent (Figures 4-5).Therefore, one potential future consideration is to reduce the LiFePO 4 electrode particle size to allow initial flaws during material processing (Figure 4) and to effectively decrease stress inside the particles (Figure 5a) in addition to reducing the diffusion path for better battery performance [42,43].
The current study provides an approach towards improving battery structural design by providing a better understanding of the micromechanics, as well as a basis for future fatigue analyses, which may incorporate a variety of other battery chemistries and additional fatigue life parameters.Thus, the current study has identified factors that play an important role in inducing fracture, and the study also provides insight into factors that can be optimized to minimize any detrimental capacity loss in a battery.This research will hopefully elucidate a relationship between micro-mechanics and battery usage so as to help in the design of a higher performance, longer lasting battery for the future.

Figure 1 .
Figure 1.(a) A representative orthotropic bi-material with variable lengths of crack sizes 0.05-0.8L/d.Expansion occurs along the a-direction and contraction occurs along the cdirection when particles undergo phase changes [52].(b) Two-dimensional finite element models are used in the analysis, and models are densely meshed around the crack tip to a size of 50 nm.Scale bar = 0.03 μm.

Figure 2 .
Figure 2. (a) The VCCT for 2-D quadrilateral elements with a unit thickness t.Energy release rates are calculated based on the crack size (Δa), reaction forces (Z and X), relative displacements (u, w, and a), and via the FEM-based VCCT from the commandcoding interface of ANSYS.(b)-(c) Anticipated Mode 1 (tension mode) and Mode 2 (sliding shear mode) fractures for a two-phase LiFePO 4 battery material under mixedmode loading conditions.

Figure 3 .
Figure 3. Crack extension calculation.The argument G T > 2γ is used to determine the required crack advancement (da) to reach mechanical equilibrium, where b is the particle size in the b-axis and  = 0.66 N/m [29,58].

Figure 4 .
Figure 4. Energy release rates for two modes of fracture at crack tips with various L/d ratios.G 1 is highly dependent on the crack size for three different particle sizes.It is observed that Mode 1 cracks are likely to occur and has higher tendency to propagate, especially when L/d ≈ 0.5-0.6.

Figure 5 .
Figure 5. (a) Mode 1 stress intensity factors at crack tips with various L/d ratios.K 1 is highly dependent on the particle and crack sizes.(b) Computational model predictions of crack advancement (Mode 1) for three different particle sizes.Note that these simulations are not based on Paris' law, but rather from systematic analyses of independent VCCT models (i.e., 16 models are solved to generate the 16 data points plotted above).It is observed that smaller particles exhibit faster crack propagation (i.e., a higher value for the exponent of ΔK 1 in the best-fit equations appearing in the plot above).