# Comprehensive modeling of lithium-ion batteries thermal runaway behavior

Sirar Chakaroun, Pierre Coste, Stephanie de Persis, Alain Bengaouer, Jérôme Cognard, Nabiha Chaumeix

► To cite this version:

Sirar Chakaroun, Pierre Coste, Stephanie de Persis, Alain Bengaouer, Jérôme Cognard, et al. Comprehensive modeling of lithium-ion batteries thermal runaway behavior. ISFEH'11 - 11th International Seminar on Fire and Explosion Hazards, Jun 2025, Rome, Italy. cea-05527865

HAL Id: cea-05527865

https://cea.hal.science/cea-05527865v1

Submitted on 25 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Comprehensive Modeling of Lithium-Ion Batteries Thermal Runaway Behavior

Chakaroun S.<sup>1,2,4,*</sup>, Coste P.<sup>1</sup>, De Persis S.<sup>2,4</sup>, Bengaouer A.<sup>3</sup>, Cognard J.<sup>3</sup>, Chaumeix N.<sup>2</sup>

<sup>1</sup> Univ. Grenoble Alpes, CEA, Liten, DTCH, 38000 Grenoble, France

<sup>2</sup> CNRS, ICARE, 45100, Orleans, France

<sup>3</sup> Univ. Grenoble Alpes, CEA, Liten, DEHT, 38000 Grenoble, France

<sup>4</sup> Université d'Orléans, 45067 Orléans Cedex 2, France

*Corresponding author email: sirar.chakaroun@cea.fr

## ABSTRACT

The thermal runaway behavior of lithium-ion batteries is a critical area of research due to the complex thermal propagation phenomena involved. These phenomena include internal heat generation from exothermic reactions leading to the ejection of substantial amounts of combustible gases and condensed particles (solid and liquid) during venting. This study develops a transient 3D model using Reynolds-Averaged Navier-Stokes (RANS) approach considering the coupling between internal heat generation and the venting process. The model incorporates internal heat generation, combustion heat, and heat transfer with the surroundings through vented material. Solid particles are modeled via a Lagrangian interaction with the gas phase. The model was validated against 18650 lithium-ion batteries with Nickel Manganese Cobalt (NMC) cathodes in two configurations: a closed vessel and an open-end vessel. For the closed vessel, the model accurately predicted energy releases between 23 kJ and 77 kJ across states of charge from 15% to 100%, with a relative error of less than 9%. For the open-end vessel, the model captured temperature distributions at various distances from the cell, showing good agreement with experimental data. These results illustrate the model's capability to predict key parameters of thermal runaway events under varying containment conditions.

**KEYWORDS:** CFD, Thermal runaway, Lithium-ion batteries, Safety, Venting process.

## NOMENCLATURE

$A_x$ pre-exponential factor

$c_p$ specific heat [J/kg K]

$E_{a,x}$ activation energy of the reaction [J/mol]

$H_{ec}$ heat of reaction [J/m³] added to the thermal abuse heat

$H_x$ heat release [J/kg]

$h$ specific enthalpy [J/kg]

$J_i$ diffusion flux of species i [kg m⁻² s⁻¹]

$k$ thermal conductivity [W·m⁻¹·K⁻¹]

$M_i$ molecular weight of species i [g/mol]

$\dot{q}$ heat source term [J·m⁻³·s⁻¹]

$R_u$ universal gas constant, $R_u = 8.314[J/mol K]$

$R_x$ reaction rate due to chemical reactions [s⁻¹]

$S_x$ heat source term [J·m⁻³·s⁻¹]

$S_{oc}$ progress variable for internal short reaction

$t$ time [s]

$T$ temperature [K]

$W_x$ specific content of x [g·m⁻³]

$X_i$ mole fraction of species i

$Y_i$ mass fraction of species i

### Greek

$\bar{\tau}$ stress tensor [Pa]

$\vec{v}$ velocity vector

$\rho$ density [kg·m⁻³]

$\nabla$ partial derivative of a quantity with respect to all directions in the chosen coordinate system

### Subscripts

sei SEI decomposition reaction

ne negative electrode-electrolyte reaction

pe positive electrode-electrolyte reaction

e electrolyte decomposition reaction

i initial value

f final value

x sei, ne, pe, or e
# INTRODUCTION

