Dimensional model

The well-known P2D model described by Doyle et al. [1] is considered herein as a baseline. As each domain is normalized with its thickness, the gradients used are:

\[\begin{split}\begin{gathered} \nabla = \begin{pmatrix} \frac{\partial}{\partial x} \\ \frac{\partial}{\partial y} \\ \frac{\partial}{\partial z} \end{pmatrix} = \begin{pmatrix} \frac{1}{L_i} \frac{\partial}{\partial \hat{x}} \\ \frac{1}{H_i} \frac{\partial}{\partial \hat{y}} \\ \frac{1}{W_i} \frac{\partial}{\partial \hat{z}} \end{pmatrix} \end{gathered}\end{split}\]

Where \(L_i\), \(H_i\) and \(W_i\) are the domain thickness, height and width. However, in the case of the P2D model, only the first row of the vector is used.

Electrochemical Model

• Mass transport in the electrolyte

  \[\begin{gathered} \varepsilon_\mathrm{e} \frac{\partial c_\mathrm{e}}{\partial t} = \nabla\cdot \left( D^\mathrm{eff}_\mathrm{e} \nabla c_\mathrm{e} \right) + \frac{1}{F}\sum_{i=0}^{n_\mathrm{mat}} a_i j_{i} \end{gathered}\]

• Charge transport in the electrolyte

  \[\begin{gathered} - \nabla\cdot \left( \kappa_\mathrm{eff} \nabla \varphi_\mathrm{s} \right) + \nabla\cdot \left( 2 \kappa_\mathrm{eff} \frac{RT}{F} (1-t_0) \left(1+\frac{\partial \ln{f_{\pm}}}{\partial \ln{c_\mathrm{e}}}\right) \nabla \ln{c_\mathrm{e}} \right) = \sum_{i=0}^{n_\mathrm{mat}} a_i j_{i} \end{gathered}\]

• Charge transport in the electrode

  The charge transport in the solid electron conductor materials is calculated as follows

  \[\begin{gathered} - \nabla\cdot \left( \sigma_\mathrm{eff} \nabla \varphi_\mathrm{s} \right) = - \sum_{i=0}^{n_\mathrm{mat}} a_i j_{i} \quad ; \quad \sigma_\mathrm{eff} \nabla \varphi_\mathrm{s} \Big|_\mathrm{tab} = I_\mathrm{app} \end{gathered}\]

  In the current colectors, a similar equation is used, but in this case, there is no exchage with the electrolyte, therefore the right hand side term equals zero.

• Mass transport in the active material (pseudodimension)

  The mass transport in the active material is calculated using Legendre polynomials.

  \[\begin{gathered} \frac{\partial c_\mathrm{s}}{\partial t} = \nabla\cdot \left( D_\mathrm{s}^\mathrm{eff} \nabla c_\mathrm{s} \right) \quad ; \quad - D_\mathrm{s}^\mathrm{eff} \nabla c_\mathrm{s} \Bigg|_{r=R_p} = \frac{j_{i}}{F} \end{gathered}\]

• Exchange current between the electrolyte and the electrode by lithium intercalation.

  \[\begin{gathered} j_i = F k_0 c_\mathrm{e}^\alpha (c_\mathrm{s}^\mathrm{surf})^{1-\alpha} (c_\mathrm{s}^\mathrm{max}-c_\mathrm{s}^\mathrm{surf})^\alpha \left(e^{\frac{-\alpha F}{RT} \eta }-e^{\frac{(1-\alpha) F}{RT} \eta} \right) \end{gathered}\]

• Overpotential

  \[\begin{gathered} \eta = \varphi_\mathrm{s} - \varphi_\mathrm{e} - U_\mathrm{eq}(c_\mathrm{s}^\mathrm{surf}) \end{gathered}\]

Thermal Model

• Energy conservation:

  The heat transfer across the cell is computed as follows

  \[\begin{gathered} \rho c_p \frac{\partial T}{\partial t} = \nabla\cdot \left( \lambda_\mathrm{eff} \nabla T \right) + q \quad ; \quad \mathbf{n} \cdot \left(-\lambda_\mathrm{eff} \nabla T \right) \Bigg|_{\Gamma}= h (T-T_\mathrm{ext}) \end{gathered}\]

• Heat generation:

  Several heat sources have been considered using formulations based on Chiew et al. [8]

  \[\begin{gathered} q = q_\mathrm{ohm}+q_\mathrm{rev}+q_\mathrm{irr} \end{gathered}\]

  • Ohmic heat source

    This corresponds to the heat generated by the transport of charge within the cell.

    \[\begin{split}\begin{align*} q_\mathrm{ohm} &= q_\mathrm{solid} + q_\mathrm{liquid}, \\ q_\mathrm{solid} &= \sigma_\mathrm{eff} \nabla \varphi_\mathrm{s} \nabla \varphi_\mathrm{s}, \\ q_\mathrm{liquid} &= \kappa_\mathrm{eff} \nabla \varphi_\mathrm{e} \nabla \varphi_\mathrm{e} - 2 \kappa_\mathrm{eff} \frac{RT}{F} (1-t_0^+) \left(1+\frac{\partial \ln{f_{\pm}}}{\partial \ln{c_\mathrm{e}}}\right) \frac{\nabla c_\mathrm{e}}{c_\mathrm{e}} \nabla \varphi_\mathrm{s}. \end{align*}\end{split}\]

  • Reversible reaction heat source

    The reversible heat caused by the reaction is proportional to the entropy change, that is approximated with the variation of Open Circuit Potential.

    \[\begin{gathered} q_\mathrm{rev} = \sum_{i=0}^{n_\mathrm{mat}} a_i j_{i} T \frac{\partial U_i(c^\mathrm{surf}_\mathrm{s})}{\partial T} \end{gathered}\]

  • Irreversible polarization heat source

    This represents the irreversible heating due to the polarization heat generated by the exchange current at the electrolyte-electrode interface.

    \[\begin{gathered} q_\mathrm{rev} = \sum_{i=0}^{n_\mathrm{mat}} a_i j_{i} \eta_{i} \end{gathered}\]

