Parameterization and Validation of an Electrochemical Thermal Model of a Lithium-Ion Battery

The key challenge in developing a physico-chemical model is the model parameterization. The paper presents a strategic model parameterization procedure, parameter values, and a developed model that allows simulating electrochemical and thermal behavior of a commercial lithium-ion battery with high accuracy. Steps taken are the analysis of geometry details by opening a battery cell under argon atmosphere, building upon reference data of similar material compositions, incorporating cell balancing by a quasi-open-circuit-voltage experiment, and adapting the battery models reaction kinetics behavior by comparing experiment and simulation of an electrochemical impedance spectroscopy and hybrid pulse power characterization. The electrochemical-thermal coupled model is established based on COMSOL Multiphysics® platform (Stockholm, Sweden) and validated via experimental methods. The parameterized model was adopted to analyze the heat dissipation sources based on the internal states of the battery at different operation modes. Simulation in the field of thermal management for lithium-ion batteries highly depends on state of charge-related thermal issues of the incorporated cell composition. The electrode balancing is an essential step to be performed in order to address the internal battery states realistically. The individual contribution of the cell components heat dissipation has significant influence on the temperature distribution pattern based on the kinetic and thermodynamic properties.


Introduction
Lithium-ion batteries have evolved to be a dependable energy storage technology for mobile and stationary applications due to their high operating voltage and high energy density. To ensure enhanced performance, whilst the applications demand is met, is an increasing challenge on their durability and safety [1,2]. The operating conditions, such as environment temperature, applied current, state of charge (SoC), and mechanical stress within the battery design induce component degradation and can individually or in combination with each other influence the lifetime performance of the battery [3][4][5][6]. Hence, the classification, prediction, and minimization of these stress factors have emerged as important but challenging topics in the field of battery research [7][8][9][10][11].
Different kinds of models were developed to simulate lithium-ion batteries at varying internal states. Simple electrical-circuit models are used to address control strategies within battery management systems. Impedance models are more complex and can be used to map physical processes, which allow for separation of internal physical phenomena and analyzing their dependence on operating conditions. Physico-chemical models are based on fundamental equations describing the migration and diffusion process, as well as the intercalation reaction kinetics within the batteries. Hence, they can be used to understand the impact of changes in material properties on the system's behavior and improve it. Several papers have been published based on the initial work of Newman and Tiedemann in 1975, and others [12][13][14][15][16]. Further approaches including temperature distributions, first in 2D and later in 3D, are reported in literature [17,18].
However, the key challenge in developing a physico-chemical model is the model parameterization. If conclusions on the battery performance should be drawn based on the internal battery states, the parameter set for the material under consideration is of major importance. Several study approaches have been witnessed on this topic, which range from simple estimations and literature-based parameterization [4,19,20], parameter determination by nonlinear-regression [21,22], to dedicated electrode experiments and full cell evaluation [23][24][25]. The problem with electrode material data is that small changes in the composition and microstructure of the coating can greatly influence the electrode behavior. Electrode performance related to the reaction kinetics is especially problematic, since referenced data is either not available or deviates within several orders of magnitude in comparison to the actual material used [24]. Therefore, several scientific gaps remain to be filled regarding the most appropriate parameterization approach and the accuracy of the model predictions with respect to parameter interdependencies and the model's application [12].
In this study, a commercially available 40 Ah nickel-manganese-cobalt (NMC) lithium-ion battery cell was selected as the research target, for which no interior construction plan or detailed information on material composition is given. We contribute with this paper by proposing a strategic model parameterization procedure to replicate the heuristic electrical and thermal behavior of the battery, while unrealistic behavior, especially induced at limiting operation modes, is avoided. Four mandatory parameter evaluation groups are identified: "Cell geometry details, thermo-physical properties, cell balancing, and reaction kinetics" and accompanied methods are introduced and executed. The model is established based on the COMSOL Multiphysics ® simulation platform and validated via experimental methods. The main advantage of the proposed model is to quickly and non-invasively identify a wide range of parameter values dependent on the choice of operating conditions during the parametrization and thereafter. Finally, the parameterized model was adopted to analyze the heat dissipation sources based on the internal states of the battery.

Theory
The conversion and transport theory behind electrochemistry and heat transfer will be introduced from an application-oriented point of view. This model is implemented in the COMSOL Multiphysics ® v5.4 environment using the "lithium-ion battery" submodule to describe the electrochemical behavior and an ordinary differential equation to respect heuristic thermal behavior of the battery.
The schematic of the model implemented is shown in Figure 1.
Batteries 2019, 5, x FOR PEER REVIEW 2 of 29 dependence on operating conditions. Physico-chemical models are based on fundamental equations describing the migration and diffusion process, as well as the intercalation reaction kinetics within the batteries. Hence, they can be used to understand the impact of changes in material properties on the system's behavior and improve it. Several papers have been published based on the initial work of Newman and Tiedemann in 1975, and others [12][13][14][15][16]. Further approaches including temperature distributions, first in 2D and later in 3D, are reported in literature [17,18]. However, the key challenge in developing a physico-chemical model is the model parameterization. If conclusions on the battery performance should be drawn based on the internal battery states, the parameter set for the material under consideration is of major importance. Several study approaches have been witnessed on this topic, which range from simple estimations and literature-based parameterization [4,19,20], parameter determination by nonlinear-regression [21,22], to dedicated electrode experiments and full cell evaluation [23][24][25]. The problem with electrode material data is that small changes in the composition and microstructure of the coating can greatly influence the electrode behavior. Electrode performance related to the reaction kinetics is especially problematic, since referenced data is either not available or deviates within several orders of magnitude in comparison to the actual material used [24]. Therefore, several scientific gaps remain to be filled regarding the most appropriate parameterization approach and the accuracy of the model predictions with respect to parameter interdependencies and the model's application [12].
In this study, a commercially available 40 Ah nickel-manganese-cobalt (NMC) lithium-ion battery cell was selected as the research target, for which no interior construction plan or detailed information on material composition is given. We contribute with this paper by proposing a strategic model parameterization procedure to replicate the heuristic electrical and thermal behavior of the battery, while unrealistic behavior, especially induced at limiting operation modes, is avoided. Four mandatory parameter evaluation groups are identified: "Cell geometry details, thermo-physical properties, cell balancing, and reaction kinetics" and accompanied methods are introduced and executed. The model is established based on the COMSOL Multiphysics ® simulation platform and validated via experimental methods. The main advantage of the proposed model is to quickly and non-invasively identify a wide range of parameter values dependent on the choice of operating conditions during the parametrization and thereafter. Finally, the parameterized model was adopted to analyze the heat dissipation sources based on the internal states of the battery.

