Electro-chemo-mechanical modelling of structural battery composite full cells

Preliminaries

The analysis considers a synthetical microscale volume element Ω = Ω_CF ∪ Ω_SBE ∪ Ω_SEP ∪ Ω_P subdivided into several domains, shown in Fig. 5. Ω_CF represents the carbon fibre domain in the negative electrode, described by the union of carbon fibres in the volume element. The carbon fibres are embedded in a structural battery electrolyte matrix denoted Ω_SBE, where Ω_N = Ω_CF ∪ Ω_SBE is referred to as the negative electrode. Redox reactions occur on the carbon fibre/SBE interfaces, Γ_CF = ∂Ω_CF ∩ ∂Ω_SBE. To prevent short circuits, the positive and negative electrodes are separated by a porous SBE-impregnated glass fibre separator, Ω_SEP. The interface between the negative electrode and the separator is denoted \({\Gamma }_{{\rm{N}}}^{\,\text{int}\,}=\partial {\Omega }_{{\rm{N}}}\cap \partial {\Omega }_{{\rm{SEP}}}\). For convenience, we define the electrolyte domain (ED) where ionic transport takes place as Ω_ED = Ω_SEP ∪ Ω_SBE. Ω_P denotes the homogenised positive electrode domain, where \({\Gamma }_{{\rm{P}}}^{{\rm{int}}}=\partial {\Omega }_{{\rm{SEP}}}\cap \partial {\Omega }_{{\rm{P}}}\) describes the interface between the positive electrode and separator domains. The cell thickness, h_cell, is roughly 390 μm, whereas the negative electrode has a thickness of h_N = 50 μm. The thickness of the Whatman GF/A separator is approximately h_SEP = 260 μm, while the positive electrode has a thickness of h_P = 80 μm^1,6. The positive electrode is conceptualised as a continuum, implying that the underlying structure is not resolved. It is acknowledged that redox reactions occur on each interface of the active LFP particles distributed at various positions inside Ω_P. To simplify the modelling approach, the redox reactions in the positive electrode are modelled on the \({\Gamma }_{{\rm{P}}}^{{\rm{int}}}\) interface. Through this simplification, the interface kinetics are clearly structural, and will vary with, e.g., the thickness of the positive electrode. Moreover, we assume that both the positive and negative electrodes function as ideal conductors. This implies that all carbon fibres are electrically connected to a current collector, establishing the global potential of the negative electrode, denoted Φ_N inside Ω_CF. The same rationale applies to the positive electrode characterised by the global potential Φ_P inside Ω_P.

The alternative text for this image may have been generated using AI.

Full size image

Schematic two-dimensional illustration of the structural battery domains and boundaries.

In addition to the interfaces already discussed, Fig. 5 illustrates an assumed periodic structure in the 2-direction. The right side \({\Gamma }^{+}={\Gamma }_{{\rm{N}}}^{+}\cup {\Gamma }_{{\rm{SEP}}}^{+}\cup {\Gamma }_{{\rm{P}}}^{+}\) has a corresponding mirror side \({\Gamma }^{-}={\Gamma }_{{\rm{N}}}^{-}\cup {\Gamma }_{{\rm{SEP}}}^{-}\cup {\Gamma }_{{\rm{P}}}^{-}\). For the subsequent introduction of periodic boundary conditions, we introduce the mapping of points φ_PER: Γ⁺ → Γ⁻ such that f⁻(x⁺) = f(φ_PER(x⁺)) and f⁺(x⁺) = f(x⁺) for any function f(x) on Γ⁺ ∪ Γ⁻.

Balance equations

In the following, we formulate the boundary value problem for the domains shown in Fig. 5. The balance of linear momentum and relevant boundary conditions applicable to the negative electrode domain is stated as

where σ is the Cauchy stress tensor, u is the displacement field and n is the (from Ω_N) outwards facing unit normal.

The ionic transport of lithium ions and the accompanying anions, here denoted X, inside Ω_ED is governed by the conservation of mass and Gauss law as follows:

