Towards Determining an Engineering Stress-Strain Curve and Damage of the Cylindrical Lithium-Ion Battery Using the Cylindrical Indentation Test

Abstract

Mechanical internal short circuit (ISC) is one of the significant safety issues in lithium-ion battery design. As a result, it is possible to subject LIB cells to thorough mechanical abuse tests to determine when and why failure may occur. The indentation test is a recommended loading condition for evaluating mechanical damage and ISC. In this study, 18,650 cylindrical battery cells underwent indentation tests and a voltage reduction following the peak force identified by the ISC. Due to the complexity of the contact surface shape between two cylinders (LIB cell and indenter), a new phenomenological analytical model is proposed to measure the projected contact area, which the FEM model confirms. Moreover, the stress-strain curve and Young’s modulus reduction were calculated from the load-depth data. In contrast to previously published models, the model developed in this paper assumes anisotropic hyperelasticity (the transversely isotropic case) and predicts the growing load-carrying capacity (scalar damage), whose variation is regulated by the Caputo-Almeida fractional derivative.

1. Introduction

In various applications, such as electric vehicles and cell phones, lithium-ion battery (LIB) cells are considered the first choice for power batteries due to their higher capacity and higher efficiency [1]. However, the safety of electric vehicles has become a crucial issue since LIB cells can be subjected to various mechanical loading conditions [2,3]. Hence, it is essential to have a better knowledge of the damage behavior of LIB cells, which has drawn extensive attention from vehicle companies and engineering groups [4]. Furthermore, sufficient understanding in this critical field is necessary to develop reliable methods to ensure the safe application of the LIB cell under different mechanical loadings [4].

However, in terms of shape, four types of LIB cells are found in the literature, including pouch [5,6,7,8,9,10,11,12], elliptical [5,13,14,15], cylindrical [16,17,18,19,20], and prismatic [21,22,23]. Furthermore, regarding the detection of failure, various mechanical abuse tests with a range of indenters, such as nail [24,25,26,27,28,29], conical [30], hemispheric [5,6,11,19,23,30,31,32,33,34,35], cylindrical [17,31,32,33,36,37,38], flat [10,16,17,31,36,39], pinch [11], and three-point bending indenters [31] were used to induce the internal short circuit (ISC) on the whole LIB cell. These tests have been developed to address the safety concerns associated with LIB cells: the risk of failure, ISC, and thermal runaway. Consequently, many attempts have been made to develop a failure model based on a function of the plastic strain to predict the failure onset and propagation in cases of mechanically abusive loadings [31,39,40,41,42]. In contrast, the current study focuses on developing a damage model to predict material degradation through mechanical deformation.

Three primary modeling methods have been employed in the literature to evaluate the LIB cell’s response and safety risks following mechanical deformation and damage. The first method considers the LIB cell as a multi-layer structure and evaluates the mechanical properties of individual components [43]. The second approach utilizes the reduced-detail model by combining several layers into only a single one [9]. Finally, the last method considers the LIB cell as a homogenized medium with a significantly lower computational cost [44]. It should be pointed out that a homogenized model was developed in this study, which boasts higher computational efficiency than others.

However, in the literature, three isotropic, macro-homogenized material models have been developed to characterize the elastoplastic response of the LIB cell under different mechanical scenarios: power-law [41,45], hyper-foam [30,46,47], and crushable foam [30,35,46,47]. Furthermore, it should be noted that strain failure criteria have been used to predict the failure of the LIB cell in the above studies. Moreover, no explicit research is available in the literature addressing the degradation of mechanical properties due to mechanical deformation. Hence, in this study, a scalar damage model for hyperplastic materials is developed to predict the mechanisms and the extent of damage to the LIB cells in the case of the indentation test.

From the continuum damage mechanics perspective, it is widely accepted that changes in the macroscopic material properties (such as stiffness and Young’s modulus) can indicate the onset, growth, and combination of microcracks and microvoids [48,49]. Consequently, characterizing Young’s modulus of LIB cell under mechanical abuse deformation is crucial to understanding how this degradation of mechanical properties is related to an appropriate damage function. Therefore, the reduction in Young’s modulus in the computational continuum damage model is evaluated to indicate the damage onset and propagation through mechanical deformation.

Ali et al. [46] proposed a homogenized model based on the hyper-foam material model, using the strain energy function to represent the LIB cell’s mechanical behavior under in-plane constrained compression tests. While such models may be used widely to predict elastoplastic behavior, they cannot predict damage in the LIB cell. To address this, the current study evaluates LIB cells’ response under mechanical abuse tests using a combination of anisotropic hyperelasticity and the time-fractional damage evolution law within a phenomenological framework. The anisotropy is introduced through the fundamental rule that the effective (undamaged) Helmholtz free-energy function ($\tilde{\mathsf{\Psi }}$) is an isotropic tensor function, including the deformation gradient ($\mathrm{F}$) and a structural tensor ($\mathrm{M}$). It should be noted that physically observed damage was modeled using a scalar variable ($\mathsf{\varphi }$), which was isotropic and whose evolution was defined based on the Perzyna model [50]. Moreover, it introduces the rate dependence and implicit space-length scale (through relaxation time). Additionally, damage evolution is governed by using the Caputo-Almeida fractional derivative operator (${}_{\mathrm{t}-{\mathrm{l}}_{\mathrm{t}}}{}^{\mathrm{C}}\mathrm{D}_{\mathrm{t}}^{\mathsf{\alpha },\mathsf{\Psi }}$) Almeida [51], which incorporates a memory (time-length scale).

Since indentation is a versatile experimental test that can easily be used to evaluate the LIB cell’s response to plastic deformation and damage propagation, an indentation test was employed in this study. Furthermore, the present study developed a novel approach to capture the reduction in the Young’s modulus due to damage in the LIB cell and convert the indentation load-depth data into the stress-strain relationship. Initially, the sink-in and elastic deflection were taken into account to measure the reduced Young’s modulus and the projected contact area. Subsequently, researchers presented an effort to modify the reported hyperelastic model in the literature to make it more suitable for the LIB cell by including softening directly in the expression for damage accumulation.

2. Materials and Methods

The LIB samples used in this study were 18,650 cylindrical LIB cells with a voltage of 3.52 V, which were 65 mm in length and 18 mm in diameter. The experimental indentation tests were performed using electromechanical MTS Insight with 100 kN load cell. For the indentation tests, a rigid cylindrical steel rod was used as the indenter with a 16 mm diameter, and LIB cells were placed on the lower plate (as seen in Figure 1). A temporary holder was used to ensure the first connection between the cell and the indenter was created, which minimized the misalignment and rotation in the battery cell. Five maximum penetration depths (2, 3, 4, 5, and 6 mm) were used to measure the reduction in the Young’s modulus through the indentation test. The LIB cell was subjected to a 2 mm displacement, the second cell to a 3 mm displacement, and so on, and the loading speed was 1 mm min⁻¹ to extract the stress-strain relation from force-depth data. In addition, it is important to note that the indenter load direction and the cell seam are parallel, which significantly influences the experimental setup and the subsequent observations reported in this work. However, the in-situ voltage of LIB cells during indentation tests was measured by an Agilent 34410A Digital Multimeter. This study uses the voltage technique to predict the internal short circuit (ISC) and detects the force drop that coincides with when the voltage drops to 0 V.