Theory
The conversion and transport theory behind electrochemistry and heat transfer will be introduced from an application-oriented point of view. This model is implemented in the COMSOL Multiphysics ® v5.4 environment using the "lithium-ion battery" submodule to describe the electrochemical behavior and an ordinary differential equation to respect heuristic thermal behavior of the battery.
The schematic of the model implemented is shown in Figure 1.

Lithium Ion Battery
The implementation rests on porous electrode theory, concentrated solution theory, and kinetics equations [12,13,16]. The partial differential equations of the electrochemical and temperature propagation are reported in similar applications [20,26].
The intercalation and de-intercalation processes of the lithium-ions in each electrode are evaluated on the surface of spherical particles (r-axis), while the transfer processes are considered unidirectional between the cell layer (x-axis). The variables solved within the electrochemical model are the solid-phase lithium-ion particle concentration (c s ), the liquid-phase lithium-ion particle concentration (c l ), the solid-phase potential (φ s ), and the liquid-phase potential (φ l ) [12].
The physical phenomena considered are current and mass conservation in the solid and liquid phase. For the electrodes, Fick's Law of diffusion for spherical particles is used to derive the solid phase lithium-ion particle concentration (c s ). The determination of the liquid phase lithium-ion particle concentration (c l ) is done using mass conservation within the electrolyte in each domain. The derivation of solid-phase potential (φ s ) in the electrodes is accounted for by Ohm's law. The calculation of the liquid-phase potential (φ l ) in the electrolyte is done by using Kirchhoff's and Ohm's laws. The calculation of the pore wall flux of the lithium-ions (j k ) into the electrodes is performed with help of the Butler-Volmer reactions kinetics equation [12,27].
A complete derivation of the model equations, including initial and continuity conditions, is beyond the scope of the paper. Derivation and discussions of the model type [28][29][30][31] and the implementation in COMSOL Multiphysics ® [20] can be found within the reference list.
The governing equations are shown in Table 1. For the sake of simplicity, whenever a parameter is modified with respect to the underlying degrees of freedom, it is denoted by (·) e f f . The phase and domain indices are denoted by l, s and k = neg, pos, sep, respectively. The parameter definition and parameterization are referred to in Section 5.5.

Reaction Kinetics
Heat dissipation sources will be quantified based on the solution variables. The total battery heat dissipation density Q tot is implemented in this study as the sum of the integrated heat dissipation sources resulting during energy transfer processes [20,32]. Heat generated by side reactions and Batteries 2019, 5, 62 4 of 29 particle mixing within the electrode domains is considered negligible [2,33]. The heat dissipation density variable is defined based on [20,31,32] and is shown in Table 2.
The specific parameter definition and parameterization is referred to in Section 5.5.

Heat Transfer
The following ordinary differential equation accounts for the heuristic temperature calculation and is defined with respect to [34]: where T is temperature, T ext is the external temperature, ρ e f f is the battery density, C e f f p is the battery heat capacity, V batt is the battery volume, Q tot is the total heat dissipation density, V stack is the electrode stack volume, A ht is the battery surface with cooling by convection, h air is the convective heat transfer coefficient, A e f f rad is the surface area-dependent effective emissivity, and σ is the Stefan-Boltzmann constant.

Coupling of the Physics
The governing equations shown in Table 1 represent the electrochemical behavioral model of the battery with respect to the operations modes: I, T, and SoC. The total heat dissipation density Q tot of the battery with respect to the aforementioned setting is defined in Table 2. In Equation (1) the thermal behavior of the battery is calculated in dependence of the heat dissipated, while simultaneously exchanging the temperature T with the electrochemical variables in the particle, electrode, and electrolyte domains.

Experimental Techniques
A model parameterization includes multiple steps: "Opening a battery cell, sample extraction and preparation, measurement execution, and evaluation [23][24][25]35,36]." Special focus was taken on an accurate pre-set of the initial SoC of the battery cell before each measurement objective. The initial SoC of 100% was approached by a 1 C constant current-constant voltage charge to 4.2 V with a cut-off current of 1/20 C. The initial SoC of 0% was approached by stepwise constant current discharge phases with 1 C, 1/2 C, 1/4 C, 1/8 C, 1/20 C, and 1/40 C until the cut-off voltage of 2.8 V was reached. Between each measurement, there was at least two hours' wait for a full thermal and electrochemical homogenization at the initial state of the test environment.
In this section, measurement techniques used for the battery under investigation are described.

Battery Description
The commercially available prismatic format (PHEV2) battery used for the model parameterization is labelled to be an NMC/Graphite type with a nominal theoretical capacity of 40 Ah and 144 Wh, measured at 1/3 C constant current rate. A voltage range between 2.8 V to 4.2 V is recommended with a nominal voltage of 3.6 V. Charging within a temperature range of 0 to 45 • C up to a maximal discharge current rate of 5 C, as well as discharging with 10 C within −10 • C to 65 • C, is approved.

Pre-Mortem Techniques
The cycling for the electrical measurements was performed by a Maccor series 4000 Potentiostat, multiple cell tester with 0 to 8 V (0.02% accuracy) for voltage and 0 to 150 A current (0.05% accuracy) range each. The climate was controlled by a Vötsch climate chamber, VC 7018 with a temperature range of −75 to 180 • C, and temperature was measured by an Agilent 34972A with PT100 sensors, −40 to 85 • C (±0.1 K accuracy) temperature range. The chamber has active air circulation to ensure accurate conditions dependent on each measurement style. The electrical impedance spectroscopy (EIS) was performed by a Solartron Analytical Modulab, XM ECS System, Potentiostat.