where c_α and j_α denote the concentration and flux, respectively, of lithium ions (α = Li) and the accompanying anion (α = X). μ_α denotes the chemical potential for the ions (α = Li, X), d is the electric flux density, F is Faraday’s constant and φ is the electric potential. Here, n denotes the outwards pointing normal (from Ω_ED). Hence, during charging, j_P denotes the flux of Li-ions entering the electrolyte from the positive electrode, and j_N denotes the flux of Li-ions moving into the negative electrode. Similarly, d_P and d_N denote the electric flux density moving from the positive and to the negative electrode, respectively.

Finally, the transport of charge-neutral lithium inside the electrodes is governed by mass conservation as

where c_Li and j_Li denote the concentration and flux, respectively, of neutral lithium atoms. Here, n is the outwards pointing unit normal (from Ω_CF/Ω_P). Hence, the conditions in Eqs. (13) and (14) combined with Eqs. (8) and (9) pertain to mass conservation over the interfaces \({\Gamma }_{{\rm{P}}}^{{\rm{int}}}\) and Γ_CF. Across the interfaces \({\Gamma }_{{\rm{p}}}^{{\rm{int}}}\) and Γ_CF, we introduce the jump operators

for the normal n_Γ associated with the convention for j_P, j_N. Hence, n_Γ points from Ω_P to Ω_SEP and from Ω_SBE to Ω_CF.

When solving the pertinent initial boundary value problem, we introduce the spatially constant electric potentials Φ_N in Ω_CF and Φ_P in Ω_P. The total current density over the positive electrode is evaluated at the interface \({\Gamma }_{{\rm{P}}}^{{\rm{int}}}\) as

where I_app is the prescribed current for a galvanostatic (current control) process. Furthermore, the initial, stress-free reference state is given at \({c}_{{\rm {Li}}}={c}_{{\rm{CF}}}^{0}\) for \({\boldsymbol{x}}\in {\Omega }_{{\rm{CF}}},{c}_{{\rm {Li}}}={c}_{{\rm{P}}}^{0}\) for x ∈ Ω_P and c_Li = c_X = c_ref for x ∈ Ω_ED at t = 0.

Time-incremental weak format of the full cell problem, current control

Upon employing the Backward–Euler rule, the galvanostatic problem for controlled current I_app(t) is stated as follows: For known values of ⁿ⁻¹c_Li, ⁿ⁻¹c_X, we seek the spatial fields at time t = t_n: \({\boldsymbol{u}},\varphi ,{\mu }_{{\rm {Li}}},{c}_{{\rm {Li}}},{\mu }_{X},{c}_{X},{\Phi }_{{\rm{P}}}\,\in {\mathbb{U}}\times {\mathbb{F}}\times {{\mathbb{M}}}_{{\rm {Li}}}\times {{\mathbb{M}}}_{X}\times {{\mathbb{L}}}_{2}({\Omega }_{{\rm{ED}}})\times {\mathbb{R}}\). We have used the notation y(x): = ⁿy(x) = y(x, t_n) and ⁿ⁻¹y(x) = y(x, t_n−1).

The relevant solution (and test) spaces are defined as

where \({{\mathbb{L}}}_{2}(w)\) and \({{\mathbb{H}}}^{1}(w)\) denotes the space of square integrable functions, and the space of functions with square integrable derivatives of order 0 and 1, respectively.

Eq. (19) corresponds to the balance of linear momentum along with boundary conditions shown in Eqs. (1)–(3). Eq. (20) correspond to Gauss law shown in Eq. (5) combined with boundary conditions in Eqs. (8)–(10). The weak representation of Li-mass balance in Eq. (21) originates from its strong form counterparts in Eqs. (4) and (12), along with the relevant boundary conditions in Eqs. (6), (8)–(10) and (13)–(16). Similarly, the weak representation of the anion mass balance in Eq. (22), is related to its strong counterpart in Eq. (4) with boundary conditions in Eqs. (6) and (11). Finally, Eqs. (23) and (24) enforce the chemical potential, μ_α, to follow the constitutive relation \({\mu }_{\alpha }^{{\rm{en}}}={\mu }_{\alpha }^{{\rm{en}}}({c}_{\alpha })\) and Eq. (25) forces the electric current over \({\Gamma }_{{\rm{P}}}^{{\rm{int}}}\) to follow the prescribed current in I_app. We recall that the negative electrode potential is set to (reference) Φ_N = 0. In the following sections, we introduce constitutive relations for \({\boldsymbol{\sigma }},{\mu }_{\alpha }^{{\rm{en}}},{{\boldsymbol{j}}}_{\alpha },{j}_{{\rm{N}}},{j}_{{\rm{P}}},{\boldsymbol{d}},{d}_{{\rm{N}}}\) and d_P.