3. Results and Discussion

3.1. Measurement of Young’s Modulus Reduction

From the homogenized model respective, the tensile or axial compression test was used to obtain Young’s modulus values for the LIB cell [31]. In the literature [33,36], a constant value, 1.5 MPa, was reported as Young’s modulus of the jellyroll. However, this definition of Young’s modulus cannot describe the elastic modulus reduction due to the damage accumulation. Therefore, this study employed the indentation test to measure the reduction of Young’s modulus through the indentation test.

A multistep indentation test was performed to measure Young’s modulus at various depths. This section performed single-loading indentation tests at five maximum penetration depths (2, 3, 4, 5, and 6 mm). The multistep indentation test provides more information about the material properties of LIB cells in a shorter time by giving information about load/reload at each step. However, as suggested by Doerner-Nix [52], a linear fit function was used for the initial 20% of the unloading segment to evaluate Young’s modulus during the indentation test.

According to Hertz’s theory of contact mechanics between two isotropic elastic solids [52], from the unloading part, the relation between the load (P) and the penetration depth (h) during the indentation test is defined as:

(1)$\mathrm{P}=\frac{4}{3}{\mathrm{E}}_{\mathrm{eff}}{\mathrm{R}}_{\mathrm{eff}}^{\frac{1}{2}}{\mathrm{h}}_{\mathrm{e}}\sqrt{{\mathrm{h}}_{\mathrm{e}}}$

where ${\mathrm{h}}_{\mathrm{e}}$ is elastic penetration depth. Effective Young’s modulus (${\mathrm{E}}_{\mathrm{eff}}$) and effective radius (${\mathrm{R}}_{\mathrm{eff}}$) are obtained as follows:

(2)${\mathrm{E}}_{\mathrm{eff}}=\frac{{\mathrm{E}}_{\mathrm{c}}{\mathrm{E}}_{\mathrm{i}}}{\left(1-{\mathrm{v}}_{\mathrm{c}}\right){\mathrm{E}}_{\mathrm{i}}+\left(1-{\mathrm{v}}_{\mathrm{i}}\right){\mathrm{E}}_{\mathrm{c}}}$

(3)${\mathrm{R}}_{\mathrm{eff}}=\frac{{\mathrm{R}}_{\mathrm{c}}{\mathrm{R}}_{\mathrm{i}}}{{\mathrm{R}}_{\mathrm{c}}+{\mathrm{R}}_{\mathrm{i}}}$

where ${\mathrm{E}}_{\mathrm{c}}$,${ \mathrm{v}}_{\mathrm{c}}$, ${\mathrm{R}}_{\mathrm{c}}$ and ${\mathrm{E}}_{\mathrm{i}}$, ${\mathrm{v}}_{\mathrm{i}}$, ${\mathrm{R}}_{\mathrm{i}}$ denote Young’s modulus, Poisson’s ratio, the LIB cell’s radius, and the corresponding indenter’s properties. Then, based on the linear approach to measuring an appropriate fit of the unloading segment, as seen in Figure 2, the contact stiffness, S, is experimentally obtained, which can also be analytically calculated by taking the first derivative of Equation (1):

(4)$\mathrm{S}=\frac{\mathrm{dP}}{\mathrm{dh}}={\mathrm{E}}_{\mathrm{eff}}{\left({\mathrm{R}}_{\mathrm{eff}}{\mathrm{h}}_{\mathrm{e}}\right)}^{\frac{1}{2}}$

Therefore, the reduced Young’s modulus can be expressed in terms of the contact stiffness and the indentation depth obtained by combing Equations (1) and (4) as:

(5)${\mathrm{E}}_{\mathrm{ff}}=\frac{\mathrm{S}}{\sqrt{{\mathrm{R}}_{\mathrm{eff}}{\mathrm{h}}_{\mathrm{e}}}}=\frac{\mathrm{S}}{{\mathrm{A}}_{\mathrm{c}}}$

where ${\mathrm{A}}_{\mathrm{c}}$ is the projected contact area, and an analytical approach will be introduced to calculate its value in the following section. The above equation confirms that Young’s modulus involves the projected contact area and the unload’s slope at the maximum depth point. Figure 2 shows that the unloading slope at the maximum depth point increased as the depth increased. However, simultaneously, the value of the projected contact area was raised. This test indicated that Young’s modulus decreased as the depth increased, which can result from damage accumulation. In the literature [35,53], a constant value used for Young’s modulus is most often considered as the LIB cell under mechanical abuse test in the computational simulations. However, in this study, a novel way was used to combine a couple of the concepts of the continuum damage mechanic and the experimental indentation test to show the degradation of Young’s modulus with increasing penetration depth, as seen in Figure 2. The rationale for using the reduction of Young’s modulus is to better describe the effect of the voids and cracks on the decrease in mechanical properties of the LIB cell. An analytical damage model will be introduced in the following sections.

3.2. Analytical Approach to Convert Force-Depth into Stress-Strain Relation

One of the challenges in the indentation test is determining the projected contact area between two crossed cylinders (the indenter and LIB cell). In this study, an analytical approach was developed to calculate the projected contact area as a function of the indenter depth. As the experimental results show in Figure 3a, the deformed area has an elliptical shape. According to Hertz’s theory, Popov [54] proposed a model to measure the contact area between two perpendicular cylinders with a radius ${\mathrm{R}}_1 {\mathrm{and} \mathrm{R}}_2$ (as seen in Figure 3b) to calculate the semi-axes of an ellipse as follows:

(6)$\mathrm{a}=\sqrt{{\mathrm{R}}_1{\mathrm{h}}_{\mathrm{c}}}$

(7)$\mathrm{b}=\sqrt{{\mathrm{R}}_2{\mathrm{h}}_{\mathrm{c}}}$

where ${\mathrm{h}}_{\mathrm{c}}$ is the contact depth being derived as a function of the maximum indenter depth (${\mathrm{h}}_{\max })$ in this study, as follows:

(8)${\mathrm{h}}_{\mathrm{c}}=0.75{\mathrm{h}}_{\max }$

Hence, the projected contact area (${\mathrm{A}}_{\mathrm{c}}$) can be obtained as:

(9)${\mathrm{A}}_{\mathrm{c}}=\mathsf{\pi }\mathrm{ab}={\mathsf{\pi }\mathrm{h}}_{\mathrm{c}}\sqrt{{\mathrm{R}}_1{\mathrm{R}}_2}$

Moreover, a finite element method (ABAQUS) is used to track deformation during the indentation test. To obtain the 3D cloud point data from the 3D mesh, a custom-made Python script was developed to read/extract data from the “odb” file and obtain the node depth/position from the contact area between the two crossed cylinders, as shown in Figure 3b. Hence, the projected contact measured from the FE model at successive indentation depth stages can verify the analytical solution for measuring the projected contact area, as shown in Figure 4.