Lithium-ion batteries (LIBs) have become essential across numerous sectors, powering devices from consumer electronics to electric vehicles, energy storage systems, and even critical transport infrastructures. However, the inherent risk of thermal runaway (TR) - a phenomenon in which exothermic reactions within the battery cell generate sufficient heat to cause self-sustained chain reactions leading to fire and potential explosion - has resulted in numerous high-impact incidents. Notable cases, including the Samsung Galaxy Note 7 recall, Tesla vehicle fires, and the fatal crashes of UPS Flight 6 and Asiana Airlines Flight 991, underscore the dangers associated with these batteries. Incidents at large-scale energy storage facilities and in maritime cargo have further highlighted the challenges in managing LIB safety. Consequently, regulatory bodies have implemented stringent measures, including restrictions on cell charge levels during transport, specialized packaging requirements, and advanced fire suppression protocols. Understanding the mechanisms behind TR and the implementation of mitigation measures are therefore essential for advancing LIB safety and ensuring their reliable integration into modern technology and energy infrastructures [1-3]. TR can lead to the release of significant amounts of energy, flammable gases, and particulate matter, posing risks such as fires and explosions. Understanding and accurately modeling these complex processes is vital for enhancing battery safety. Previous research has largely focused on specific aspects of TR: (i) chemical kinetics involved in exothermic reactions, (ii) thermal propagation behavior within the cell. One commonly used approach is the Arrhenius model [4-6], which relates temperature, chemical reaction rates, and the concentration of reactants to describe the energy released during TR events. While these models capture some essential aspects of the TR process, they often lack the ability to comprehensively simulate the coupling between internal heat generation and the external venting of gases and particles, which can significantly affect the total energy release and the potential hazard level. Recent work has aimed to address these limitations by coupling internal cell reactions with external phenomena like gas venting. The need for a model that considers internal heat generation within the cell and the energy released during venting was highlighted in [7]. A model [8], focusing on the interaction between heat generation and the jet flow during TR, led to promising results in terms of capturing the venting gases during thermal runaway behavior. However, it did not address the injection of solid particles, a critical component of TR events, nor did it explore the influence of key parameters such as state of charge (SOC) and total energy release from the battery cell. Furthermore, while some multi-scale models have been developed to capture the dynamic behavior of gases and particles during TR, their validation has often been limited to specific cell types, reducing their generalizability [9]. Consequently, there is a strong need for a more comprehensive approach that can simulate the combined effects of internal heat generation, gas venting, and solid particle release in various LIB systems. This study aims to fill these gaps by developing a transient 3D computational fluid dynamics (CFD) model using the Reynolds-Averaged Navier-Stokes (RANS) approach. In this work, simulations were performed using the commercial software ANSYS Fluent® 2023R1 [10]. The model incorporates key TR phenomena, including internal heat generation, combustion of vented gases, and the interaction between vented gases and solid particles, modeled using a Lagrangian approach. A critical feature of the model is the coupling between internal and external processes, ensuring mass conservation and accurately capturing the temperature, gas venting, and particle release behavior. The model is compared to experimental data obtained for 18650 lithium-ion cells with Nickel Manganese Cobalt (NMC) cathodes, across a range of SOCs from 15 to 100%, corresponding to energy releases between 23 and 77 kJ.

# MATERIAL AND METHODS

## Selection of an experiment

A series of tests are conducted at CEA-Liten/DEHT on various commercially available LIBs. The primary objective of these tests is to measure the energy released by the battery cells during TR at
different SOC levels. A commercial LG HG2 18650 Li-ion cell is selected. A resistance heater is attached to the external shell of the cell, gradually increasing its temperature from ambient to the critical temperature required to initiate TR. The resistance heater generates a constant heat rate of 6°C/min, ensuring a controlled and steady increase in temperature until the onset of TR. Two types of experiments, closed and open-end vessels, are performed to analyze the cells behavior under thermal abuse conditions. For the closed configuration (Fig. 1), tests are conducted in a closed vessel with a 0.54 L volume, with variations in the cell SOC level. The key measurements recorded include the initial and final mass of the cell, pressure inside the vessel, and temperature at various locations, such as within the vessel, near the outer surface of the cell, and along the external shell of the vessel. These measurements are used to validate the model's predictions of pressure buildup, heat generation, and temperature distribution during TR. To minimize heat loss to the surroundings, the vessel is placed inside a closed box made of vermiculite, providing thermal insulation. For the open-end configuration (Fig. 2), the vessel is modified with a 2 cm diameter opening near the end to allow for venting and to study gas ejection behavior during TR. The key measurements for this configuration include the temperatures recorded by three thermocouples (0.5 mm type-K): TC0, TC1, and TC2, placed 1.5 cm, 20 cm, and 40 cm from the opening, respectively. These thermocouples track the temperature changes at different distances from the cell to evaluate the thermal impact of vented gases. The vessel is placed inside a closed box of 1 m³, which captures the vented gases and enables the analysis of their effects. The temperature profiles and pressure changes observed are compared with the simulation results to validate the model's ability to replicate the experimental conditions.

<figure><img src="image_2.png" /><figcaption>Fig. 1. View of the LIB cell inside the stainless steel vessel.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Fig. 2. The open-end vessel.</figcaption></figure>

## Model description

The numerical model developed in this study aims at simulating the TR behavior of lithium-ion batteries by capturing both the internal heat generation mechanisms and the external venting of gases and particles. This section details the three primary components of the model: (1) the modeling of heat generation inside the battery cell, (2) the coupling between the ejected gases and particles during the venting process, (3) and the reaction of the gas mixture outside the cell.

## Modelling of heat generation inside the cell

The energy transport within a solid can be expressed as well as the heat source term $\dot{q}$:

$$ \rho c_p \frac{\partial T}{\partial t} = \nabla \cdot (k \nabla T) + \dot{q} \quad (1) $$

$$ \dot{q} = S_{abuse\_chem} + S_{ec} \quad (2) $$

Where: k is the thermal conductivity of the cell is considered orthotropic, $K_z=30.4 \text{ W.m}^{-1}.\text{K}^{-1}$, $K_i=0.2 \text{ W.m}^{-1}.\text{K}^{-1}$ for the axial and radial directions, respectively [8]; $S_{abuse\_chem}$ is a model consisting of 4 key reactions contributing to the thermal energy released inside the cell during TR. They include the exothermic decomposition of: the Solid Electrolyte Interface (SEI), the negative electrode-electrolyte reaction (ne), the positive electrode-electrolyte reaction (pe), and the electrolyte decomposition reaction (ele). This can be written as [11]:
$$S_{abuse\_chem} = S_{sei} + S_{ne} + S_{pe} + S_{ele} \quad (3)$$

SEI decomposition reaction, $R_{sei}$, negative-solvent reaction, $R_{ne}$, positive-solvent reaction, $R_{pe}$ and electrolyte decomposition reaction, $R_e$, are defined as shown below [11-12]:

$$R_{sei}(T, c_{sei}) = A_{sei} \times \exp\left(-\frac{E_{a,sei}}{R_u T}\right) \times c_{sei}^{m_{sei}} \quad (4)$$

$$S_{sei} = H_{sei} \times W_c \times R_{sei} \quad (5)$$

$$\frac{dc_{sei}}{dt} = -R_{sei} \quad (6)$$

$$R_{ne}(T, c_{ne}, t_{sei}) = A_{ne} \times \exp\left(-\frac{t_{sei}}{t_{sei,ref}}\right) \times c_{ne}^{m_{ne,n}} \times \exp\left(-\frac{E_{a,ne}}{R_u T}\right) \quad (7)$$

$$S_{ne} = H_{ne} \times W_c \times R_{ne} \quad (8)$$

$$\frac{dc_{neg}}{dt} = -R_{ne} \quad (9)$$

$$\frac{dt_{sei}}{dt} = R_{ne} \quad (10)$$

$$R_{pe}(T, \alpha) = A_{pe} \times \alpha^{m_{pe,p1}} \times (1 - \alpha)^{m_{pe,p2}} \times \exp\left(-\frac{E_{a,pe}}{R_u T}\right) \quad (11)$$

$$S_{pe} = H_{pe} \times W_c \times R_{pe} \quad (12)$$

$$\frac{d\alpha}{dt} = R_{pe} \quad (13)$$

$$R_e(T, c_e) = A_e \times \exp\left(-\frac{E_{a,e}}{R_u T}\right) \times c_e^{m_e} \quad (14)$$

$$S_{ele} = H_e \times W_e \times R_e \quad (15)$$

$$\frac{dc_e}{dt} = -R_e \quad (16)$$

$c_{sei}$, $c_{ne}$, $\alpha$, and $c_e$ are dimensionless variables representing the fraction of the remaining reactant attributed to each type of reaction occurring in the medium as the reactions progress. $S_{ec}$ is the heat due to the internal short reaction and it can be determined [13]:

$$S_{ec} = H_{ec} \left| \frac{dS_{oc}}{dt} \right| \quad (17)$$

$$\frac{dS_{oc}}{dt} = -ISC \times Soc \times A_{ec} \times \exp(-E_{ec} / RT) \quad (18)$$

where ISC is the internal short-circuit factor used to make sure that the reaction is activated only after reaching the onset temperature (0 for $T < T_{trigger}$; 1 for $T \ge T_{trigger}$). Since the cell used in this study is an LG HG2 with an NCM-type cathode. The constants and parameters used for modeling the heat
generation are based on data from similar cathods [14]. Physically, the onset temperature corresponds to the melting temperature of the separator.

## Coupling gas and particles ejected from the cell

Gas and particles vented by the cell during TR are coupled with the internal thermal propagation within the cell. This coupling is modeled, at each time step, using a User-Defined Function (UDF), which dynamically updates the quantities of ejected gases and particles. Their mass flow rates are computed as:

$$ \dot{m}_{gas} = R_e \times m_{tot(gas)} \qquad (19) $$

$$ \dot{m}_{particles} = R_e \times m_{tot(particles)} \qquad (20) $$

Additionally, particles and vented gas are initially injected through the external surface of the cell at the average cell temperature, which is calculated and reported at each time step.

During TR, the cell loses approximately half of its mass due to the gases and particles ejection. To account for this mass loss while keeping the volume constant, the cell density must be updated dynamically. The calculations are performed through the UDF as follows:

$$ \rho = \left( \frac{m_{i(cell)} - \int_{t_i}^{t_f} \left( \frac{\partial m_{i(gas)}}{\partial t} \right) dt - \int_{t_i}^{t_f} \left( \frac{\partial m_{i(particles)}}{\partial t} \right) dt}{m_{i(cell)}} \right) * \rho_i \qquad (21) $$

## Reaction of the gas mixture outside the cell