Degradation Models

• SEI Models: Solvent-diffusion

  This model is implemented inside the SEI class. The model considers that the SEI is originated by the electrochemical reaction between a EC solvent molecule, two lithium ions and two electrons at the electrode surface:

  \[\begin{gathered} \rm EC + 2 Li^+ + 2 e^- \rightarrow V_\mathrm{\scriptstyle SEI}. \end{gathered}\]

  Therefore, the reaction equation reads:

  \[\begin{gathered} j_\mathrm{\scriptscriptstyle SEI} = F k_\mathrm{\scriptscriptstyle SEI} c_\mathrm{\scriptscriptstyle EC}\big|_{r=R_\mathrm{s}} \exp{\left(\frac{-\beta F}{RT}\left(\eta + U_\mathrm{eq}-U_\mathrm{eq,\scriptscriptstyle SEI}\right)\right)}, \end{gathered}\]

  where the concentration of EC solvent at the SEI is modelled according to the transport equation:

  \[\begin{gathered} \frac{\partial c_\mathrm{\scriptscriptstyle EC}}{\partial t} = \nabla\cdot \left( D_\mathrm{\scriptscriptstyle EC} \nabla c_\mathrm{\scriptscriptstyle EC} - \frac{ \partial \delta_\mathrm{\scriptscriptstyle SEI}}{\partial t} c_\mathrm{\scriptscriptstyle EC} \right), \end{gathered}\]

  with the following boundary conditions:

  \[\begin{gathered} \left( D_\mathrm{\scriptscriptstyle EC} \nabla c_\mathrm{\scriptscriptstyle EC} - \frac{ \partial \delta_\mathrm{\scriptscriptstyle SEI}}{\partial t} c_\mathrm{\scriptscriptstyle EC} \right) \Bigg|_{r=R_\mathrm{s}} = \frac{j_\mathrm{\scriptscriptstyle SEI}}{F} \quad ; \quad c_\mathrm{\scriptscriptstyle EC} \big|_{r=R_\mathrm{s}+\delta_\mathrm{\scriptscriptstyle SEI}} = \epsilon_\mathrm{\scriptscriptstyle SEI} c_\mathrm{\scriptscriptstyle{EC},\scriptstyle{sln}}. \end{gathered}\]

  The SEI growth can be calculated from the reaction rate and SEI components properties:

  \[\begin{gathered} \frac{\partial \delta_\mathrm{\scriptscriptstyle SEI}}{\partial t} = - \frac{M_\mathrm{\scriptscriptstyle SEI}}{2 F \rho_\mathrm{\scriptscriptstyle SEI}} j_\mathrm{\scriptscriptstyle SEI}. \end{gathered}\]

  Thus, the total exchange current has two components:

  \[\begin{gathered} j_\mathrm{tot} = j_\mathrm{int} + j_\mathrm{\scriptscriptstyle SEI}. \end{gathered}\]

  And the overpotential has now an additional component corresponding to the voltage drop caused by SEI resistance:

  \[\begin{gathered} \eta = \varphi_\mathrm{s} - \varphi_\mathrm{e} - U_\mathrm{eq} - G_\mathrm{\scriptscriptstyle SEI} j_\mathrm{tot}, \end{gathered}\]

  where the resistance of the SEI is calculated as:

  \[\begin{gathered} G_\mathrm{\scriptscriptstyle SEI} = R_{0,\mathrm{\scriptscriptstyle SEI}} + \frac{\delta_\mathrm{\scriptscriptstyle SEI}}{\kappa_\mathrm{\scriptscriptstyle SEI}}. \end{gathered}\]

• LAM model

  This model is implemented inside the SEI class. The model computes the lost of active material due to particle cracking driven by stresses. Therefore, the decrease of the volume fraction of active material is computed as

  \[\begin{gathered} \sigma_\mathrm{h}=\frac{\sigma_\mathrm{r}+2\sigma_\mathrm{t} }{3}=\frac{2\Omega E}{9\left ( 1-\nu \right )}\left ( 3\int_{0}^{R}\tilde{c}r^2dr-\tilde{c} \right ). \end{gathered}\]

  And the hydrostatic stress is computed from the equilibrium of stresses of a spherical electrode particle

  \[\begin{gathered} \frac{\partial \varepsilon_\mathrm{s}}{\partial t}=-\beta \left ( \frac{\sigma_\mathrm{h}}{\sigma_\mathrm{cr}} \right )^m \qquad \sigma_\mathrm{h}>0 \end{gathered}\]