Figure 5 shows that the cell voltage drops at a depth value of approximately 6.2 mm in indentation testing. Noting that the drop in voltage and peak force occur simultaneously confirms the presence of ISC. Moreover, the main objective of this study is to characterize the elastoplastic response of the LIB cell under the indentation test. Many attempts have been made in the literature to determine stress-strain curves of the LIB cell from axial and radial compression tests [13,16,41,44,55]. Since the indentation test is one of the typical mechanical abuse tests to trigger ISC, it is challenging to calculate the stress-strain curve from the force depth data from the indentation test rather than the axial compression test. Hence, in this section the projected contact area, measured in the previous section, is used to propose a new analytical method presented in this section to calculate the stress-strain relation. From the analytical point of view, one defines effective stress (${\mathsf{\sigma }}_{\mathrm{eff}}$) as:

(10)${\mathsf{\sigma }}_{\mathrm{eff}}=\frac{1}{\mathsf{\alpha }}\frac{\mathrm{P}}{{\mathrm{A}}_{\mathrm{c}}}=\frac{1}{\mathsf{\alpha }}\frac{\mathrm{P}}{{\mathsf{\pi }\mathrm{h}}_{\mathrm{c}}{\mathrm{R}}_{\mathrm{eff}}}$

where $\mathsf{\alpha }$ is the stress constraint factor.

In contrast to the previous studies in the literature [56,57], due to the difficulties involved in calculating the contact area, it came to the attention of the authors to define the indentation strain in novel ways as follows:

(11)${\mathsf{\epsilon }}_{\mathrm{xx}}=\frac{\mathrm{a}}{{\mathrm{R}}_{\mathrm{eff}}}$

(12)${\mathsf{\epsilon }}_{\mathrm{yy}}=\frac{\mathrm{b}}{{\mathrm{R}}_{\mathrm{eff}}}$

(13)${\mathsf{\epsilon }}_{\mathrm{zz}}=\frac{{\mathrm{h}}_{\mathrm{c}}}{{\mathrm{R}}_{\mathrm{eff}}}$

The effective indentation strain can be calculated as follows:

(14)${\mathsf{\epsilon }}_{\mathrm{eff}}=\frac{1}{\mathsf{\beta }}\sqrt{\frac{3}{2}\left({\mathsf{\epsilon }}_{\mathrm{xx}}^2+{\mathsf{\epsilon }}_{\mathrm{yy}}^2+{\mathsf{\epsilon }}_{\mathrm{zz}}^2\right)}$

where $\mathsf{\beta }$ is the strain constraint factor. Figure 6 shows the effective stress-strain curve obtained from Equations (10) and (14).

3.3. Fundamental Assumptions

In the simulation of the LIB cell under the indentation test, it is assumed that the conservation of mass, the balance of momentum, the balance of moment of momentum, and the balance of energy and entropy production hold. The first assumption states that the free energy function was modelled as a scalar function of the tensorial argument [58]:

(15)$\mathsf{\Psi }=\mathsf{\Psi }\left(F, M;\mathsf{\mu }\right)$

where $F$ is the deformation gradient, $M$ denotes the structural tensor (symmetric and positive definite), and $\mathsf{\mu }$ represents a set of internal state variables, which consists of the damage variable $\left(\mathsf{\varphi }\right)$. Hence, the degradation of mechanical properties of the LIB cell is described using a single scalar damage parameter, as follows:

(16)$\mathsf{\mu }=\left\{\mathsf{\varphi }\right\}$

As usual in damage mechanics, the variable $\mathsf{\varphi }$ takes the value within limits $\left[0,1\right]$. Moreover, it should be noted that the damage variable has a phenomenological meaning, so that $\mathsf{\varphi }=0$ is equivalent to the undamaged material, and $\mathsf{\varphi }=1$ means the complete loss of load-carrying capacity (loss of continuity). Next, one can rewrite the free energy function Equation (16) as:

(17)$\mathsf{\Psi }=\mathsf{\Psi }\left(F, M;\mathsf{\mu }\right)=\left(1-\mathsf{\varphi }\right)\tilde{\mathsf{\Psi }}\left(F,M\right)$

where $\tilde{\mathsf{\Psi }}\left(F,M\right)$ is the effective (undamaged) Helmholtz free-energy function. Furthermore, based on LIB geometry and its internal structure, we assume that anisotropic behavior is of transverse isotropic type [58], therefore:

(18)$\tilde{\mathsf{\Psi }}=\widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}=\widehat{{\mathsf{\Psi }}^{\mathrm{ti}}}\left({\mathrm{I}}_1,{\mathrm{I}}_2,{\mathrm{I}}_3,{\mathrm{J}}_4,{\mathrm{J}}_5\right)={\mathsf{\Psi }}^{\mathrm{iso}}\left({\mathrm{I}}_1,{\mathrm{I}}_2,{\mathrm{I}}_3\right)+{\mathsf{\Psi }}^{\mathrm{ti}}\left({\mathrm{J}}_4,{\mathrm{J}}_5\right)$

where

(19)${\mathrm{I}}_1=\mathrm{tr}C,{\mathrm{I}}_2=\mathrm{tr}\left[\mathrm{Cof}C\right],{\mathrm{I}}_3=\det C,{\mathrm{J}}_4=\mathrm{tr}\left[CM\right],{\mathrm{J}}_5=\mathrm{tr}\left[\mathrm{Cof}\left[C\right]M\right]$

where ${\mathsf{\Psi }}^{\mathrm{iso}}$ stands for the isotropic part, ${\mathsf{\Psi }}^{\mathrm{ti}}$ stands for transverse isotropy term, $\mathrm{C}$ denotes the right Cauchy-Green tensor, $\mathrm{Cof}\left[C\right]=\det \left[C\right]C^{-1}$ and $\mathrm{M}=\mathrm{diag}\left(\mathrm{a},\mathrm{a},\mathrm{b}\right)$ (${\mathrm{X}}_3$ spatial direction as chosen as a preferred one—axonal structure direction).

Under the above assumptions, the constitutive relation for the HTFD model is, in general

(20)$S=2\left(1-\mathsf{\varphi }\right)\left[\left(\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{I}}_1}+\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{I}}_2}{\mathrm{I}}_1\right)1-\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{I}}_2}C+\left(\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{I}}_3}{\mathrm{I}}_3+\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{J}}_5}{\mathrm{J}}_5\right){\mathrm{C}}^{-1}+\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{J}}_4}M-\frac{\partial \widetilde{{\mathsf{\Psi }}^{\mathrm{ti}}}}{\partial {\mathrm{J}}_5}{\mathrm{I}}_3{C}^{-1}M{C}^{-1}\right]$

where $S$ denotes the 2nd Piola-Kirchhoff stress tensor, an explicit definition of $S$ requires a specific form of the effective (undamaged) Helmholtz free-energy function, and the meaning of damage evolution—experimental data inspire both [59].

3.4. Free-Energy Function and Damage Evolution