Constitutive relations for the carbon fibre domain in the negative electrode

Lithiation of carbon fibres are accompanied by a change in moduli and an anisotropic expansion characterised by α_CF, containing longitudinal and transverse expansion coefficients α_L and α_T, respectively^6,15. The Helmholtz free energy for the carbon fibres is expressed as

where E_CF(c_Li) represents the fourth order (anisotropic) stiffness tensor of the carbon fibres described in Larsson et al.¹⁶. The relevant material parameters related to E_CF(c_Li) are the following moduli and Poisson ratios E_L(c_Li), E_T(c_Li), ν_LT, ν_TT, G_LT and G_TT. Furthermore, ϵ(u) is the (small) strain tensor. The maximum (reversible) concentration can be defined as the sum of inserted lithium concentration and the concentration of vacant sites, \({\widetilde{c}}_{{\rm{CF}}}={c}_{v}+{c}_{{\rm {Li}}}\), where c_v is the concentration of vacant sites for possible lithium insertion. In this work, we consider conditioned carbon fibres, assuming the maximum concentration is constant. Consequently, no Li-degradation mechanisms in the electrode such as lithium trapping, dendrite growth, or SEI formation are considered. The chemical contribution to the free energy is expressed as

where R is the ideal gas constant and T is the temperature. Λ(c_Li, c_v) represents the deviation from ideal conditions caused by interaction of inserted lithium and the vacant sites. We adopt a formulation similar to Landstorfer et al.¹⁷, with the solubility parameter k_CF so that

We can now express the following constitutive relations

where the carbon fibre reference chemical potential for the carbon fibres is obtained as \({\mu }_{{\rm{CF}}}^{0}={\mu }_{{\rm {Li}},{\rm{CF}}}^{0}-{\mu }_{v,{\rm{CF}}}^{0}\). \({\mu }_{{\rm {Li}},{\rm{CF}}}^{0}\) and \({\mu }_{{\rm{v}},{\rm{CF}}}^{0}\) are the reference chemical potentials for lithium and vacancies, respectively. The contribution from the concentration dependent stiffness is ignored. The mass flux of neutral lithium within the carbon fibres is governed by the gradient of chemical potential scaled by an isotropic mobility

The mobility along the fibre (1-dir) is not equal to the mobility in the transversely isotropic plane (2-3-dir). Here, we assign the isotropic mobility equal to the mobility in the 2-3-dir.

Constitutive relations for the structural battery electrolyte domain

The porous, bi-phasic structural battery electrolyte matrix enables ion transport as well as mechanical load transfer. Recent studies carried out by Duan et al.¹⁸ show that the pores are homogeneously distributed and that the SBE is isotropic. We model the SBE using linear isotropic elasticity and introduce the Helmholtz free energy for the SBE domain as

where E_SBE is the homogeneous isotropic fourth-order stiffness tensor with related shear and bulk moduli G and K. c_ref is the salt concentration in the liquid electrolyte phase, ε is the electric permittivity, ∇φ the electric field and \({\mu }_{\alpha }^{0}\) the reference chemical potential of lithium ions (α = Li) and the anions (α = X). The following (decoupled) constitutive relations can be derived

The mass flux of Li ions and accompanying anion are expressed as

where η_Li and η_X are the isotropic mobility coefficients of lithium ions and anions, respectively. We note that both chemical and electric potential gradients contribute to ionic transport, i.e. diffusion and migration. The electric field affects the ionic transport of the cation and anion in opposite directions due to the different charge of the ions. The mobilities of Eqs. (41) and (42) follow from the assumption of constant diffusion coefficients in the electrolyte.

Constitutive relations for the separator domain