Quasi-Open-Circuit Voltage (qOCV)-Low Current Method
The battery was measured during a separated charge and discharge phase with a 1/40 C rate and a constant temperature environment of 25 • C.

Open-Circuit Voltage (OCV)-Rest Method
The battery was measured during a discharge and a charge phase. Each stepwise discharge of 10% nominal capacity was followed by a break until a voltage change of smaller than 2 mV occurred or 1 h duration was reached. The initial state during the charge phase was the last state of charge during the discharge phase, after relaxation, when a cut-off voltage of 2.8 V was reached.

Electrochemical Impedance Spectroscopy (EIS)
The potentiostatic measurement was performed at the battery at 50% SoC, roughly at 22.5 • C, with 10 mV amplitude and a maximal current of 2 A within a frequency range of 10 mHz to 1 kHz.

Internal Resistance Measurement by Hybrid Pulse Power Characterization (HPPC)
The resistance of the battery within a 3 C discharge pulse for 30 s, measured after 2 s, 10 s, and 18 s, and within a 2 C charge pulse, measured after 10 s, was evaluated at nine different SoCs (10% to 90%), roughly at a temperature of 22.5 • C. After the first and second pulse, a relaxation time of 40 s was performed. Each stepwise discharge of 10% nominal capacity was followed by a break until a voltage change of smaller than 2 mV occurred, or 0.5 h duration was reached.

Electrical/Thermal Validation
Three batteries were measured during separate charge and discharge phases at 0.5 C, 1 C, and 2 C constant current profiles with an initial temperature set to 25 • C. While the battery was charged/discharged, the chambers climate control was turned off to avoid undesired fan cooling. The exterior environment temperature of 22.5 • C and imposed reflective radiation by the chambers interior metal surfaces are assumed to have a negligible effect on the measurement results.

Cell Opening and Sample Preparation
The battery was opened at a SoC of 0% and an ambient temperature of 22 • C. Due to the possibility of undesired material modification by contact to air, a full disassembly of the battery cell is recommended to be done under argon atmosphere within a glovebox [37]. However, to be able to separate the electrode stack and terminal from the metal casing of the prismatic battery cell, a cut at the top edge of the battery cell was performed first with a Dremel ® tool (Racine, WI, USA) under a fume hood [38]. During this process, the batteries' temperature and voltage were carefully measured. Fast changes would have indicated the activation of electrode stack internal chemical heating processes due to the reaction with air or due to the heating by the cutting process. A voltage drop above 0.1 V or a high temperature increase above 5 K difference did not happen. Thereafter, the lower part of the metal casing could be detached from the terminal and electrode stack carefully with no excessive force.
The final disassembly step of the electrode stack and terminal was done within a glovebox under an inert argon atmosphere on non-conducting ground, atmospheric pressure, and with water and oxygen below 5 ppm. Due to the electrode stack design, a clean cut with a sharp ceramic knife between terminal and stack was executed. This step has to be done with care, however, to avoid cross-contamination or short circuits by metal dust and swarfs between the electrode types. The double-sided electrodes and the separator are detached from each other and stored and sealed in separate containers. In these conditions, the electrode materials quality is not stable, since the material's fraction of active lithium might still react with small amounts of oxygen or water and electrolyte might evaporate even when the bags are sealed airtight. Therefore, further measurements should be performed within days to avoid the effect of uncertain sample dissipation.
For measurement methods outside the glovebox or when material is intended to be stored for a longer duration, however, a washing process with solvent is recommended to reduce the impact of sample dissipation. The procedure chosen within this study consists of two bathing steps per material sheet with 10 mL pure dimethyl carbonate (DMC) and 1 min duration each, always performed in exactly the same way. The procedure chosen within this study consist of two bathing steps per material sheet with 10 mL pure DMC and 1 min duration each, always performed in exactly the same way [37].

Geometrical Data Measurement
The total thickness of the electrode sheets and the separator were measured with a micrometer-screw at five different points of each material and an average was taken.
X-ray micro-computed (µ-CT, SkyScan 1172, Bruker, Belgium) and nano-computed tomography (n-CT, SkyScan 2011, Bruker, Belgium) was used to access the electrodes layer and coating configuration. All tests were done at room temperature with a set voltage of 80 kV for the µ-CT and 40 kV for the n-CT. The recorded data sets were reconstructed and selected volume of interest(s) (VOI) and region of interest(s) (ROI) were quantitatively analyzed using associated software products (SkyScan, Bruker, Belgium): • NRecon: The obtained datasets of recorded images are reconstructed to a 3D cross-Section image stack. The ring artefact and beam hardening values are kept constant for all the samples. • Data Viewer: The reconstructed images are adjusted parallel to their respective viewing plane, coronal and sagittal, and the VOI is defined. • CTAn: For the 3D morphometric analysis, the shape and sheet number of ROI are defined.
For the determination of the current collector thickness, the electrode sheets were measured by the µ-CT at a resolution of 1.17 µm/pixel with a rotation step of 0.75 • per projection, each projection being an average of four image frames. The 3D reconstructed data set was viewed at coronal and sagittal planes, the thickness was determined at a minimum of ten different points on each plane, and an average was taken. For the analysis of particle size and porosity percentage, the negative and positive electrode coatings were detached from the current collector foils. The negative electrode coating was measured by the µ-CT with a resolution of 1 µm/pixel and a rotation step of 0.18 • per projection, with each projection being an average of six frames, while the positive electrode coating was measured by the n-CT with a resolution of 0.2 µm/pixel and a rotation step of 0.2 • per projection, with each projection being an average of four frames. The morphometric analysis was done on the defined VOI and square shaped ROI of 500 µm × 500 µm for the negative electrode coating and ROI of 100 µm × 100 µm for the positive electrode coating.
Scanning electron microscopy (SEM) was utilized to visualize and confirm the particle sizes and homogenization at the electrode's surfaces. An InLens electron detector and an accelerator voltage of 20 kV were used for both electrode measurements.