Mass, momentum, energy, and species transport equations are solved in order to simulate the combustion of the gas vented by the cell:

### Mass conservation

$$ \frac{\partial \rho}{\partial t} + \nabla \cdot (\rho \vec{v}) = S_{mass} \qquad (22) $$

Mass conservation equation includes a source term $S_{mass}$, which represents the rate of mass injected into the domain per unit volume. This term accounts for the amount of the vented gases ejected by the cell.

### Momentum conservation

$$ \frac{\partial (\rho \vec{v})}{\partial t} + \nabla \cdot (\rho \vec{v} \vec{v}) = -\nabla p + \nabla \cdot \tau + S_{momentum} \qquad (23) $$

Boundary conditions are implemented using standard wall functions to account for near-wall turbulence, while the momentum equation's source term ($S_{momentum}$) includes contributions from body forces, pressure gradients, and viscous forces.

### Energy conservation

$$ \frac{\partial (\rho e)}{\partial t} + \nabla \cdot (\rho e \vec{v}) = \nabla \cdot (k \nabla T - \sum_{i}^{k} h_i J_i) + S_{energy} \qquad (24) $$

Where $S_{energy}$ is the source term in the energy equation, $e$ denotes the energy, $\nabla \cdot (\sum_{i}^{k} h_i J_i)$ is the transport of enthalpy due to species diffusion.
## Species Transport

$$ \frac{\partial(\rho_i)}{\partial t} + \nabla \cdot (\rho_i \vec{v}) = \nabla \cdot J_i + \dot{\omega}_i \quad (25) $$

The reaction source term $\dot{\omega}_i$, which accounts for the production of species i due to chemical reactions. The k-ω SST turbulence model is employed, where the stress tensor in the momentum equation is calculated using the turbulent viscosity derived from the model. This approach combines the k-ω model near the walls with the k-ε model in the free stream for enhanced accuracy.

## Combustion model

The Finite-Rate/Eddy-Dissipation (FR/ED) combustion model is employed to simulate the combustion process during TR. This combined model computes the reaction rates for both the Finite-Rate (FR) and Eddy-Dissipation models (EDM). The FR approach considers chemical kinetics and assumes that reactions occur at finite rates, whereas the EDM assumes that the reaction rates are controlled by turbulence and that chemical reactions occur instantly under turbulent mixing. To capture the complete combustion dynamics, the minimum reaction rate from both models is used, ensuring that the lower bound dominates in controlling the reaction kinetics. The FR/ED model has been applied in a wide range of scenarios to estimate the heat of reaction and combustion rate [15]. In the present study, five main species are considered based on previous studies on the composition of gases released during TR [16-17] These species include: carbon monoxide (CO), methane (CH$_4$), ethylene (C$_2$H$_4$), hydrogen (H$_2$), and carbon dioxide (CO$_2$). The following four combustion reactions represent the gas composition vented from an NMC-cathode battery.

$$ \begin{aligned} & CO + 0.5O_2 \rightarrow CO_2 \\ & H_2 + 0.5O_2 \rightarrow H_2O \end{aligned} \qquad \begin{aligned} & CH_4 + 1.5O_2 \rightarrow CO + 2H_2O \\ & C_2H_4 + 2O_2 \rightarrow 2CO + 2H_2O \end{aligned} $$

## Solid particles

In addition to gas ejection, the interaction of particles expelled from the cell is considered in this model. The particles are coupled with the continuous phase of vented gases and are modeled using the Discrete Phase Model (DPM), which tracks the particles' trajectory over time. The particles are treated using unsteady particle tracking, where their movement is calculated during each time step of the fluid flow simulation. This approach enables a mechanical description of particle behavior as they are expelled from the cell during TR. In terms of physical modeling, the drag law (for spherical particles) is used to account for momentum exchange between the particles and the surrounding gas. The Ranz-Marshall model is applied for heat exchange between the particles and gases. Although particle sizes can vary significantly in the literature [18] this study assumes a constant particle diameter to simplify the validation process. Importantly, this assumption does not affect the primary goal of this study, which is to validate the energy released by the cell during TR. The energy released is largely governed by the heat exchange between the particles and the gas and is less sensitive to variations in particle size for the range of conditions studied here. Future work will explore particle size variations to better capture the diverse behaviors observed during TR. The focus of this work is to establish a validated baseline model, with simplifying assumptions such as constant particle diameter applied to prioritize accuracy in energy validation while deferring more complex modeling to subsequent studies.

## Geometry and mesh

Two 3D geometries were created to validate the developed model: a closed vessel (Fig. 3a) and a vessel with an opening (Fig. 4a). The closed vessel geometry consists of a 0.54-liter stainless steel volume with an internal diameter of 46 mm (Fig. 3b). An 18650 Li-ion cell, a commercial cylindrical cell with a diameter of 18 mm and a length of 65 mm, is fixed centrally near the bottom of the vessel on a PTFE plug. The cell casing, also made of stainless steel, has a shell thickness of 3 mm.
The outside wall of the closed vessel is modeled as adiabatic to account for thermal insulation. The upper part of the vessel includes a safety valve set at 200 bars to prevent rupture (Fig. 3a). However, due to its small diameter, the valve is neglected in the simulation, while the vessel total volume is conserved. The open-end vessel is placed inside a 1 m³ enclosure of air at atmospheric pressure (Fig. 4a). To validate temperature measurements experimentally, 3 thermocouples are modeled as cylinders positioned at the same level as the vessel opening, capturing temperature data at various locations.