Based on experimental evidence, the isotropic part of the effective (undamaged) Helmholtz free-energy function is assumed to be governed by the compressible Mooney- Rivlin model, namely

(21)${\mathsf{\Psi }}^{\mathrm{iso}}\left({\mathrm{I}}_1,{ \mathrm{I}}_2,{ \mathrm{I}}_3\right)={ \mathsf{\alpha }}_1{\mathrm{I}}_1+{\mathsf{\alpha }}_2{\mathrm{I}}_2+{\mathsf{\delta }}_1{\mathrm{I}}_3-{\mathsf{\delta }}_2\ln \left(\sqrt{{\mathrm{I}}_3}\right)$

where ${\mathsf{\alpha }}_1,{\mathsf{\alpha }}_2,{\mathsf{\delta }}_1,{\mathsf{\delta }}_2\ge 0$ denotes material parameters. Next, the transverse isotropic term has the form proposed by Schroder et al., (2008)

(22)${\mathsf{\Psi }}^{\mathrm{ti}}=\frac{1}{{\mathsf{\alpha }}_4{\left(\mathrm{trM}\right)}^{{\mathsf{\alpha }}_4}}{\mathsf{\eta }}_1\left({\mathrm{J}}_4^{{\mathsf{\alpha }}_4}+{\mathrm{J}}_{45}^{{\mathsf{\alpha }}_4}\right)$

where $M=\mathrm{diag}\left({\mathsf{\gamma }}_1^2,\frac{1}{{\mathsf{\gamma }}_1},\frac{1}{{\mathsf{\gamma }}_1}\right), \mathrm{and}$ ${\mathsf{\gamma }}_1,{\mathsf{\eta }}_1\ge 0$ and ${\mathsf{\alpha }}_4\ge 1$ are material parameters. Furthermore, to fulfill the stress-free initial configuration, the above material parameters should satisfy the requirement

(23)${\mathsf{\delta }}_2=2{\mathsf{\alpha }}_1+4{\mathsf{\alpha }}_2+2{\mathsf{\delta }}_1+2{\mathsf{\eta }}_1$

The above assumptions allow one to rewrite the constitutive relation to the form

(24)$S=2\left(1-\mathsf{\varphi }\right)\left[\left({\mathsf{\alpha }}_1+{\mathsf{\alpha }}_2{\mathrm{I}}_1\right)1-{\mathsf{\alpha }}_{2}C+\left({\mathsf{\delta }}_1{\mathrm{I}}_3-\frac{{\mathsf{\delta }}_2}{2}+\frac{{\mathsf{\eta }}_1}{{\left(\mathrm{trM}\right)}^{{\mathsf{\alpha }}_4}}{\mathrm{J}}_5^{{\mathsf{\alpha }}_4}\right)C^{-1}+\frac{{\mathsf{\eta }}_1}{{\left(\mathrm{trM}\right)}^{{\mathsf{\alpha }}_4}}{\mathrm{J}}_4^{{\mathsf{\alpha }}_4-1}M-\frac{{\mathsf{\eta }}_1}{{\left(\mathrm{trM}\right)}^{{\mathsf{\alpha }}_4}}{\mathrm{J}}_5^{{\mathsf{\alpha }}_4-1}{\mathrm{I}}_3{C}^{-1}M{C}^{-1}\right]$

As mentioned, the Caput-Almeida fractional derivative is used to describe the evolution of damage in time, which simultaneously introduces memory to the general constitutive law. One has then the following form for damage evolution

(25)${}_{\mathrm{t}-{\mathrm{l}}_{\mathrm{t}}}{}^{\mathrm{C}}\mathrm{D}_{\mathrm{t}}^{\mathsf{\alpha },\mathsf{\Psi }}\mathsf{\varphi }=\frac{1}{{\mathrm{T}}^{■}}\mathsf{\Phi }\frac{{\mathrm{I}}_{\mathsf{\varphi }}}{{\mathsf{\tau }}_{\mathsf{\varphi }}}-1$

where ${}_{\mathrm{t}-{\mathrm{l}}_{\mathrm{t}}}{}^{\mathrm{C}}\mathrm{D}_{\mathrm{t}}^{\mathsf{\alpha },\mathsf{\Psi }}$ denotes the Caputo-Almeida fractional derivative operator (Almeida, 2017)

(26)${}_{\mathrm{t}-{\mathrm{l}}_{\mathrm{t}}}{}^{\mathrm{C}}\mathrm{D}_{\mathrm{t}}^{\mathsf{\alpha },\mathsf{\Psi }}\mathrm{f}\left(\mathrm{t}\right)=\frac{1}{\mathsf{\Gamma }\left(\mathrm{n}-\mathsf{\alpha }\right)}{{\displaystyle \int }}_{\mathrm{t}-{\mathrm{l}}_{\mathrm{t}}}^{\mathrm{t}} {\mathsf{\Psi }}^{\prime }\left(\mathsf{\tau }\right){\left(\mathsf{\Psi }\left(\mathrm{t}\right)-\mathsf{\Psi }\left(\mathsf{\tau }\right)\right)}^{\mathrm{n}-\mathsf{\alpha }-1}{\left(\frac{1}{ {\mathsf{\Psi }}^{\prime }\left(\mathrm{t}\right)}\frac{\mathrm{d}}{\mathrm{dt}}\right)}^{\mathrm{n}}\mathrm{f}\left(\mathrm{t}\right)\mathrm{d}\mathsf{\tau }$

where $\mathsf{\alpha }$ denotes the order of fractional velocity of damage, ${\mathrm{l}}_{\mathrm{t}}$ is the time-length scale, $\mathsf{\Psi }\left(\mathrm{t}\right)\in {\mathrm{C}}^{\mathrm{n}}\in \left[\mathrm{t}-{\mathrm{l}}_{\mathrm{t}},\mathrm{t}\right]$ is an increasing function that fulfills the requirement ${\mathsf{\Psi }}^{\prime }\left(\mathrm{t}\right)\ne 0$, $\mathrm{n}=\left[\mathsf{\alpha }\right]+1$ and $\mathrm{n}=\mathsf{\alpha }$ for $\mathsf{\alpha }\in \mathrm{N}$, $\mathrm{T}$ is the characteristic time, $\left(■\right)$ defines power which depends on the selected $\mathsf{\Psi }\left(\mathrm{t}\right)$, $\mathsf{\Phi }$ is the overstress function, ${\mathrm{I}}_{\mathsf{\varphi }}$ is the stress intensity invariant, ${\mathsf{\tau }}_{\mathsf{\varphi }}$ is the static threshold stress for damage evolution, and the bracket $⟨\cdot ⟩$ defines the ramp function. It is vital to notice that the classical damage mechanical model (without memory effects) is a particular case of the proposed evolution law and is obtained for $\mathsf{\Psi }\left(\mathrm{t}\right)=\mathrm{t}$ and $\mathsf{\alpha }=1$.

Finally, based on results obtained by Voyiadjis and Sumelka [58]