Thermal Parameter Measurement
The cell components double coated graphite anode on copper foil, double coated NMC cathode on aluminum foil, and separator of the electrode stack were characterized using differential scanning calorimetry (DSC 204 F1, Netzsch, Germany) to determine their heat capacities C p .
The DSC device setup and measurement fundamentals for the evaluation of C p are described in [39]. The anode and cathode coating were scrapped from the copper and aluminum foil using a scalpel. The scrapped samples were loaded and sealed in hermetic aluminum pans. The measurements parameters were defined on the DSC 204 F1 phoenix software (Netzsch, Germany). The measurements were performed from a temperature of 25 • C to 60 • C with a heating rate of 10 K/min. The baseline correction measurement was performed with an empty pan and the calibration measurement with a sapphire standard of 12.73 mg, a thickness of 0.25 mm, and a diameter of 4 mm. The obtained data was loaded onto the proteus analysis software (Netzsch, Germany) for C p evaluation. A similar application of the methods for the characterization of lithium-ion battery components is shown in [1].

Simulation Techniques
The focus in this section is to assess the types of simulation methods that are needed to analyze the design key parameters that define the correlation of the electrode and cell voltage. The determination of the contributing factors on the cell voltage due to the current flow is separated in the following standardized way, as shown in [40]: Therein V cell represents the cell voltage, V OCV is the open circuit voltage without an applied load, µ ohmic is a pure ohmic voltage drop, µ f is a voltage drop due to passivation film or layers on each electrode respectively, µ pol are the polarization voltage drops due to interfacial charge-transfer reactions at each electrode, and µ con is a concentration-related effect like particle diffusion, migration, or convection [40]. Each of the different phenomena is reported to occur on different time scales within the battery's application and deviate dependent on the considered chemistry [40].
At the start, the behavior of the electrode vs. cell balancing ratio is analyzed in Section 5.3, when no current applied. Hence, V cell in Equation (2) is solely defined by the electrode's contributions: Afterward, reaction kinetics methods are considered that are able to separate contributing factors on the battery voltage in Section 4.2, when a current is applied. In both sections, a set of subsequent equations is provided that point out modifications for the application in the model.

Cell and Electrode Examination
The aim within this examination is to find a suitable and realistic representative function of the negative electrode equilibrium potential U neg under the assumption of a known positive electrode equilibrium potential U re f pos and a measured quasi-equilibrium cell voltage V qOCV . A similar approach is reported in previous studies [3,25]. The stoichiometry of both electrodes at 100% and 0% SoC (four parameters: SoC min pos , SoC max pos , SoC min neg , and SoC max neg ) are initially chosen to be the same as in [25]. The following steps one to three are performed to examine the electrode potentials based on the investigated cell SoC: 1.
The qOCV potential for the cell is measured during a charge and discharge phase to investigate the actual capacity Q cell . The OCV potential of the cell is approximated by the average of the qOCV profiles at normalized time, since the overvoltage error is canceled with respect to the direction of the energy transfer: The overvoltage error cancels out with respect to the direction of each energy transfer.

2.
The positive electrodes potential in dependence of the electrode SoC is chosen to be taken from literature [25]: The positive electrode potential profile is known to show no extreme characteristic slopes in comparison to the negative electrode potential within the mean SoC range.

3.
By applying Equations (4) and (5) to Equation (3), the negative electrodes potential is calculated as follows: By shifting and rescaling parts of the electrode potential curves to fit to the difference of the cell potential the parameters SoC min pos , SoC max pos , SoC min neg and SoC max neg are derived that uniquely define the electrode SoC windows. A shift of the function U neg is happening, when the threshold values for U re f pos are transformed by the same value, while rescaling of the function U neg is happening, when the difference ∆SoC k = SoC max k − SoC min k is altered. The characteristic potential slope S2, which is localized at a theoretical electrode SoC of 50%, is considered as fit criterion to uniquely define the final thresholds values [41]. A discussion of this step and the results is referred to in Section 5.3. The data processing and implementation of the potential profiles is referred to in Appendix A.
Based on the SoC thresholds of the electrodes and with respect to the maximal lithium particle concentrations of each electrode, the following equations are established with respect to ideas of [28,36] and can be used to define the initial lithium particle concentration of the electrodes at the initial state of charge SoC 0 : The upcoming set of equations are derived and used to reproduce the design specification of the full cell based on the electrode's properties (k = neg, pos) and are used similarly by [23][24][25]: • The inactive fraction of the electrodes SoC windows is calculated by: Batteries 2019, 5, 62 9 of 29 • The actual electrode capacities Q neg and Q pos are calculated based on the fact that the usable capacity Q cell represents the active material loading, while the counterparts in each electrode remain inactive: • The theoretical electrode capacities are calculated as follows: The electrodes active material fraction is received by calculating the ratio of the actual electrode capacity of the theoretical electrode capacity: • The specific surface area is calculated as follows: • Similarly, the actual surface area is calculated as follows: The final results are reported in Table 6 in Section 5.5.