Similar to the SBE, the same equations can be applied over the separator domain. The separator consists of a mixture of SBE and isotropic glass fibre separator. The volume fraction of glass fibre separator is estimated as 34% using available data from Asp et al.². The effective mobilities and the elasticity tensor are obtained using volume averaging, where we assume that the glass fibres block ionic transport. The effective mobilities thus become \({\bar{\eta }}_{{\rm {Li}}}=(1-{V}_{{\rm{GF}}}){\eta }_{{\rm {Li}}},{\bar{\eta }}_{X}=(1-{V}_{{\rm{GF}}}){\eta }_{X}\) and E_SEP = (1− V_GF)E_SBE + V_GFE_GF, where V_GF is the volume fraction of glass fibres in the separator. The Helmholtz free energy applicable to Ω_SEP is expressed as

The constitutive relations are expressed as

Constitutive relations for the positive electrode domain

We adopt the same formulation as introduced in the section “Introduction” for the positive electrode. It is important to note that data pertaining to the homogeneous expansion of the positive electrode is unavailable. This choice effectively disables the coupling of chemo-mechanical fields in the positive electrode, and we consider the positive electrode as stress-free upon lithium insertion. It is acknowledged that α_P can be obtained through measurements or by employing computational homogenisation^16,19. The elastic properties of a positive electrode slurry were measured by Gupta et al.²⁰. They reported a tensile modulus of 0.90 GPa for the NMC based positive electrode continuum. We base the elastic properties of the homogenised positive electrode on these measurements. Additionally, the maximum concentration and solubility parameter of the positive electrode are denoted as \({\widetilde{c}}_{{\rm {P}}}\) and k_P, respectively. The Helmholtz free energy of the positive electrode is expressed as

where

and

The constitutive relations in the positive electrode can be expressed as

The mass flux in the positive electrode is characterised by the effective diffusion coefficient D_Li,P, representing the homogenised underlying transport mechanisms in the positive electrode slurry. The mobility of neutral lithium in the positive electrode is derived assuming a constant diffusion coefficient such that the mass flux becomes

where \(\frac{\partial {\mu }_{{\rm{Li}}}}{\partial c}\) is computed from Eq. (54).

Constitutive relations for electrode–electrolyte interfaces

Figure 6 shows a generic electrode/electrolyte interface where Γ describes the transition from the electrode (− side) to the electrolyte (+ side). The interface mass flux of lithium ions is denoted j_n = j_Li ⋅ n_Γ, where n_Γ is the unit normal pointing from (−) to (+). We consider the case where no lithium ions accumulates on Γ, whereby \({j}_{{\rm{n}}}^{-}={j}_{{\rm{n}}}^{+}={j}_{\rm{{n}}}\). Furthermore, we do not allow accumulation of free charge on the boundary, whereby the current \({i}_{{\rm{n}}}^{-}={i}_{{\rm{n}}}^{+}={i}_{{\rm{n}}}\) also becomes continuous. Since there is no flux of anions across Γ, we conclude that

Assuming no net free charge on the interface, we also obtain continuity in the electric flux densities \({d}_{{\rm{n}}}^{-}={d}_{{\rm{n}}}^{+}={d}_{\rm{{n}}}\). Finally, for an isotropic ideal conductor, the electric flux density can be expressed as

at the electrode side, where ε_c and κ_c are the electric permittivity and conductivity of the homogeneous electrode. Using Eqs. (56) and (57) for the continuous fluxes, we arrive at

Hence, j_n is the only remaining constitutive relation needed on Γ.

The alternative text for this image may have been generated using AI.

Full size image

Illustration of a generic electrode/electrolyte interface.

The interface mass flux over the positive and negative electrode/electrolyte interfaces is described by the Butler–Volmer relation, where the mass flux is continuous over the interface while the electric and chemical potential fields are discontinuous. The interface relation over positive electrode and separator interface is expressed as

where the overpotential

and the interface electric flux density is expressed as