(27)$\mathsf{\Phi }⟨\frac{{\mathrm{I}}_{\mathsf{\varphi }}}{{\mathsf{\tau }}_{\mathsf{\varphi }}}-1⟩={⟨\left(1-\mathsf{\varphi }\right){\left(\frac{\widetilde{{\mathrm{I}}_{\mathsf{\varphi }}}}{{\mathsf{\tau }}_{\mathsf{\varphi }}}\right)}^{{\mathrm{n}}_1}-1⟩}^{{\mathrm{n}}_2}$

where ${\mathrm{n}}_1,{\mathrm{n}}_2$ are material parameters, and

(28)$\widetilde{{\mathrm{I}}_{\mathsf{\varphi }}}={\mathrm{b}}^{\ast }\left({\mathrm{a}}^{\ast }+{\left(\frac{\overline{\mathsf{\Psi }}}{{\mathrm{c}}^{\ast }}\right)}^{{\mathrm{n}}^{\ast }+1}{{\displaystyle \sum }}_{{\mathrm{k}}^{\ast }=0}^{{\mathrm{n}}^{\ast }}\left(\begin{array}{c}{\mathrm{n}}^{\ast }+{\mathrm{k}}^{\ast }\\ {\mathrm{k}}^{\ast }\end{array}\right)\left(\begin{array}{c}2{\mathrm{n}}^{\ast }+1\\ {\mathrm{n}}^{\ast }-{\mathrm{k}}^{\ast }\end{array}\right){\left(-\frac{\overline{\mathsf{\Psi }}}{{\mathrm{c}}^{\ast }}\right)}^{{\mathrm{k}}^{\ast }}\right)+{\mathrm{f}}^{\ast }\exp \left({\mathrm{d}}^{\ast }\left(\overline{\mathsf{\Psi }}-{\mathrm{e}}^{\ast }\right)\right)$

In above equation ${\mathrm{a}}^{\ast },{\mathrm{b}}^{\ast },{\mathrm{c}}^{\ast },{\mathrm{d}}^{\ast },{\mathrm{e}}^{\ast },{\mathrm{f}}^{\ast },{\mathrm{n}}^{\ast }$ are material parameters, $\overline{\mathsf{\Psi }}=\overline{{\mathrm{I}}_{\mathsf{\varphi }}}\left(C,M\right)-\overline{{\mathrm{I}}_{\mathsf{\varphi }0}}$, $\overline{{\mathrm{I}}_{\mathsf{\varphi }0}}=$ $\overline{{\mathrm{I}}_{\mathsf{\varphi }}}\left(C,M\right){|}_{\mathrm{C}=0}$, and

(29)$\overline{{\mathrm{I}}_{\mathsf{\varphi }}}={\mathsf{\alpha }}_1^{\ast }{\mathrm{I}}_1+{\mathsf{\alpha }}_2^{\ast }{\mathrm{I}}_2+{\mathsf{\delta }}_1^{\ast }{\mathrm{I}}_3-{\mathsf{\delta }}_2^{\ast }\ln \left(\sqrt{{\mathrm{I}}_3}\right)+\frac{1}{{\mathsf{\alpha }}_4^{\ast }{\left(\mathrm{trM}\right)}^{{\mathsf{\alpha }}_4^{\ast }}}{\mathsf{\eta }}_1^{\ast }\left({\mathrm{J}}_4^{{\mathsf{\alpha }}_4^{\ast }}+{\mathrm{J}}_5^{{\mathsf{\alpha }}_4^{\ast }}\right)$

The material parameters are such that ${\mathsf{\alpha }}_1^{\ast },{\mathsf{\alpha }}_2^{\ast },{\mathsf{\delta }}_1^{\ast },{\mathsf{\delta }}_2^{\ast },{\mathsf{\eta }}_1^{\ast }\ge 0,{\mathsf{\alpha }}_4^{\ast }\ge 1$, and ${\mathsf{\delta }}_2^{\ast }=2{\mathsf{\alpha }}_1^{\ast }+4{\mathsf{\alpha }}_2^{\ast }+2{\mathsf{\delta }}_1^{\ast }+2{\mathsf{\eta }}_1^{\ast }$. Altogether we have 23 material parameters:

• ➢. 4 defining anisotropic hyperelastic part ${\mathsf{\alpha }}_1,{\mathsf{\alpha }}_2,{\mathsf{\delta }}_1,{\mathsf{\alpha }}_4,{\mathsf{\eta }}_1$

• ➢. 19 defining damage evolution $\mathsf{\alpha },{\mathrm{l}}_{\mathrm{t}},\mathrm{T},{\mathrm{n}}_1,{\mathrm{n}}_2,{\mathsf{\tau }}_{\mathsf{\varphi }},{\mathsf{\alpha }}_1^{\ast },{\mathsf{\alpha }}_2^{\ast },{\mathsf{\delta }}_1^{\ast },{\mathsf{\alpha }}_4^{\ast },{\mathsf{\eta }}_1^{\ast },{\mathrm{a}}^{\ast },{\mathrm{b}}^{\ast },{\mathrm{c}}^{\ast },{\mathrm{d}}^{\ast },{\mathrm{e}}^{\ast },{ \mathrm{f}}^{\ast },{\mathrm{n}}^{\ast }$ plus the modulation function $\mathsf{\Psi }$.

For model validation, the number of material parameters is initially reduced to 17 under the assumption that ${\mathsf{\alpha }}_1={\mathsf{\alpha }}_1^{\ast },{\mathsf{\alpha }}_2={\mathsf{\alpha }}_2^{\ast },{\mathsf{\delta }}_1={\mathsf{\delta }}_1^{\ast },{\mathsf{\alpha }}_4={\mathsf{\alpha }}_4^{\ast },{\mathsf{\eta }}_1={\mathsf{\eta }}_1^{\ast }$, and that full-memory assumption holds, therefore ${\mathrm{l}}_{\mathrm{t}}$ could be reduced as well.

3.5. Model Validation

The hyperelastic time-fractional damage (HTFD) model has been implemented in the Mathematica software and identified based on experimental evidence on the LIB specimen. Calibration has been prepared under the assumption of best fitting curve utilizing an optimization task (FindFit function from the Mathematica software [60], utilizing “Quasi Newton” optimization, i.e., the quasi-Newton BFGS approximation to the Hessian was applied) for selected modulation functions $\mathsf{\Psi }$. The best results were obtained for $\mathsf{\Psi }\left(\mathrm{t}\right)={\mathrm{t}}^3+1$, therefore, in consequence in damage evolution law $■=3\mathsf{\alpha }$. Material parameters are collected in Table 1, whereas the comparison of experimental data with HTFD model predictions is presented in Figure 7.