The mesh was generated using ANSYS Fluent Meshing. A mesh sensitivity study was conducted for the closed vessel, focusing on the fully charged cell during TR, involving the creation of three meshes with different resolutions (Table 1): M1 (fine, 169,490 elements), M2 (medium, 45,995 elements), and M3 (coarse, 25,725 elements). The results demonstrated excellent agreement in energy released by the cell during TR, with differences of only 0.27% between M1 and M2, and 1.16% between M1 and M3. Similarly, the maximum pressure reached inside the vessel showed minimal variation, with differences of 0.53% (M1 vs. M2) and 6.48% (M1 vs. M3). While larger differences were observed in pressure decay rates (15.38% between M1 and M2, and 7.64% between M1 and M3), these results primarily highlight the sensitivity of transient behaviors to mesh resolution. Importantly, M2 (Fig. 3c) demonstrated sufficient accuracy for capturing key parameters like energy release and peak pressure, making it a reliable compromise between computational efficiency and result accuracy. Future studies will focus on refining mesh strategies to further enhance the transient behavior predictions, ensuring robust performance across all key metrics. These initial findings, however, provide a strong foundation for understanding the dynamics of the system. For the open-end vessel, 281,018 elements was used, incorporating localized refinement around the thermocouples for increased accuracy (Fig. 4b).

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** The closed vessel (a); view of the LIB cell inside the computational domain representing the closed vessel (b); Mesh M2 and Refinement around vessel wall (c).</figcaption></figure>

## Model input

The pressure-velocity coupling is done using a coupled scheme. Discretizing all the solution variables is done using a second-order upwind. The first order implicit is used for the transient formulation. The targeted residuals for energy, momentum, species, and continuity are set to be $10^{-7}$, $10^{-3}$, $10^{-3}$, and $10^{-5}$ respectively. A standard $k-\omega$ SST turbulence model is used to model turbulence. As mentioned previously, the thermal runaway heat generation is described using four equation abuse model. CFD is solved using the finite volume method using the commercial CFD software, ANSYS FLUENT® 2023R1. For the purposes of this simulation, particles were assigned a constant diameter of $10^{-5}$ m.
Simulations were initially conducted at ambient temperature (20°C) and atmospheric pressure, the external surface of the cell is adjusted to achieve a temperature increase at a rate of 6 °C.min⁻¹ until the triggering temperature is reached. Triggering temperatures for different SOCs were determined experimentally. Simulations were performed using 20 processors on an Intel(R) Xeon(R) Gold 6242R. The total computing time was approximately one hour for the closed vessel, and four hours for the open-end vessel. To optimize computing time, initially, only the energy equation was applied. Species and flow equations were introduced a few time steps before reaching the triggering temperature. Resistive heating was deactivated when reaching the triggering temperature. Several simplifying assumptions were made to enhance computational efficiency, including neglecting the thickness of the resistive heating wire and disregarding changes in the cell volume during TR. Experimentally, the vessel was securely positioned within a box surrounded by vermiculite, a recognized insulating material. Therefore, the simulated vessel was supposed as adiabatic with no heat exchange with its surroundings.

The key inputs of the model include the geometric and thermal parameters of the cell (including the venting area), the kinetic and thermodynamic parameters of the degradation reactions (including ISC), the total mass and composition of the gas, as well as the total mass, average diameter, and density of the particles, and the environmental setup (e.g., thermal insulation for the closed configuration). The outputs include temperature profiles, pressure evolution, gas flow dynamics, and energy release during TR. For the closed configuration, the model provides insights into pressure buildup, heat transfer, and gas ejection behavior under contained conditions. For the open configuration, it focuses on venting dynamics and the spatial distribution of expelled gas and heat. These outputs are compared directly with experimental data such as thermocouple measurements. Moreover, the model allows exploration of phenomena that are challenging to measure experimentally, such as detailed gas flow dynamics, localized heat transfer mechanisms, and species distribution during TR.

<figure><img src="image_5.png" /><figcaption>Fig. 4. Geometry and mesh of the open vessel, (a) and (b), respectively.</figcaption></figure>

<table><caption>Table 1 Mesh parameters of the computational model for the closed vessel.</caption><thead><tr><td>Parameter</td><td>Fine mesh (M1)</td><td>Medium mesh (M2)</td><td>Coarse mesh (M3)</td></tr></thead><tbody><tr><td>Cells</td><td>169490</td><td>45995</td><td>25725</td></tr><tr><td>Aspect ratio</td><td>&lt; 30.7</td><td>&lt; 31,24</td><td>&lt; 5.8</td></tr><tr><td>Orthogonal quality</td><td>&gt; 0.21</td><td>&gt; 0,33</td><td>&gt; 0.45</td></tr></tbody></table>