Here, i_0,P(c_Li) is the exchange current density parametrised in the electrode concentration c_Li, pertinent to the positive electrode/separator interface expressed as a third degree polynomial with coefficients q_i so that \({i}_{0,{\rm{P}}}({c}_{{\rm{Li}}})={q}_{0}+{q}_{1}{c}_{{\rm{Li}}}+{q}_{2}{c}_{{\rm{Li}}}^{2}+{q}_{3}{c}_{{\rm{Li}}}^{3}\). Furthermore, assuming a near-constant concentration in the electrolyte. Φ_P is the electric potential, ε_P is the electric permittivity and κ_P is the electric conductivity of the positive electrode. Similarly, the interface mass flux and electric charge flux density over the carbon fibres in the negative electrode and SBE interfaces are expressed as

with the overpotential

i_0,N(c_Li) is the exchange current density of the SBE/carbon fibre interface parametrised in local carbon fibre Li-concentration, Φ_N is the electric potential in the negative electrode, ε_N is the electric permittivity and κ_N is the electric conductivity of the negative electrode.

Experimental overview

This section addresses the determination of material parameters aimed at minimising the difference between the experimental charge/discharge curves and simulated response voltage profiles for the same charge rates. Two unique cells are considered in the calibration and validation process. The first cell is manufactured in-house following the process described by Siraj et al.¹ and contains approximately 12,000 T800 carbon fibres in the negative electrode. Furthermore, the cell extends 33 mm the fibre direction and is cycled at C/20 and C/10, including a 30 min rest time. The second cell was manufactured at a different point in time by Siraj et al.¹ and contains ~24,000 fibres with a length of 50 mm in the fibre direction. The charge–discharge data from the second cell, cycled at C/2, is used to validate the calibrated parameters. The applied current related to each cell and C-rate, is normalised with the calculated weight of carbon fibres in the cell where ρ_CF = 1750 kg/m³ is the density of the carbon fibres. These currents are applied to the model and given as a function of time through the variable I_app in Eq. (25).

Identification of carbon fibre electrode properties

Kjell et al.²¹ conducted measurements of the diffusion coefficient and exchange current density of a single IMS65 carbon fibre in relation to the state of charge. To reduce the number of calibration parameters, the diffusion coefficient and exchange current density data is utilised to determine expressions for the equivalent mobility and exchange current density for the negative electrode. The measurements on the carbon fibre was performed in liquid electrolyte, which implies stress free expansion of the carbon fibre as lithium inserts.

By identifying \({M}_{{\rm {Li}},{\rm{CF}}}({c}_{{\rm {Li}}})={(\frac{\partial {\mu }^{{\rm{en}}}}{\partial {c}_{{\rm {Li}}}})}^{-1}{D}_{{\rm {Li}},{\rm{CF}}}({c}_{{\rm {Li}}})\), where \({D}_{{\rm {Li}},{\rm{CF}}}^{{\rm{M}}}\) contains discrete data and \({(\frac{\partial {\mu }^{{\rm{en}}}}{\partial {c}_{{\rm {Li}}}})}^{-1}\) is obtained from Eq. (34), we define \({D}_{{\rm {Li}},{\rm{CF}}}^{* }({c}_{{\rm {Li}}})\) as a second-degree polynomial determined by

Thereby, the mass flux in the carbon fibres is expressed as

Similarly, the exchange current density data is approximated by a third-degree polynomial determined by minimising the function

The resulting polynomial fits are shown in Fig. 7 where the solid curves represent the discrete data points and the dashed curves are the continuous functions obtained using Eqs. (67) and (69).

The alternative text for this image may have been generated using AI.

Full size image

a Diffusion coefficient discrete data from Kjell et al.²¹ and the polynomial fit as a function of normalised concentration. b Exchange current density discrete data from Kjell et al.²¹ and the related polynomial fit as a function of normalised concentration.

Furthermore, the reference carbon fibre chemical potential \({\mu }_{{\rm{CF}}}^{0}\) is obtained from measurements by Carlstedt et al.⁶, where the carbon fibre negative electrode is analysed in a half-cell configuration vs. lithium metal. The equilibrium potential of the half-cell was measured across various concentrations. We utilise the experimental measurements to uniquely determine the magnitude of the reference chemical potential \({\mu }_{{\rm{CF}}}^{0}\). In a half-cell configuration, the carbon fibres serve as the positive electrode denoted •_CF, while lithium metal acts as the negative electrode denoted •_ref. Furthermore, quantities related to the electrolyte are denoted without subscript. At net zero current, the following conditions hold:

Here, η_CF = F[φ−Φ_CF] + μ−μ_CF and η_ref = F[Φ_ref−φ] + μ_ref−μ are the overpotentials at carbon fibre / electrolyte and lithium metal/electrolyte interfaces, respectively. After relaxation, at zero current, a stationary electric potential is obtained. Additionally, the concentration of excess lithium in the electrolyte inserts into the carbon fibres, leading to c_Li = c_ref in the electrolyte. Moreover, the reference chemical potential in the electrolyte is assumed to be zero. Setting the electric potential of the negative electrode as the reference, the overpotentials at each interface can be expressed as

The chemical potential of lithium metal is not significantly affected by concentration variations, i.e., the activity of the solid phase is one, \({\mu }_{{\rm{ref}}}={\mu }_{{\rm{ref}}}^{0}\)²². From Eq. (72), η_ref = 0 gives the electrolyte potential

We adopt Eq. (34) as the model for the chemical potential in the carbon fibre. Assuming a a stress free fibre together with Eq. (73) together with η^CF = 0 in Eq. (72) gives

where c_Li is the Li-concentration in the fibre. The reference chemical potential of lithium metal phase approximated as \({\mu }_{{\rm{ref}}}^{0}=0\), see for instance Mayur et al. or work carried out by Lai et al.^23,24. At \({c}_{{\rm {Li}}}/{\widetilde{c}}_{{\rm{CF}}}=0.5\) the measured equilibrium potential was \({\Phi }_{{\rm{CF}}}^{0.5}=0.2\) V⁶.

thereby, the reference chemical potential of the fibre can be expressed as

Hence, \({\mu }_{{\rm{CF}}}^{0}\) is determined uniquely from the coefficient k_CF.

SBE conductivity

The ionic conductivity of the SBE denoted κ_SBE, can be related to the mobility coefficients of the cations and anions under the assumption that no concentration gradients are present, and c_Li ≈ c_X ≈ c_ref, is expressed as follows

by identification, the ionic conductivity in the SBE can be expressed as

where z_α is the charge number of the ion. Additionally, the transference number t⁺ relates the mobility coefficients of the cation and anion to each other, we express

In this way, both the mobility of lithium and the related anion can be uniquely determined. Cattaruzza et al.¹³ conducted measurements regarding the ionic conductivity, κ_SBE, and transference number, t⁺, for various salt concentrations and electrolyte volume fractions. They report a variation in conductivity, ranging from 0.037 mS/cm at 40 wt% electrolyte content to 0.29 mS/cm at 50 wt% electrolyte content. Meanwhile, the transference number exhibits a slight increase, moving from 0.34 to 0.43 with an increase in electrolyte content. The reference salt concentration of the SBE liquid phase remained constant at 1000 mol/m³ throughout these measurements. For the cells manufactured by Siraj et al.¹, the SBE contained 50 wt% liquid electrolyte.

Parameter identification

To identify the remaining material parameters we utilise the experimental charge-rest-discharge profiles from the first cell described in section “Introduction”. The optimal material parameters, θ^*, are determined as the following minimiser

Here, Φ_exp,C/20(t) and Φ_exp,C/10(t) is the experimental cell potentials over time, containing N_C/20 and N_C/10 data points, respectively. The material parameters, θ, considered in the calibration are collected in Table 4.

Although the ionic conductivity and transference number in the SBE are measured parameters, it is evident that they are sensitive to measurement errors based on the spread of the parameters with respect to the weight percent of solid phase¹³. Therefore, we include κ_SBE and t⁺ in the parameter set, where variations within the measured range is allowed. Once θ^* is obtained, we assess the sensitivity of the simulated cell potential to perturbations in each material parameter. The purpose of the sensitivity analysis is to quantify how the model parameters affect the obtained, simulated cell potential over time. The sensitivity of parameter i is computed as

where h_i is the magnitude of the perturbation and \({({{\boldsymbol{e}}}_{i})}_{j}={\delta }_{ij}\).