Reaction Kinetics Examination
Several studies reported on EIS experiments of lithium-ion batteries and half-cells to evaluate cell performance, as shown in Equation (2), and employed the method as reliable to separate the dynamic effects based on the measured impedance response [25,40,[42][43][44][45][46]. A small sinusoidal current signal I is imposed on the battery at a given frequency ω and induces a sinusoidal voltage V response. Using complex numbers to represent the sinusoidal current and voltage, the complex impedance Z of the battery at a given frequency can be calculated as follows: In fact, two characteristic points and corresponding resistance values will be assessed that roughly describe the dynamic behavior of the battery [46]: The impedance spectrum crosses the x-axis at the frequency, when the inductive behavior is compensated by capacitive behavior of battery components. Hence, the battery impedance at this point is nearly the pure ohmic resistance R 0 of the battery. • R Ct at min −Im Z : The negative imaginary part of the impedance spectrum shows a local minimum at the end of a large semicircle at low frequencies. The semicircle is usually associated with interfacial charge transfer reaction combined with double layer capacitance.
The entire impedance spectrum of the battery is evaluated based on the frequencies in the given frequency range 1000 Hz to 0.01 Hz, when a sinusoidal voltage response with 10 mV is produced. The COMSOL Multiphysics ® physico-chemical modeling approach [45,47] will be introduced based on the model study one needs to perform to evaluate equation (15). The model study is a stationary frequency analysis and requires modifications of the governing equations: The time-dependent model variables, here represented by n, are modified by a sinusoidal perturbation and solved in the complex-valued domain for predefined frequencies ω :

2.
Where n re f is the initial time-dependent value at time zero and n is the complex-valued extension of n. The particle concentration time derivative mentioned in Table 1 is exchanged by a frequency mode as follows: where j k incorporates the double layer capacitance C dl,k at each electrode: Correspondingly, j k instead of j e f f k in the governing equations in Table 1.

3.
Thus, the cell impedance Z is calculated dependent on the perturbed boundary flux variables φ s and I s and incorporated conditions defined at the positive current collector: Therein R 0 and L 0 are a resistance and inductance component, i s is equal to the applied current density in the solid domain, φ s is the potential within the solid domain, φ s,re f and i s,re f are the initial time-dependent values at time zero, n is the normal boundary flux vector, ∆φ s is the amplitude of the potential perturbation and ∆i s is the amplitude of the current perturbation induced by ∆φ s .
The model results are initiated based on reference values [25] and refined with respect to the EIS experiment at cell SoC 50% and temperature at 25 • C. The final fitting results are shown in Table 6 in Section 5.5.
For the characterization of low frequency battery dynamics, a HPPC experiment and simulation is considered, as suggested in literature [48,49]. The pulse methodology is not able to separate dynamic contributing factors to the internal resistance value, but it is able to show its evolution and is considered to be the more standardized way in present battery diagnostics and scientific evaluations [40]. To analyze the internal resistance progression at multiple time constants and different SoC levels, the internal resistance is evaluated as the battery's voltage response of a constant current pulse during a discharge and charge phase: where t is time, V(t) is the voltage measured at time t, V(0) is the voltage measured the pulse initiation, and I is the applied current. The following set of equations were chosen and partially refined to reproduce the reaction kinetics of the cell with respect to the governing equations and underlying domains k = neg, pos, sep defined in Table 1: • The temperature dependence of model variables is considered by applying the Arrhenius relation as considered in [25]: where n is the variable, n re f is the variables value at the reference temperature T re f , E act n is the activation energy, R is the universal gas constant, and T is the temperature variable.
• Each electrode's exchange current density is defined temperature-dependent with respect to [25]: • Each electrode's diffusion coefficient is dependent on the electrodes SoC and temperature. Therefore, the following equations are defined with respect to [25]: where • The electrolyte diffusion coefficient is defined as follows, as taken from [50]: where D re f l,k c l,k = 8.794 · 10 −11 c 2 l,k − 3.972 · 10 −10 c l,k + 4.862 · 10 −10 (29) • The Li + transference number t + 0 is defined as follows, as taken from [50]: • The effective thermodynamic activity f e f f cl coefficient is defined as follows, as taken from [50,51]: • The electrolyte conductivity is defined as follows, as taken from [3,51]: where σ re f l,k c l,k = 0.1297 · c 2 l,k − 2.51 · c 1.5 l,k − 0.3329 · c l,k • The effective electrode conductivity is defined as follows, similar to [3,25]: • The effective diffusional electrolyte conductivity is defined as follows, similar to [3]: The final results are reported in Table 6 in Section 5.5.

Thermal Behaviour Examination
Effective parameter sets for the density ρ and heat capacity C p are considered with respect to the geometry of the cell and are calculated based on volume, averaging similar to [1]. Therefore, the geometry was reconstructed based on volumetric measurements of each material component of the battery cell: where is the volume fraction, p i is the density, and C p,i is the heat capacity of the i-th material component of the battery cell.
Convective cooling by air is considered material-independent and the surface area A ht is simply the sum of each material surfaces exposed to air. Cooling by radiation, however, is dependent on the materials surface emissivity. Therefore, the effective emissivity is calculated as follows: where A e f f rad is the surface-area dependent emissivity, A rad,i is a surface with applied cooling by radiation, and i is the i-th material emissivity.
The final results are reported in Table 6 in Section 5.5.

Results and Discussion
In this section, the experimental characterization results and evaluated key parameters are gathered and incorporated into the simulation model of the commercially available 40 Ah battery cell. It is investigated how the electrodes balancing and reaction kinetics related phenomena can be derived from full cell measurements and used to define the electrodes characteristics. The simulations are compared with validation experiments to evaluate the electrical and thermal model behavior under the influence of simulated heat dissipation.

Electrode Geometry Evaluation
The geometry parameters presented in this section are analyzed by post-mortem analysis of the commercial battery cell. The experimental procedures undertaken for this step are explained in Section 3.3.
Based on the material analysis by SEM and data sheet description of the manufacturer, the negative electrode is considered to be graphite, while the material of the positive electrode is supposed to be of the NMC class. The negative electrode consists of flake-shaped particles and the positive electrode consists of agglomerates of smaller particle types as shown in SEM pictures in Figure 2: under the influence of simulated heat dissipation.

Electrode Geometry Evaluation
The geometry parameters presented in this section are analyzed by post-mortem analysis of the commercial battery cell. The experimental procedures undertaken for this step are explained in Section 3.3.
Based on the material analysis by SEM and data sheet description of the manufacturer, the negative electrode is considered to be graphite, while the material of the positive electrode is supposed to be of the NMC class. The negative electrode consists of flake-shaped particles and the positive electrode consists of agglomerates of smaller particle types as shown in SEM pictures in Figure 2: The electrode stack of the battery cell is constructed based on a z-folding of the separator with each electrode contained within the folds in alternate order. All electrode sheets are double-coated, even the exterior negative electrode sheets, although no counter electrode is present. Only the interior electrode area is considered to be electrochemically active, therefore the electrode surface area with counter electrode is chosen to be used for the parameterization. Electrodes material properties considered for the parameterization are the thicknesses, particles sizes, and the porosities. The The electrode stack of the battery cell is constructed based on a z-folding of the separator with each electrode contained within the folds in alternate order. All electrode sheets are double-coated, even the exterior negative electrode sheets, although no counter electrode is present. Only the interior electrode area is considered to be electrochemically active, therefore the electrode surface area with counter electrode is chosen to be used for the parameterization. Electrodes material properties considered for the parameterization are the thicknesses, particles sizes, and the porosities. The coating thickness is recalculated based on the current collector thickness. The evaluated geometry properties are listed in Table 3.

Battery Thermal Parameter Examination
Within each component of the battery, thermo-physical properties are different; hence, each component has to be taken into account. The thermal parameters for the battery are listed in Table 4. The obtained values for the single components are in good agreement with data reported for a NMC/Graphite battery [2]. Temperature dependence of the heat capacity is not considered for the model parameterization, since the analysis of the cell sheet components have shown almost constant behavior within the range from 25 • C to 60 • C. Therefore, the mean value is taken instead.
The radiation related surface properties of the battery materials, shown in Table 5, are taken from reference [53].

Cell Balancing Evaluation
The purpose of the cell balancing evaluation is to access electrode SoC and cell SoC vice versa for SoC-dependent model application and evaluation.
The capacity in the qOCV-experiment during the discharge-phase is Q cell, dch = 37.3 Ah and during the charge-phase is Q cell, ch = 36.8 Ah. The measured capacity difference of 0.5 Ah might be attributed to one of the following reasons: • Due to the difficulty to set up the initial cell SoC at 0%, as described in Section 3.2, an off-set within the experiment is realistic.

•
An overhang effect at the negative electrode during a prior charge-phase lead to an increased discharge capacity after relaxation [54].
Since this effect is not reproducible with this model type, the cell capacity is defined based on the discharge phase result.
The existence of characteristic potential slopes in the interim SoC range of the full cell are reported to occur within the Graphite OCV profile rather than the positive electrode of different Graphite/Li-Metal-Oxides [23,25,41]. A typical graphite electrodes OCV profile shows several potential slopes dependent of the SoC stages. This is caused by the structure of graphite, where the lithium intercalation between the graphene layers occurs in an ordered pattern [41,55]. The potential profiles for the cell balancing evaluation, a simulation variant, and an experimentally determined OCV profile by rest method are shown for comparison with the suggested lithiation stages within the calculated graphite electrode OCV profile in Figure 3: Graphite/Li-Metal-Oxides [23,25,41]. A typical graphite electrodes OCV profile shows several potential slopes dependent of the SoC stages. This is caused by the structure of graphite, where the lithium intercalation between the graphene layers occurs in an ordered pattern [41,55]. The potential profiles for the cell balancing evaluation, a simulation variant, and an experimentally determined OCV profile by rest method are shown for comparison with the suggested lithiation stages within the calculated graphite electrode OCV profile in Figure 3: The lithium distribution is numerated by a stage number, which indicates the number of graphene layers between the layers of intercalated lithium. Areas of coexisting intercalation phases are seen as plateaus and single-phase areas as change in potential. The intercalation stages are named The lithium distribution is numerated by a stage number, which indicates the number of graphene layers between the layers of intercalated lithium. Areas of coexisting intercalation phases are seen as plateaus and single-phase areas as change in potential. The intercalation stages are named with an increase in lithium-ion concentration as follows: dilute stage-1 (S1'), stage-4 (S4), stage-3 (S3), liquid-like-stage 2 (S2L), ordered-stage-2 (S2), and ordered stage-1 (S1) [41,56]. S2 is considered as characteristic potential change for the electrode SoC threshold evaluation, since it is located at the theoretical electrode SoC of 50% [25,56,57].

Reaction Kinetics Evaluation
By fitting the experimental EIS results at cell SoC 50% and temperature 25 • C, the initial electrode exchange current densities taken from [25] are updated.
The Nyquist plot obtained from the experiment and the cell impedance model are shown together, with each contributing impedance part in the frequency domain in Figure 4.
Three different processes with the identified frequency ranges are annotated in Figure 4a: region I: 100 Hz up to 1000 Hz, region II: 1 Hz up to 100 Hz, and region III: 0.01 Hz up to 1 Hz. Each region is reported to coincide with a separate underlying physical phenomenon: experimental test-setup (region I), charge-transfer processes (Region II), and diffusion processes (region III) [58]. In region I, an inductive component is existent. A deviation between the simulation and the experiment results is visible. The combination of an inductive component L 0 connected in series with a resistance R 0 is not able to reproduce the slope of the experimental pattern. The effect however can be interpreted as measurement artifacts that originate from connections at high frequencies and might be attributed to the wiring of the experimental test setup and not necessarily to the battery cell [59].
(region I), charge-transfer processes (Region II), and diffusion processes (region III) [58]. In region I, an inductive component is existent. A deviation between the simulation and the experiment results is visible. The combination of an inductive component 0 connected in series with a resistance 0 is not able to reproduce the slope of the experimental pattern. The effect however can be interpreted as measurement artifacts that originate from connections at high frequencies and might be attributed to the wiring of the experimental test setup and not necessarily to the battery cell [59]. The semicircle of the charge transfer region II is fully resolved and can be fitted to derive values for the reaction rate constants and their double layer capacitances. , [43]. The diffusion process within region III cannot be fully analyzed due to the limitations in the frequency range of the experimental setup below 0.01 Hz. Hence, a diffusion-induced contribution on the cells internal resistance will be evaluated based on the pulse methodology. The absolute mean errors for the fitting procedure are shown in Figure 4b,c. The refined values are denoted in Table 6 in Section 5.5.
The evaluation of the HPPC test of the internal resistance is performed based on Equation (22). In Figure 5, mean experimental results and deviations of four different battery cells of the same kind are compared against the simulation results of the battery that is considered for the model's parameterization. The semicircle of the charge transfer region II is fully resolved and can be fitted to derive values for the reaction rate constants K k and their double layer capacitances. C dl,k [43].
The diffusion process within region III cannot be fully analyzed due to the limitations in the frequency range of the experimental setup below 0.01 Hz. Hence, a diffusion-induced contribution on the cells internal resistance will be evaluated based on the pulse methodology. The absolute mean errors for the fitting procedure are shown in Figure 4b,c. The refined values are denoted in Table 6 in Section 5.5.
The evaluation of the HPPC test of the internal resistance is performed based on Equation (22). In Figure 5, mean experimental results and deviations of four different battery cells of the same kind are compared against the simulation results of the battery that is considered for the model's parameterization. The general trend visible is an increase of the internal resistance ( ) with an increasing time constant = 2 s to = 18 s at all SoC levels during the discharge and charge pulse in the simulation and in the experiment. By comparing the magnitude of ( ) within the EIS experiment in Figure 4 at the SoC level 0.5 and temperature 25 °C, it can be concluded that the internal resistance value (2 ) already comprises the ohmic resistance, the complete charge-transfer-resistance, and parts of The general trend visible is an increase of the internal resistance R(t) with an increasing time constant t = 2 s to t = 18 s at all SoC levels during the discharge and charge pulse in the simulation and in the experiment. By comparing the magnitude of Re(Z) within the EIS experiment in Figure 4 at the SoC level 0.5 and temperature 25 • C, it can be concluded that the internal resistance value R(2s) already comprises the ohmic resistance, the complete charge-transfer-resistance, and parts of the diffusion resistance. Hence, the resistance increases between R(2s) and R(18s) within the considered temperature regime can mainly be attributed to the effect of electrolyte salt diffusion and electrode diffusion [40,59,60]. Table 6. Summary of all the material and cell properties resembling the analyzed battery cell.  Table 3 Electrode Plate Area A (m 2 ) 2.1024 Table 3 Particle Radius r (µm) 9.89 1.72 Table 3 Actual  At the SoC 0.1, the experimental resistance values are increased at all time constants within comparison to the SoC window 0.2 < SoC < 0.9. Similarly, the mean error within the SoC window 0.2 < SoC < 0.9 is increased in comparison to the SoC 0.1.The simulation results are located within the error margin of the experimental results in Figure 5a and remain close to or within the error margin in Figure 5b,c. Two distinct peaks are visible within the simulation results that are located at each time constant within the SoC range 0.1 to 0.4 and 0.6 to 0.8, called peak one and peak two, respectively. Their location coincides with the potential slopes due to the structural changes of S4 and S2 within the negatives electrode, as seen in Figure 3. A shift of the peak might be attributed to a SoC change during the applied pulse profile.
The peaks, although heavily shaved, are also visible within the experimental results, around the peaks SoC windows. There are two major reasons that would indicate a peak shaving behavior within the experimental result:

•
Each nominal SoC level of the batteries could comprise a deviation in their real capacity that cannot be captured by the experiment setup.

•
Each battery contains ninety cells gathered in a parallelized electrode stack. Non-uniformity of the cell SoCs within the stack would indicate a variation in duration and time until the slopes S2 or S4 are passed. Since the batteries' internal resistance comprises the contributions of each of the cell's contribution, a slight shift within the peak's potential would shave the peak magnitude and move its location.
The model is not capable of reproducing the peak shaving behavior, since only one cell sheet is considered; however, the model accurately predicts the increase of the internal resistance dependent on the Cell SoC within comparable magnitude.

Model Parameterization
In the following, the model's setup is demonstrated.

Cell Voltage Validation
The model validation considers both electrochemical and thermal properties. The experimental test procedures considered electrochemical and thermal constant current tests, as described in Section 3.2. For electrical validation, both the voltage curves during a 0.5 C, 1 C, and 2 C discharge and charge phases were taken into account, as shown in Figure 6a-f. The simulated curves match well with the experimental curve and show a relative error within 1% deviation.
More pronounced deviations are visible at the end of the discharge phase and at the beginning of the charge phase. These regions are correlated with a low SoC of the negative electrode and a high SoC of the positive electrode.

Cell Voltage Validation
The model validation considers both electrochemical and thermal properties. The experimental test procedures considered electrochemical and thermal constant current tests, as described in Section 3.2. For electrical validation, both the voltage curves during a 0.5 C, 1 C, and 2 C discharge and charge phases were taken into account, as shown in Figure 6a-f. The simulated curves match well with the experimental curve and show a relative error within 1% deviation.
More pronounced deviations are visible at the end of the discharge phase and at the beginning of the charge phase. These regions are correlated with a low SoC of the negative electrode and a high SoC of the positive electrode.

Temperature Validation
For the thermal validation, both the temperature curves during a 0.5 C, 1 C, and 2 C discharge and charge phases were taken into account. The experiments were conducted as mentioned in Section 3.2. The battery temperature was measured by fourteen thermocouples, positioned at different locations on the batteries surface. The average battery temperature and maximal temperature variability of the experiment were compared with the simulation as shown in Figure 7a-f:

Temperature Validation
For the thermal validation, both the temperature curves during a 0.5 C, 1 C, and 2 C discharge and charge phases were taken into account. The experiments were conducted as mentioned in Section 3.2. The battery temperature was measured by fourteen thermocouples, positioned at different locations on the batteries surface. The average battery temperature and maximal temperature variability of the experiment were compared with the simulation as shown in Figure 7a The range within the variability of the experimental values and the simulation are basically the same at each constant current rate. However, the relative error increases from 0.5 C to 2 C and the model tends to overestimate the experimental results at elevated current rates. The deviations between the experimental and the simulated temperature values might be attributed to one of the three upcoming reasons:  No spatial resolution of the battery cell is considered. Hence, no thermal gradient based on the The range within the variability of the experimental values and the simulation are basically the same at each constant current rate. However, the relative error increases from 0.5 C to 2 C and the model tends to overestimate the experimental results at elevated current rates. The deviations between the experimental and the simulated temperature values might be attributed to one of the three upcoming reasons: • No spatial resolution of the battery cell is considered. Hence, no thermal gradient based on the interior cell component connection is created.

•
The active material within the electrode stack is considered to be homogenous and stable according to the electrochemical model definition, but side reactions or the effect of heat of mixing could influence the experimental results, although considered negligible by others as well [2,32,33].

•
The environment conditions of the battery could affect the experimental results. The climate chamber was turned off, when the experiment was initialized at a constant temperature to prevent weak forced cooling, but fluctuations could still arise by the buoyancy effect within air [64]. The wires and fastener connected to the battery cell could influence the experimental thermal behavior of the battery cell by inducing a thermal outflow [2].

The Effect of Heat Dissipation on Thermal Behaviour
The heat dissipation sources resulting from effects from each layer and the current collectors within the electrochemical model are categorized into two categories: reversible heat Q rev , originated from entropy change of the electrode materials, and irreversible heat, divided into polarization heat Q irr,pol and ohmic heat. Q irr,ohmic The simulated heat dissipation densities and temperature curves during 0.5 C, 1 C, and 2 C discharge and charge phases are simulated and shown in Figure 8a The total heat dissipation is visible as sum of the several contributions. The hatched area depicts compensated dissipation by the mix of endothermic and exothermic heat sources.
The irreversible ohmic heat , ℎ is almost constant at each C-rate, independent of the charge or discharge direction, and increases from 0.5 C to 2 C practically 15 times.
The irreversible polarization heat , has two characteristic peaks within the interims SoC and increases at SoC below 10% significantly at all current rates. The peaks are shifted from their original location at 0.5 C dependent on the charge and discharge direction. They grow and widen with increased current rate up to 15 times their original magnitude. The existence and origin of the peaks is discussed in Section 5.4 based on the resistance evaluation.
The reversible heat dissipation comprises the electrodes entropy profiles [32,63]. The SoCdependency in this study, however, solely originated from the negative electrode, since the NMCbased electrodes contribution is assumed to be constant [63]. The reversible heat source acts as a coolant at low current SoC during the charge phase and at the interims SoC range during the discharge phase. The first case is reported to correlate to a transition between ordered and disordered states, while the second case is correlated to the transition between the stages S2L and S2 within the lithiation of the negative electrode, as mentioned in Section 5.3.
With increasing current rate, irreversible heat dissipation is dominant. The result is represented by the fact that the characteristic thermal pattern, which is originated from reversible heat dissipation, almost vanished at the 2 C rate; both visible in the simulation and in the experiments.

Conclusions
The strategic model parameterization procedure and the developed model established in this study allow simulating electrochemical and thermal behavior of a commercial lithium-ion battery cell with high accuracy.
This work shows that commercially available battery cell models can be parametrized by a sophisticated mixture of material characterization, literature review, and model refinement. This has been confirmed by multiple experimental characterization approaches focused on the battery cell at hand. The hierarchical model parameterization procedure starting from geometry details, incorporating cell balancing and adapting the model to reaction kinetics behavior might be applied by other model developers as well, when data of the considered electrodes are not available. The total heat dissipation Q tot is visible as sum of the several contributions. The hatched area depicts compensated dissipation by the mix of endothermic and exothermic heat sources.
The irreversible ohmic heat Q irr,ohmic is almost constant at each C-rate, independent of the charge or discharge direction, and increases from 0.5 C to 2 C practically 15 times.
The irreversible polarization heat Q irr,pol has two characteristic peaks within the interims SoC and increases at SoC below 10% significantly at all current rates. The peaks are shifted from their original location at 0.5 C dependent on the charge and discharge direction. They grow and widen with increased current rate up to 15 times their original magnitude. The existence and origin of the peaks is discussed in Section 5.4 based on the resistance evaluation.
The reversible heat dissipation Q rev comprises the electrodes entropy profiles [32,63]. The SoC-dependency in this study, however, solely originated from the negative electrode, since the NMC-based electrodes contribution is assumed to be constant [63]. The reversible heat source acts as a coolant at low current SoC during the charge phase and at the interims SoC range during the discharge phase. The first case is reported to correlate to a transition between ordered and disordered states, while the second case is correlated to the transition between the stages S2L and S2 within the lithiation of the negative electrode, as mentioned in Section 5.3.
With increasing current rate, irreversible heat dissipation is dominant. The result is represented by the fact that the characteristic thermal pattern, which is originated from reversible heat dissipation, almost vanished at the 2 C rate; both visible in the simulation and in the experiments.

Conclusions
The strategic model parameterization procedure and the developed model established in this study allow simulating electrochemical and thermal behavior of a commercial lithium-ion battery cell with high accuracy.
This work shows that commercially available battery cell models can be parametrized by a sophisticated mixture of material characterization, literature review, and model refinement. This has been confirmed by multiple experimental characterization approaches focused on the battery cell at hand. The hierarchical model parameterization procedure starting from geometry details, incorporating cell balancing and adapting the model to reaction kinetics behavior might be applied by other model developers as well, when data of the considered electrodes are not available.
Under consideration of the heuristic battery thermal properties, this model has proven to be a trustworthy tool to investigate the heat dissipation during diverse battery operation modes. It could be confirmed that the thermal behavior of the battery cell at lower current rates is dominated by reversible heat dissipation. With increasing current rates, a shift towards the dominance of irreversible heat dissipation has been profound.
Funding: This research "BaSyMo" (FKZ: 03ET6087I) was funded by the Federal Ministry of Economic Affairs and Energy (BMWi).

Acknowledgments:
The authors gratefully acknowledge the financial support provided by BMWi, the German Federal Ministry of Economic Affairs and Energy, for funding of the project BaSyMo (FKZ: 03ET6087I). We like to thank our project partner for their excellent collaboration on the topic of battery research: ElringKlinger AG.

Conflicts of Interest:
The authors declare no conflict of interest. The founding sponsors had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, and in the decision to publish the results.

OCV Profiles
For the use of interpolated functions within COMSOL Multiphysics ® , a reduced set of nodes and corresponding function values are required. Within the simulation platform, the piecewise cubic interpolation method is considered within the node interpolation functions. The method generates a Hermite polynomial h(x; h(x i ) = y i ∀i ∈ I), I index set, which preserves the shape of the data, respects monotonicity, and has a continuous first derivative. The interpolation values of the positive electrode profile are created based on the coefficient data taken from literature [25]. The derivation of the negative electrode profile is described in Section 4.1.