# RESULTS AND DISCUSSION

## Closed vessel

The total number of moles of gas released during TR was measured experimentally in the case of the closed vessel for different SOC levels (Table 3) and used as an input in the calculation. Using the initial conditions of the closed vessel—0.54 L of air at ambient temperature and atmospheric
pressure—the number of moles inside the closed vessel is calculated using the ideal gas law, yielding 0.022 moles. Each gas species fraction for different SOC levels is taken from [16]. **Table 2** presents the mole (Xi) and mass fractions of each gas species released by the cell for a fully charged state of charge. The same conditions were applied for various SOC levels (15% to 100%) considering the total mass of gaseous species reported from experiments in **Table 3**. The data obtained under these conditions and used for the simulations are summarized in **Table 3**.

<table><caption>**Table 2** Calculated mass of each species ejected by the cell, mole numbers obtained from experiments, and gas fractions assumed from Xu et al. [16].</caption><thead><tr><th>SOC100</th><th>CO2</th><th>CO</th><th>H2</th><th>C2H4</th><th>CH4</th><th>O2</th><th>Total</th></tr></thead><tbody><tr><td>Mole fraction Xi [16]</td><td>0.3962</td><td>0.2435</td><td>0.2163</td><td>0.0908</td><td>0.0329</td><td>0.0204</td><td>1</td></tr><tr><td>Molar mass (g.mol<sup>-1</sup>)</td><td>44.01</td><td>28.01</td><td>2.02</td><td>28.05</td><td>16.04</td><td>32.00</td><td>nn</td></tr><tr><td>Mass fraction Yi</td><td>0.61</td><td>0.24</td><td>0.02</td><td>0.09</td><td>0.02</td><td>0.02</td><td>1</td></tr><tr><td>Mass of gaseous species [g]</td><td>4.53</td><td>1.77</td><td>0.11</td><td>0.66</td><td>0.14</td><td>0.17</td><td>7.39</td></tr></tbody></table>

<table><caption>**Table 3** Cell mass data and mass loss fractions utilized for conducting simulations at various SOC levels.</caption><thead><tr><th>SOC (%)</th><th>15</th><th>20</th><th>30</th><th>50</th><th>75</th><th>100</th></tr></thead><tbody><tr><td>Experimental moles of gas released</td><td>0.04</td><td>0.0495</td><td>0.0468</td><td>0.08</td><td>0.194</td><td>0.26</td></tr><tr><td>Cell initial mass (g) (exp.)</td><td></td><td></td><td></td><td>46</td><td></td><td></td></tr><tr><td>Cell final mass (g) (exp.)</td><td>10.80</td><td>18.48</td><td>12.96</td><td>11.00</td><td>24.43</td><td>22.20</td></tr><tr><td>Total mass of gaseous species (g)</td><td>1.37</td><td>1.70</td><td>1.53</td><td>2.09</td><td>5.45</td><td>7.39</td></tr><tr><td>Total mass of particles (g)</td><td>33.82</td><td>25.81</td><td>31.50</td><td>32.90</td><td>16.12</td><td>16.41</td></tr></tbody></table>

The energy and other key parameters associated with thermal runaway are evaluated. For the closed vessel, **Fig. 6** shows a comparison between the experimental and calculated energy released by the cell for different SOC. E_max represents the maximum energy released by the cell, including the energy supplied by the thermal heater. E_TR is the difference between E_max and the energy supplied by the thermal heater, representing the energy released due to TR. The simulated energy released during thermal runaway events ranged from 23.4 kJ to 78.5 kJ which is in good agreements with the experimental data, moreover the non-linear variation of the energy release observed experimentally is reproduced by the simulation. These results demonstrate that incorporating the full set of relevant phenomena, along with their interactions, is sufficient to accurately predict the energy release during thermal runaway, with a relative error, between experiments and simulations, of less than 9%. Kinetic energy can be disregarded when compared to the thermal energy released by the cell. Consequently, the shape of the cell cap was neglected, and injection was assumed to occur from the upper external surface of the cell, with a diameter of 12 mm. **Fig. 5** shows the velocity contour of the gas ejected during TR from a cross-sectional side view of the vessel, with SOC=100.

The pressure dynamics inside the closed vessel during TR were analyzed for various SOC (Fig. 7) (100%, 75%, 50%, and 30%) using both experimental measurements and CFD simulations. The results indicate that the peak pressure increases significantly with SOC, with the highest pressure observed at SOC 100% and the lowest at SOC 30%. This trend is attributed to the greater availability of reactive material at higher SOCs, leading to more vigorous exothermic reactions and gas generation during TR. Post-peak pressure showed a gradual decline, with faster dissipation observed at lower SOCs due to reduced gas production. CFD simulations closely aligned with experimental data, particularly in capturing peak pressures and overall trends, though minor discrepancies highlight potential areas for refinement, such as accounting for secondary reactions and material heterogeneities. These findings underscore the critical role of SOC in determining the severity of TR events and emphasize the importance of designing robust containment systems to handle high-pressure scenarios, particularly for cells at high SOC.
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Velocity profile (m/s) - SOC=100%; Cross-Sectional View – without particles.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Fig. 6.** A Comparison between the calculated and measured (exp.) energy released by the LIB cell during TR for different SOC.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Calculated vs measured pressure inside the vessel for different SOC levels.</figcaption></figure>