Figure 7 presents the stretch (${\mathsf{\lambda }}_1$) vs. stress (${\mathrm{P}}_{11}-$1st Piola-Kirchhoff stress component) for 3 cases: experimental data, undamaged answer (purely elastic case; ${\mathsf{\tau }}_{\mathsf{\varphi }}\to \infty$), and full-range HTFD model prediction. Since the main novelty of the proposed description is the application of fractional calculus for the definition of internal state variable evolution, one should notice that the fractional-velocity of damage of order $\mathsf{\alpha }=0.72$ (with the physical dimension ${\mathrm{s}}^{-3\ast 0.72}$) provided the best approximation of experimental data. One can observe that the first stage of softening appears for ${\mathsf{\lambda }}_1\approx 1.08$ and rapidly grows up to ${\mathsf{\lambda }}_1\approx 1.18$. In the next step, damage grows slowly and at the limit ${\mathsf{\lambda }}_1\approx 1.48$ jumps drastically up to its extreme value $\mathsf{\varphi }\to 1^{-}$.

Finally, it is clear that the HTFD model provides proper predictions in the phenomenological sense. Complex damage phenomena in real specimens can be correctly described in a homogenized sense through a single scalar variable whose evolution follows a time-factional velocity curve. One should point out that such flexibility is obtained through one parameter ($\mathsf{\alpha }$) and a single modulation function ($\mathsf{\Psi }$), compared to the classical damage mechanics approach.

4. Conclusions

The indentation test was performed using the cylindrical indenter. This research addresses the difficulty of converting force depth data into engineering stress-strain relations. Measuring the projected contact area between two cylinders is one of the issues posed by this mechanical test, which this study addresses by proposing a few analytical approaches and the FEM model. However, another innovation of this research is a multistep indentation test that was performed to measure the reduction of Young’s modulus. In contrast to previously published models, which can predict elastoplastic behavior but cannot predict LIB cell damage, an analytical damage model has been developed in this work. Hence, a definition is provided for an anisotropic hyperelastic theory that accounts for fractional scalar damage development regulated by the Caputo-Almeida fractional derivative. The refined damage model is then implemented for the constrained indentation configuration of the experiment. This study demonstrates that the numerical findings qualitatively and quantitatively closely match the experimental data.

Author Contributions

G.Z.V. and W.S. Towards: Review and editing; E.A.: Writing—reviewing and editing, visualization, simulation, experiment, Analytical; B.Ł.: Analytical. All authors have read and agreed to the published version of the manuscript.

Informed Consent Statement

Informed consent was obtained from all subjects involved in the study.

Data Availability Statement

The data presented in this study are available on request from the corresponding author. The data are not publicly available due to long-time experiment and rich information.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Table

Figure 1. Experimental setup of indentation of a cylindrical LIB cell with a cylindrical indenter.

Figure 2. Five single-loading indentation tests at five different penetration depths (a), unloading segments (b), and degradation of Young’s modulus (c).

Figure 3. (a) Deformed LIB cell through the indentation test, (b) Schematic of a cylindrical LIB cell being indented by a rigid cylindrical indenter, (c) Graphical representation of fitting an ellipsoidal cap to the deformed surface.

Figure 4. Comparison between the calculated contact area from the analytical and FE model.

Figure 5. Measured indentation load–depth data.

Figure 6. The effective stress-strain curve.

Figure 7. Comparison of experimental data (orange line) with HTFD model predictions (blue line) is contrasted with purely elastic (undamaged; [Forumla omitted. See PDF.]) case (green line).

Table 1

HTFD material parameters for LIB.

${\mathsf{\alpha }}_1={\mathsf{\alpha }}_1^{\ast }=10.50 \left[\mathrm{MPa}\right]$	$\mathrm{T}=0.257 \left[\mathrm{s}\right]$	${\mathrm{c}}^{\ast }=1.72\left[\mathrm{MPa}\right]$
${\mathsf{\alpha }}_2={\mathsf{\alpha }}_2^{\ast }=7.60 \left[\mathrm{MPa}\right]$	${\mathrm{n}}_1=1.76\left[-\right]$	${\mathrm{d}}^{\ast }=9.78\left[{\mathrm{MPa}}^{-1}\right]$
${\mathsf{\delta }}_1={\mathsf{\delta }}_1^{\ast }=7.00 \left[\mathrm{MPa}\right]$	${\mathrm{n}}_2=1.83 \left[-\right]$	${\mathrm{e}}^{\ast }=15.29\left[\mathrm{MPa}\right]$
${\mathsf{\alpha }}_4={\mathsf{\alpha }}_4^{\ast }=18.9 \left[-\right]$	${\mathsf{\tau }}_{\mathsf{\varphi }}=0.55 \left[\mathrm{MPa}\right]$	${\mathrm{f}}^{\ast }=1.01\left[\mathrm{MPa}\right]$
${\mathsf{\eta }}_1={\mathsf{\eta }}_1^{\ast }=1.90\left[\mathrm{MPa}\right]$	${\mathrm{a}}^{\ast }=0.55\left[-\right]$	${\mathrm{n}}^{\ast }=2\left[-\right]$
$\mathsf{\alpha }=0.72 \left[-\right]$	${\mathrm{b}}^{\ast }=1.01\left[\mathrm{MPa}\right]$

References