## Open-end vessel

For the open-end vessel experiments, two repeated tests were performed: Reference 1 and Reference 2. **Fig. 8** illustrates the temperature evolution for the 3 thermocouples, TC0, TC1, and TC2, shown in panels (a), (b), and (c), respectively. The first thermocouple, TC0, positioned 1.5 cm away from the vessel opening, recorded a maximum temperature elevation of 1055°C in Reference 2, more than twice the value obtained in the CFD simulation, which was 580°C. This discrepancy indicates that there is still an unaccounted amount of energy that could be better represented, specifically in terms of heat exchange between the ejected matters and the surrounding environment. One potential explanation lies in the assumption that no $O_2$ remains in the ejecta, as none is measured experimentally after TR. This assumption leads to a non-premixed flame with possible detachment, which differs from the experimental observations where no detachment is present. This suggests that the assumption may be incorrect and could be revisited. In this study, solid particles ejected by the cell were assumed to have initially the average temperature of the cell. Consequently, the energy released due to the oxidation of these particles was neglected, resulting in a lower energy prediction. Future work will focus on studying this phenomenon to improve the modeling of thermal heat exchange around the cell. For the second thermocouple, TC1, located 0.2 m from the vessel opening, the numerical and experimental results (shown in **Fig. 8b**) exhibit a correct agreement, with a negligible relative error compared to Reference 1. This suggests that the combustion of gases accounted for in the simulations is sufficiently accurate to predict temperatures at this distance from the vessel cap. In addition, for TC2, positioned 0.68 m from the cap, a relative error of 6% was observed between the peak CFD and measured temperatures in Reference 2. The temperature rise measured by TC1 and TC2 aligns well with the simulation data, indicating acceptable agreement. However, it was observed that the experimental cooling rate of the thermocouples after flame exposure was slower than in the CFD simulations. This discrepancy may be attributed to several factors, including the thermal mass of the thermocouples, differences in heat transfer mechanisms (such as convection and radiation), contact
resistance, and the response time of the thermocouples. Additionally, interactions between incandescent particles ejected during thermal runaway and the thermocouples, as observed in high-speed camera footage, could have influenced the measurements. Environmental surroundings and limitations of the measurement equipment could have contributed to these differences. Refining the boundary conditions and incorporating more accurate material properties in the simulations could further improve the agreement between the experimental and numerical results.

<figure><img src="image_9.png" /><figcaption>Fig. 8. Plots of the temperature of the thermocouples for TC0 (a), TC1 (b), and TC2 (c) versus time.</figcaption></figure>

## CONCLUSIONS AND PERSPECTIVES

This study investigated thermal runaway events in 18650 lithium-ion batteries with NMC cathodes under two different geometries: a closed and an open-end vessel. For the closed vessel, which had a volume of 0.54 L, six simulations were performed, demonstrating a relative error of less than 9% across a wide range of states of charge (SOC) from 15% to 100% on the prediction of the energy release. The energy release predictions, ranging from 23.4 kJ to 78.5 kJ, highlighted the model capacity to accurately capture the dynamics of these hazardous events when critical thermal and chemical phenomena were fully integrated and considered. The model's primary inputs include the cell physical dimensions, material properties (thermal and chemical), initial and boundary conditions such as vent opening geometry and triggering temperature, and the SOC. The outputs generated by the model include the energy release during TR, pressure evolution, and temperature profiles at various locations. These outputs are directly comparable with experimental data such as thermocouple measurements and pressure readings, allowing for a thorough validation of the model. Key assumptions made in the model include neglecting particle oxidation and assuming uniform particle diameter during venting, which were critical in optimizing computational efficiency but could influence certain aspects of the results. The second geometry studied was the open-end configuration. In this setup, experimental data were collected using thermocouples at various distances from the cell to analyze the temperature evolution. The results revealed that the first thermocouple (TC0), positioned 0.28 m from the cap of the cell, exhibited a significant temperature discrepancy between experimental and numerical results, possibly due to the underestimation of the gas combustion or to the particles oxidation. In contrast, the second thermocouple (TC1), located 0.48 m from the cap of the cell, showed good agreement with negligible relative error, validating the model's combustion predictions. The third thermocouple (TC2), situated 0.68 m from the opening, demonstrated a relative error of 6%, indicating reasonable predictive accuracy. Overall, while the agreement between simulations and experimental data was acceptable, discrepancies in the cooling rates were observed, likely due to differences in thermal mass, heat transfer mechanisms, and environmental factors affecting the thermocouples.

Future work will focus on refining heat exchange modeling, particularly by incorporating the energy contribution from particle oxidation, to further enhance the accuracy of the simulations and improve the understanding of thermal runaway dynamics in various containment scenarios. Further exploration of parameters such as the composition of the vented gases and the efficiency of the combustion model
could be further improved by incorporating oxygen release, these improvements could provide deeper insights into the dynamics of thermal runaway and improve the model's applicability to various cell systems. Additionally, extending the model to full battery modules will enable risk assessments at a larger scale.

## REFERENCES

[1] ‘Tesla Vehicle Safety Report | Tesla’. Accessed: Nov. 05, 2024. [Online]. Available: https://www.tesla.com/VehicleSafetyReport

[2] ‘Samsung Recalls Galaxy Note7 Smartphones Due to Serious Fire and Burn Hazards’, U.S. Consumer Product Safety Commission. Accessed: Nov. 05, 2024. [Online]. Available: https://www.cpsc.gov/Recalls/2016/Samsung-Recalls-Galaxy-Note7-Smartphones

[3] ‘Lithium Battery Incidents | Federal Aviation Administration’. Accessed: Nov. 05, 2024. [Online]. Available: https://www.faa.gov/hazmat/resources/lithium_batteries/incidents

[4] S. Abada, M. Petit, A. Lecocq, G. Marlair, V. Sauvant-Moynot, and F. Huet, ‘Combined experimental and modeling approaches of the thermal runaway of fresh and aged lithium-ion batteries’, *J. Power Sources*, vol. 399, pp. 264–273, 2018.

[5] P. Huang, H. Chen, A. Verma, Q. Wang, P. Mukherjee, and J. Sun, ‘Non-dimensional analysis of the criticality of Li-ion battery thermal runaway behavior’, *J. Hazard. Mater.*, vol. 369, pp. 268–278, May 2019.

[6] T. Katrašnik, I. Mele, and K. Zelič, ‘Multi-scale modelling of Lithium-ion batteries: From transport phenomena to the outbreak of thermal runaway’, *Energy Convers. Manag.*, vol. 236, p. 114036, May 2021.

[7] D. Kong, G. Wang, P. Ping, and J. Wen, ‘A coupled conjugate heat transfer and CFD model for the thermal runaway evolution and jet fire of 18650 lithium-ion battery under thermal abuse’, *eTransportation*, vol. 12, p. 100157, May 2022.

[8] R. Zhao, Z. Lai, W. Li, M. Ye, and S. Yu, ‘Development of a coupled model of heat generation and jet flow of lithium-ion batteries during thermal runaway’, *J. Energy Storage*, vol. 63, p. 107048, Jul. 2023.

[9] G. Wang et al., ‘Revealing particle venting of lithium-ion batteries during thermal runaway: A multi-scale model toward multiphase process’, *eTransportation*, vol. 16, p. 100237, Apr. 2023.

[10] ‘Ansys Fluent | Fluid Simulation Software’. Accessed: Aug. 05, 2024. [Online]. Available: https://www.ansys.com/products/fluids/ansys-fluent

[11] T. D. Hatchard, D. D. MacNeil, A. Basu, and J. R. Dahn, ‘Thermal Model of Cylindrical and Prismatic Lithium-Ion Cells’, *J. Electrochem. Soc.*, vol. 148, no. 7, p. A755, 2001.

[12] G.-H. Kim, A. Pesaran, and R. Spotnitz, ‘A three-dimensional thermal abuse model for lithium-ion cells’, *J. Power Sources*, vol. 170, no. 2, pp. 476–489, Jul. 2007.

[13] P. T. Coman, E. C. Darcy, C. T. Veje, and R. E. White, ‘Modelling Li-Ion Cell Thermal Runaway Triggered by an Internal Short Circuit Device Using an Efficiency Factor and Arrhenius Formulations’, *J. Electrochem. Soc.*, vol. 164, no. 4, p. A587, Jan. 2017.

[14] P. Peng and F. Jiang, ‘Thermal safety of lithium-ion batteries with various cathode materials: A numerical study’, *Int. J. Heat Mass Transf.*, vol. 103, pp. 1008–1016, Dec. 2016.

[15] B. F. Magnussen and B. H. Hjertager, ‘On mathematical modeling of turbulent combustion with special emphasis on soot formation and combustion’, *Symp. Int. Combust.*, vol. 16, no. 1, pp. 719–729, Jan. 1977.

[16] L. Xu et al., ‘Thermal runaway propagation behavior and gas production characteristics of NCM622 battery modules at different state of charge’, *Process Saf. Environ. Prot.*, vol. 185, pp. 267–276, May 2024.

[17] M. Qiu, J. Liu, B. Cong, and Y. Cui, ‘Research Progress in Thermal Runaway Vent Gas Characteristics of Li-Ion Battery’, *Batteries*, vol. 9, no. 8, p. 411, Aug. 2023.

[18] W. Li et al., ‘Enhancing understanding of particle emissions from lithium-ion traction batteries during thermal runaway: An overview and challenges’, *eTransportation*, vol. 22, p. 100354.