1. Masias, A.; Felten, N.; Sakamoto, J. Characterizing the mechanical behavior of lithium in compression. J. Mater. Res.; 2021; 36, pp. 729-739. [DOI: https://dx.doi.org/10.1557/s43578-020-00028-x]

2. Kermani, G.; Keshavarzi, M.M.; Sahraei, E. Deformation of lithium-ion batteries under axial loading: Analytical model and Representative Volume Element. Energy Rep.; 2021; 7, pp. 2849-2861. [DOI: https://dx.doi.org/10.1016/j.egyr.2021.05.015]

3. Spielbauer, M.; Berg, P.; Soellner, J.; Peters, J.; Schaeufl, F.; Rosenmüller, C.; Bohlen, O.; Jossen, A. Experimental investigation of the failure mechanism of 18650 lithium-ion batteries due to shock and drop. J. Energy Storage; 2021; 43, 103213. [DOI: https://dx.doi.org/10.1016/j.est.2021.103213]

4. Zhu, F.; Zhou, R.; Sypeck, D.J. Numerical modeling and safety design for lithium-ion vehicle battery modules subject to crush loading. Energies; 2020; 14, 118. [DOI: https://dx.doi.org/10.3390/en14010118]

5. Kisters, T.; Sahraei, E.; Wierzbicki, T. Dynamic impact tests on lithium-ion cells. Int. J. Impact Eng.; 2017; 108, pp. 205-216. [DOI: https://dx.doi.org/10.1016/j.ijimpeng.2017.04.025]

6. Luo, H.; Xia, Y.; Zhou, Q. Mechanical damage in a lithium-ion pouch cell under indentation loads. J. Power Source; 2017; 357, pp. 61-70. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2017.04.101]

7. Keshavarzi, M.M.; Gilaki, M.; Sahraei, E. Characterization of in-situ material properties of pouch lithium-ion batteries in tension from three-point bending tests. Int. J. Mech. Sci.; 2022; 219, 107090. [DOI: https://dx.doi.org/10.1016/j.ijmecsci.2022.107090]

8. Zhu, X.; Wang, H.; Wang, X.; Gao, Y.; Allu, S.; Cakmak, E.; Wang, Z. Internal short circuit and failure mechanisms of lithium-ion pouch cells under mechanical indentation abuse conditions: An experimental study. J. Power Source; 2020; 455, 227939. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2020.227939]

9. Lian, J.; Koch, M.; Li, W.; Wierzbicki, T.; Zhu, J. Mechanical deformation of Lithium-ion pouch cells under in-plane loads—Part II: Computational modeling. J. Electrochem. Soc.; 2020; 167, 090556. [DOI: https://dx.doi.org/10.1149/1945-7111/ab9eee]

10. Gao, Z.; Zhang, X.; Xiao, Y.; Gao, H.; Wang, H.; Piao, C. Influence of low-temperature charge on the mechanical integrity behavior of 18650 lithium-ion battery cells subject to lateral compression. Energies; 2019; 12, 797. [DOI: https://dx.doi.org/10.3390/en12050797]

11. Ren, F.; Cox, T.; Wang, H. Thermal runaway risk evaluation of Li-ion cells using a pinch–torsion test. J. Power Source; 2014; 249, pp. 156-162. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2013.10.058]

12. Ellersdorfer, C.; Höschele, P.; Heider, E.; Kovachev, G.; Gstrein, G. Safety Assessment of High Dynamic Pre-Loaded Lithium Ion Pouch Cells. Batteries; 2023; 9, 71. [DOI: https://dx.doi.org/10.3390/batteries9020071]

13. Gilaki, M.; Song, Y.; Sahraei, E. Homogenized characterization of cylindrical Li-ion battery cells using elliptical approximation. Int. J. Energy Res.; 2022; 46, pp. 5908-5923. [DOI: https://dx.doi.org/10.1002/er.7531]

14. Kermani, G.; Dixon, B.; Sahraei, E. Elliptical lithium-ion batteries: Transverse and axial loadings under wet/dry conditions. Energy Sci. Eng.; 2019; 7, pp. 890-898. [DOI: https://dx.doi.org/10.1002/ese3.318]

15. Aphale, S.; Murugkar, M.; Lulla, S.; Date, P. Deformation and Failure Behavior of Cylindrical Lithium-Ion Batteries Subjected to Mechanical Loading; SAE Technical Paper: Warrendale, PA, USA, 2020; pp. 0148-7191.

16. Avdeev, I.; Gilaki, M. Structural analysis and experimental characterization of cylindrical lithium-ion battery cells subject to lateral impact. J. Power Source; 2014; 271, pp. 382-391. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2014.08.014]

17. Raffler, M.; Sevarin, A.; Ellersdorfer, C.; Heindl, S.F.; Breitfuss, C.; Sinz, W. Finite element model approach of a cylindrical lithium ion battery cell with a focus on minimization of the computational effort and short circuit prediction. J. Power Source; 2017; 360, pp. 605-617. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2017.06.028]

18. Sheikh, M.; Elmarakbi, A.; Elkady, M. Thermal runaway detection of cylindrical 18650 lithium-ion battery under quasi-static loading conditions. J. Power Source; 2017; 370, pp. 61-70. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2017.10.013]

19. Chung, S.H.; Tancogne-Dejean, T.; Zhu, J.; Luo, H.; Wierzbicki, T. Failure in lithium-ion batteries under transverse indentation loading. J. Power Source; 2018; 389, pp. 148-159. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2018.04.003]

20. Muresanu, A.D.; Dudescu, M.C. Numerical and Experimental Evaluation of a Battery Cell under Impact Load. Batteries; 2022; 8, 48. [DOI: https://dx.doi.org/10.3390/batteries8050048]

21. Zhu, F.; Zhou, R.; Sypeck, D.; Deng, J.; Bae, C. Failure behavior of prismatic Li-ion battery cells under abuse loading condition-A combined experimental and computational study. J. Energy Storage; 2022; 48, 103969. [DOI: https://dx.doi.org/10.1016/j.est.2022.103969]

22. Li, W.; Xing, B.; Watkins, T.R.; Xia, Y.; Wang, H.; Zhu, J. Damage of prismatic lithium-ion cells subject to bending: Test, model, and detection. EcoMat; 2022; 4, e12257. [DOI: https://dx.doi.org/10.1002/eom2.12257]

23. Wang, H.; Simunovic, S.; Maleki, H.; Howard, J.N.; Hallmark, J.A. Internal configuration of prismatic lithium-ion cells at the onset of mechanically induced short circuit. J. Power Source; 2016; 306, pp. 424-430. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2015.12.026]

24. Lamb, J.; Orendorff, C.J. Evaluation of mechanical abuse techniques in lithium ion batteries. J. Power Source; 2014; 247, pp. 189-196. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2013.08.066]

25. Mao, B.; Chen, H.; Cui, Z.; Wu, T.; Wang, Q. Failure mechanism of the lithium ion battery during nail penetration. Int. J. Heat Mass Transf.; 2018; 122, pp. 1103-1115. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2018.02.036]

26. Feng, X.; Sun, J.; Ouyang, M.; Wang, F.; He, X.; Lu, L.; Peng, H. Characterization of penetration induced thermal runaway propagation process within a large format lithium ion battery module. J. Power Source; 2015; 275, pp. 261-273. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2014.11.017]

27. Ma, T.; Chen, L.; Liu, S.; Zhang, Z.; Xiao, S.; Fan, B.; Liu, L.; Lin, C.; Ren, S.; Wang, F. Mechanics-morphologic coupling studies of commercialized lithium-ion batteries under nail penetration test. J. Power Source; 2019; 437, 226928. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2019.226928]

28. Aiello, L.; Gstrein, G.; Erker, S.; Kaltenegger, B.; Ellersdorfer, C.; Sinz, W. Optimized Nail for Penetration Test on Lithium-Ion Cells and Its Utilization for the Validation of a Multilayer Electro-Thermal Model. Batteries; 2022; 8, 32. [DOI: https://dx.doi.org/10.3390/batteries8040032]

29. Doose, S.; Haselrieder, W.; Kwade, A. Effects of the Nail Geometry and Humidity on the Nail Penetration of High-Energy Density Lithium Ion Batteries. Batteries; 2021; 7, 6. [DOI: https://dx.doi.org/10.3390/batteries7010006]

30. Sahraei, E.; Meier, J.; Wierzbicki, T. Characterizing and modeling mechanical properties and onset of short circuit for three types of lithium-ion pouch cells. J. Power Source; 2014; 247, pp. 503-516. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2013.08.056]

31. Sahraei, E.; Campbell, J.; Wierzbicki, T. Modeling and short circuit detection of 18650 Li-ion cells under mechanical abuse conditions. J. Power Source; 2012; 220, pp. 360-372. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2012.07.057]

32. Sahraei, E.; Kahn, M.; Meier, J.; Wierzbicki, T. Modelling of cracks developed in lithium-ion cells under mechanical loading. RSC Adv.; 2015; 5, pp. 80369-80380. [DOI: https://dx.doi.org/10.1039/C5RA17865G]

33. Wang, W.; Yang, S.; Lin, C. Clay-like mechanical properties for the jellyroll of cylindrical Lithium-ion cells. Appl. Energy; 2017; 196, pp. 249-258. [DOI: https://dx.doi.org/10.1016/j.apenergy.2017.01.062]

34. Cai, W.; Wang, H.; Maleki, H.; Howard, J.; Lara-Curzio, E. Experimental simulation of internal short circuit in Li-ion and Li-ion-polymer cells. J. Power Source; 2011; 196, pp. 7779-7783. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2011.04.024]

35. Sahraei, E.; Hill, R.; Wierzbicki, T. Calibration and finite element simulation of pouch lithium-ion batteries for mechanical integrity. J. Power Source; 2012; 201, pp. 307-321. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2011.10.094]

36. Greve, L.; Fehrenbach, C. Mechanical testing and macro-mechanical finite element simulation of the deformation, fracture, and short circuit initiation of cylindrical Lithium ion battery cells. J. Power Source; 2012; 214, pp. 377-385. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2012.04.055]

37. Wang, Q.; Mao, B.; Stoliarov, S.I.; Sun, J. A review of lithium ion battery failure mechanisms and fire prevention strategies. Prog. Energy Combust. Sci.; 2019; 73, pp. 95-131. [DOI: https://dx.doi.org/10.1016/j.pecs.2019.03.002]

38. Xu, J.; Liu, B.; Hu, D. State of Charge Dependent Mechanical Integrity Behavior of 18650 Lithium-ion Batteries. Sci. Rep.; 2016; 6, 21829. [DOI: https://dx.doi.org/10.1038/srep21829] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/26911922]

39. Zhu, J.; Zhang, X.; Sahraei, E.; Wierzbicki, T. Deformation and failure mechanisms of 18650 battery cells under axial compression. J. Power Source; 2016; 336, pp. 332-340. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2016.10.064]

40. Song, Y.; Gilaki, M.; Keshavarzi, M.M.; Sahraei, E. A universal anisotropic model for a lithium-ion cylindrical cell validated under axial, lateral, and bending loads. Energy Sci. Eng.; 2022; 10, pp. 1431-1448. [DOI: https://dx.doi.org/10.1002/ese3.1111]

41. Wierzbicki, T.; Sahraei, E. Homogenized mechanical properties for the jellyroll of cylindrical Lithium-ion cells. J. Power Source; 2013; 241, pp. 467-476. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2013.04.135]

42. Wang, L.; Yin, S.; Xu, J. A detailed computational model for cylindrical lithium-ion batteries under mechanical loading: From cell deformation to short-circuit onset. J. Power Source; 2019; 413, pp. 284-292. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2018.12.059]

43. Zhu, J.; Li, W.; Wierzbicki, T.; Xia, Y.; Harding, J. Deformation and failure of lithium-ion batteries treated as a discrete layered structure. Int. J. Plast.; 2019; 121, pp. 293-311. [DOI: https://dx.doi.org/10.1016/j.ijplas.2019.06.011]

44. Li, W.; Zhu, J. A Large Deformation and Fracture Model of Lithium-Ion Battery Cells Treated as a Homogenized Medium. J. Electrochem. Soc.; 2020; 167, 120504. [DOI: https://dx.doi.org/10.1149/1945-7111/aba936]

45. Voyiadjis, G.Z.; Akbari, E.; Kattan, P.I. Damage model for lithium-ion batteries with experiments and simulations. J. Energy Storage; 2023; 57, 106285. [DOI: https://dx.doi.org/10.1016/j.est.2022.106285]

46. Ali, M.Y.; Lai, W.-J.; Pan, J. Computational models for simulation of a lithium-ion battery module specimen under punch indentation. J. Power Source; 2015; 273, pp. 448-459. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2014.09.072]

47. Lai, W.-J.; Ali, M.Y.; Pan, J. Mechanical behavior of representative volume elements of lithium-ion battery modules under various loading conditions. J. Power Source; 2014; 248, pp. 789-808. [DOI: https://dx.doi.org/10.1016/j.jpowsour.2013.09.128]

48. Murakami, S. Continuum Damage Mechanics: A Continuum Mechanics Approach to the Analysis of Damage and Fracture; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2012; Volume 185.

49. Chow, C.; Wang, J. An anisotropic theory of continuum damage mechanics for ductile fracture. Eng. Fract. Mech.; 1987; 27, pp. 547-558. [DOI: https://dx.doi.org/10.1016/0013-7944(87)90108-1]

50. Perzyna, P. The constitutive equations for rate sensitive plastic materials. Q. Appl. Math.; 1963; 20, pp. 321-332. [DOI: https://dx.doi.org/10.1090/qam/144536]

51. Almeida, R. A Caputo fractional derivative of a function with respect to another function. Commun. Nonlinear Sci. Numer. Simul.; 2017; 44, pp. 460-481. [DOI: https://dx.doi.org/10.1016/j.cnsns.2016.09.006]

52. Hertz, H. Miscellaneous Papers; Macmillan: Stuttgart, Germany, 1896.

53. Wang, H.; Duan, X.; Liu, B. The Anisotropic Homogenized Model for Pouch Type Lithium-Ion Battery Under Various Abuse Loadings. J. Electrochem. Energy Convers. Storage; 2021; 18, 021015. [DOI: https://dx.doi.org/10.1115/1.4049239]

54. Popov, V.L. Contact Mechanics and Friction; Springer: Berlin/Heidelberg, Germany, 2010.

55. Yang, S.; Wang, W.; Lin, C.; Shen, W.; Li, Y. Improved constitutive model of the jellyroll for cylindrical lithium ion batteries considering microscopic damage. Energy; 2019; 185, pp. 202-212. [DOI: https://dx.doi.org/10.1016/j.energy.2019.07.012]

56. Xu, B.; Chen, X. Determining engineering stress–strain curve directly from the load–depth curve of spherical indentation test. J. Mater. Res.; 2010; 25, pp. 2297-2307. [DOI: https://dx.doi.org/10.1557/jmr.2010.0310]

57. Kalidindi, S.R.; Pathak, S. Determination of the effective zero-point and the extraction of spherical nanoindentation stress–strain curves. Acta Mater.; 2008; 56, pp. 3523-3532. [DOI: https://dx.doi.org/10.1016/j.actamat.2008.03.036]

58. Voyiadjis, G.Z.; Sumelka, W. Brain modelling in the framework of anisotropic hyperelasticity with time fractional damage evolution governed by the Caputo-Almeida fractional derivative. J. Mech. Behav. Biomed. Mater.; 2019; 89, pp. 209-216. [DOI: https://dx.doi.org/10.1016/j.jmbbm.2018.09.029] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30292967]

59. Sumelka, W.; Łuczak, B.; Gajewski, T.; Voyiadjis, G.Z. Modelling of AAA in the framework of time-fractional damage hyperelasticity. Int. J. Solids Struct.; 2020; 206, pp. 30-42. [DOI: https://dx.doi.org/10.1016/j.ijsolstr.2020.08.015]

60. FindFit. Available online: https://reference.wolfram.com/language/ref/FindFit.html (accessed on 20 March 2023).