A semi-empirical, electrochemistry-based model for Li-ion battery performance prediction over lifetime

Abstract Predicting the performance of Li-ion batteries over lifetime is necessary for design and optimal operation of integrated energy systems, as electric vehicles and energy grids. For prediction purposes, several models have been suggested in the literature, with different levels of complexity and predictability. In particular, electrochemical models suffer of high computational costs, while empirical models are deprived of physical meaning. In the present work, a semi-empirical model is suggested, holding the computational efficiency of empirical approaches (low number of fitting parameters, low-order algebraic equations), while providing insights on the processes occurring in the battery during operation. The proposed model is successfully validated on experimental battery cycles: specifically, in conditions of capacity fade >20%, and dynamic cycling at different temperatures. A comparable performance to up-to-date empirical models is achieved both in terms of computational time, and correlation coefficient R2. In addition, analyzing the evolution of fitting parameters as a function of cycle number allows to identify the limiting processes in the overall battery degradation for all the protocols considered. The model suggested is thus suitable for implementation in system modelling, and it can be employed as an informative tool for improved design and operational strategies.

Predicting the performance of Li-ion batteries over lifetime is necessary for design and optimal operation of integrated energy systems, as electric vehicles and energy grids. For prediction purposes, several models have been suggested in the literature, with different levels of complexity and predictability. In particular, electrochemical models suffer of high computational costs, while empirical models are deprived of physical meaning. In the present work, a semi-empirical model is suggested, holding the computational efficiency of empirical approaches (low number of fitting parameters, low-order algebraic equations), while providing insights on the processes occurring in the battery during operation. The proposed model is successfully validated on experimental battery cycles: specifically, in conditions of capacity fade > 20%, and dynamic cycling at different temperatures. A comparable performance to up-to-date empirical models is achieved both in terms of computational time, and correlation coefficient R 2 . In addition, analyzing the evolution of fitting parameters as a function of cycle number allows to identify the limiting processes in the overall battery degradation for all the protocols considered. The model suggested is thus suitable for implementation in system modelling, and it can be employed as an informative tool for improved design and operational strategies.

Introduction
Due to increasing concerns about global warming, greenhouse gas emissions, and depletion of fossil fuels, increasing efforts have been put in reducing environmental issues on a global scale. In the automotive industry, electric vehicles (EVs) have been widely accepted and promise a significant reduction in CO 2 emissions [1]. Moreover, the energy conversion is shifting from conventional, fossil fuels-based to distributed and based on renewable sources, especially wind and photovoltaic. However, due to the intermittent nature of renewable sources, their integration in the grid might challenge its reliable and stable operation. Among other solutions, battery electrochemical energy storage (BEES) has been proposed for power and energy buffering, and for (primary) frequency regulation [2].
Li-ion batteries (LIBs) are considered an interesting option for BEES, and are already evaluated in several demonstration projects worldwide [2]. This is due to their exceptional energy-to-weight ratio, lack of memory effect [3], fast response, high power capability, long-cycle lifetime at partial cycles, low self-discharge rates [2], and decreasing manufacturing costs [4]. Plus, BEESs pose less stringent performance requirements compared to EVs [5], allowing to re-use traction LIBs in the grid either as "second-life" [6], or as vehicle-to-grid [7]. Both options present potential technological, economic and environmental opportunities for improving energy systems [6], while mitigating the negative impact of high EVs penetration [7].
In designing and operating systems (either electric drive trains, or BEESs in grids of any size), system costs need to be evaluated [8,9] in early stages, and safe and optimal operation guaranteed [10]. In order to do so, several parameters describing the state of the batteries during operation need to be known, especially state of charge (SOC), state of health (SOH), end of discharge (EOD) and end of lifetime (EOL) [11,12]. In most cases, these parameters cannot be measured, since only on-line voltage, current and temperature readings are generally available [10]; therefore, they have to be inferred from model-based estimation [11]. The chosen model should hold its validity across the entire lifetime of the batteries, i.e. take into account the impact of nonlinear and usage-dependent ageing effects [13] on the above mentioned parameters, either directly (through following the time evolution of ageing-related parameters [12]) or indirectly. At the same time, computational efficiency of such model is paramount, due to the need of real-time results [11], and to the small computational power generally allocated in EVs [11] and in simulating BEES-containing energy grids [14].
Four different classes of models can be found in the literature for modelling LIBs: (1) empirical models (mathematical equations and equivalent circuit models (ECMs)), (2) electrochemical models (among others, the single particle model, and Doyle-Fuller-Newman model [15]), (3) multi-physics models (as for example, electrochemical thermal modelling [16,17]) and (4) molecular and atomistic models [3,12]. The increasing level of detail brings to a higher computational cost [3]: as a result, multi-physics models and molecular models are generally not an option for real-time predictions. Even electrochemical models are considered computationally inefficient [18], due to their complexity (further increased when taking ageing into account [19]), the need of solving partial differential equations, and the large number of unknown input parameters [11,20]; as a result, they are seldom employed for real-time predictions [12]. Furthermore, the required parameterization is complex and time-consuming [8,19]. In order to improve the computational efficiency of electrochemical models, some simplified solutions have been suggested, either by transforming the partial differential equations into ordinary differential equations [12,20], or by using lumped parameters [12]. Although the computational efficiency is successfully increased, differential equations still need to be solved, and most of all, the number of input parameters to be determined for each individual battery remains high [20].
On the other hand, empirical models allow for very fast computation [3]. ECMs, although effective in describing the dynamic performance of LIBs [11], often underestimate the influence of SOC on the battery behaviour (particularly Thevenin-based and impedance-based models [11]). This is problematic when different cycling regimes, and ageing effects need to be described [21]. In order to mitigate this drawback, ECMs require the input of open circuit potential (OCP) and internal resistance of the battery as a function of different SOC. These values are generally not provided by the manufacturer, and collecting them experimentally often implies interrupting the cycling protocol (in some cases, as for OCP curve building, for prolonged time [22]).
Alternatively, low-order algebraic equations are the simplest, yet effective approach for describing the terminal voltage of a battery as a function of SOC and current [21]. These models are characterized by a low number of fitting and input parameters; in addition, the latter can be obtained directly from the supplier's specifications [23]. Classical mathematical models are Shepherd, Unnewehr universal and Nernst models [24]. Although these models are designed to describe different battery chemistries [25] at different stages of lifetime, they suffer of poor accuracy (5-20% error [21]), which can be improved by the addition of extra parameters [21]. A modified Shepherd-based model is suggested in [21] for lead-acid batteries discharge curves, and proved effective in [23] for Ni-metal hydride batteries as well. A similar Shepherd-based model is used for LIBs, and implemented in the modelling of energy grids in [14,26].
The state-of-the-art literature thus offers several mathematical models, which can successfully describe LIBs performances. These models are inexpensive and computationally light, easily implemented, and able to mathematically reproduce the voltage profile during cycling; as such, they can be considered good predictive models, in the sense that they can successfully replicate the experimental terminal voltage response with a small error. However, due to its empirical nature, the mathematical function describing the terminal voltage holds little to no physical meaning, making these tools poor descriptive models; as a result, its fitting parameters are hard to generalize, and no information can be retrieved on the impact of the operating conditions on the cell's lifetime. Such information is precious for quick failure detection, battery sizing and design, and optimized operational strategies targeting a longer battery lifetime. Therefore, what is missing is a model with the computational efficiency of an empirical approach (thus easily implementable in the description of larger systems), and the added value of providing physically meaningful information about the processes in the battery.
Addressing this gap in the literature, the present work aims to: description of the physico-chemical phenomena in LIBs. As a consequence, also its fitting parameters would hold a physical meaning, and could be used to track LIBs performances. The degree of simplification is chosen following a bottom-up approach typical of system engineering, where the simplest, fundamentally strong model that produces accurate enough predictions to address the objectives is chosen [3].
• The performance of the proposed model should at least equal the most up-to-date empirical models for LIBs, in order to maintain the high computational efficiency and low number of input parameters.
• A thorough validation should be performed under different conditions of current, temperature and SOC. This is due to the fact that all models, at any level of description, require a proper training before being able to perform accurate predictions.
• Finally, the fitting parameters in different operational conditions should be physically meaningful, and follow the trends reported in the literature. This would allow to use such parameters to decouple the principal ageing phenomena, and to pinpoint the limiting ones.
The article is organized as follows: Section 1 is the introduction; Section 2 describes the models in use; Section 3 presents the experimental datasets chosen for this work; Section 4 contains the results and the related discussion; finally, Section 5 draws the conclusions of the work.

Modified Shepherd model
The model chosen, described in [23], is one of the most recent Shepherd modifications proposed in the literature, and was successfully employed for predicting LIBs performances in [14,26]. It is formulated as follows: being E the predicted voltage (in V), I the current (in A), V 0 the OCP value provided by the manufacturer (in V), Q nom the nominal capacity of the cell (in Ah), and Qt = ∫ It the actual battery capacity (also in Ah). R ohm (in Ω) is the ohmic resistance (or direct current (DC) resistance) of the battery, due to current collectors, separator, and electrolyte ionic conductivity, and it is not altered significantly during cycling. The alternate current (AC) resistance (more commonly known as internal resistance [27]) is qualitatively pictured by the polarization resistance K (in Ω), one of the fitting parameters. The other two fitting parameters are the exponential zone amplitude A e (in V), and the exponential zone time constant inverse B (in A −1 h −1 ); they are contained in an exponential term for improved fitting at high depth of discharge (DOD).
The polarization resistance term in Eq. (1), similarly to the classical Shepherd formulation, is a linearized Butler-Volmer equation, valid for overpotentials η ≤ 0.03 V [25]. The polarization voltage term, an addition to the classical Shepherd model, is introduced for describing the OCP-SOC relation, initially represented by a constant value (V 0 ).
Among the approximations taken in Eq. (1), the following are worth mentioning: • The whole polarization in LIBs is described through linearized Butler-Volmer equation, thus neglecting the contributions due to mass limitations. Although this approximation is often taken for experiments performed at modest C-rates (being C-rate the necessary current for discharging the nominal capacity of the battery in 1 h) [28], limitations in Li + solid-state transport are characterized by a strong SOC dependency [29], and can occur at different stages of lifetime [30,31] due to the overall small values of Li + diffusion coefficients in layered oxides [32].
• In addition, the polarization voltage and the exponential term in Eq.
(1) have a strong empirical character. This feature makes modified Shepherd a poor descriptive model, since it implies that the fitting parameters K, A e and B cannot be easily connected to specific processes occurring in the battery. As a consequence, their evolution as a function of cycle number cannot be used to investigate the degradation of LIBs performance.
• Finally, there are no terms which are explicitly dependent on temperature.
Due to these limitations, Eq. (1) is suitable for mathematically predicting the terminal voltage of LIBs, but no further analysis is possible, due to the empirical nature of its formulation.

Proposed model
The proposed model is formulated as follows: with E 0 being the steady-state voltage (in V), IR ohm the ohmic drop (in A Ω = V), η kin the overpotential towards charge transfer (in V), and η diff the overpotential towards mass limitation (also in V). The separated mathematical equations for each term of Eq. (2) are the following: • Steady-state voltage, E 0 : Eq. (3) describes the OCP-SOC relation in terms of a modified gas lattice isotherm (one of the simplest models to describe the thermodynamic behaviour of LIBs [33]). The last term, based on mean-field theory, is a correction to account for mutual interactions of intercalated Li + ions in the host electrode materials. J is the interaction parameter (dimensionless), and it is strongly impacted by the presence (or absence) of phase change(s) in the materials considered [33]. Since a cumulative J is used, instead of separate ones for the positive (which, in case of layered oxides, display no phase change) and negative (graphite, with multiple phase changes) electrodes, their partial contribution in terms of lattice interactions cannot be decoupled.
• Kinetic overpotential, η kin : Eq. (4) describes the charge transfer overpotential through a linearized Butler-Volmer equation, whose derivation is reported in Appendix A. K kin (in Ω) is a lumped parameter, called kinetic resistance, which represents the hindrance towards charge transfer.
• Diffusive overpotential, η diff : Eq. (5) describes the overpotential towards mass transfer as an (approximated) semi-infinite linear diffusion in a current-controlled experiment, with quasi-reversible one-electron process [34]. The approximation taken consists in using an average Li + concentration instead of separate bulk and surface Li + concentrations, which are hard to quantify experimentally. The assumption of semi-infinite linear diffusion in LIBs is experimentally supported by their electrochemical impedance spectroscopy (EIS) response: specifically, by the occurrence of a Warburg tail in the low frequency interval [27,35]. Since pure electrochemical reversibility in LIBs only occurs in limited C-rate intervals [36,37], and it worsens with ageing, quasi-reversibility is assumed for a more generalized description of the process. K diff (in Ω) is a lumped parameter, called diffusive resistance, which represents the hindrance towards mass transport.K kin and K diff together represent the internal resistance of the battery, divided into its kinetic and diffusive contributions. This is due to the fact that the internal resistance is actually an impedance (Z(Im) ≠ 0) in the mid-low frequency range (∼1 Hz) [27], thus positioned in a region of mixed kinetic control where both charge transfer (occurring in the mid frequency range) and mass limitations (occurring in the low frequency interval) are influential [35].
The proposed model directly addresses the limitations previously pinpointed for the modified Shepherd model: • The polarization in LIBs is now described by two terms, specifically η kin and η diff , which separate and address the kinetic and diffusive components of the total overpotential η.
• All the terms contained in Eq. (2) are described by (simplified) equations derived from electrochemical theory. Consequently, all fitting parameters hold a physical meaning, and can be linked to different processes in LIBs. For this reason, the proposed model can be considered a better descriptive tool than modified Shepherd.
• A simple, explicated temperature dependency is present for all the terms in Eq. (2), in the form of an RT/F coefficient.
The proposed model has the same number of fitting parameters as the modified Shepherd (three), and it requires the same input parameters (V 0 , R ohm , Q nom ), which are inferred from the supplier's specifications. The evolution of the fitting parameters J, K kin and K diff is assumed to capture LIBs performance over cycling, without the need of additional terms describing specific degradation phenomena.

Modelling tools
Both modified Shepherd and the proposed model were implemented in Matlab® (The MathWorks Inc.), version 2017a. During validation, the fitting parameters for both models were extracted by means of least square minimization method, using Levenberg-Marquardt algorithm as solver. The default values for maximum iterations (400), function tolerance (1e−6), optimality tolerance (1e−6), step tolerance (1e−6), and maximum number of function evaluation (100 × number of variables) were chosen. Each charge/discharge was fitted separately, and returned a set of fitting parameters as output. The parameters thus obtained were employed for dynamic simulation, and the goodness of the fit was reported in the form of the correlation coefficient R 2 . The average run time for both models was calculated on a set of 10 separate runs.

Experimental datasets
The chosen datasets were obtained from open-access on-line data repositories of LIBs from different manufacturers, which implies different chemistry and formulation of electrodes and electrolyte. All the batteries had the most up-to-date positive electrodes on the market [38]: layered oxides (LiMO 2 , where M is a transition metal, or a combination of) [39,40], characterized by similar structural [41] and electrochemical properties [42]. Consequently, the voltage curves of the considered LIBs were all characterized by a sloped profile, typical of these materials [43].
The datasets were selected in order to represent common testing conditions as cycle testing, and load testing. The first consists of repeated charge/discharge cycles until the capacity drops below 80% of the initially rated value, and it is performed to assess battery lifetime and reliability. The second can be performed according to several standardized protocols, all apt to determine if the battery can deliver the required power for a certain application.

Cycle testing
The experimental data were obtained from the National Aeronautics and Space Administration (NASA) Open Data Portal, Li-ion Battery Ageing Datasets. 1 The dataset was composed by experiments performed on commercial 18,650 LIBs with 2 Ah nominal capacity (Q nom ). The values for R ohm = 0.047 Ω and V 0 = +3.6 V vs. Li were also provided. Charging was carried out at room temperature in a constant current (CC) mode at 1.5 A (0.75 C) until the battery voltage reached 4.2 V, and then continued in a constant voltage (CV) mode until the charge current dropped to 20 mA. Discharge was carried out in CC (2 A, i.e. 1 C) until the battery voltage fell to 2.7 V. The cycling was repeated until 30% capacity loss was reached (from 2 A h to 1.4 A h).

Load testing at different temperatures
The experimental data were obtained from the website of the Center for Advanced Life Cycle Engineering (CALCE) at University of Maryland. 2 The batteries employed in the dataset were commercial INR 18650-20R LIBs with 2 Ah rated capacity (Q nom ). The values for R ohm = 0.06 Ω, and for V 0 = +3.6 V vs. Li were also specified. The batteries were cycled according to the Dynamic Stress Test (DST) protocol, a simplified driving cycle standardized by the United States Advanced Battery Consortium (USABC). Three temperatures were considered: 0, 25 and 45°C, maintained constant over cycling by means of a temperature chamber [39]. Since beginning of life (BOL) characterization was available only at 25°C, R ohm , V 0 and Q nom were assumed constant for all temperatures.

Cycle testing
In Fig. 1, the dataset is fitted with both the proposed model and modified Shepherd. An accurate fitting is required for correctly estimating SOC, SOH and EOD, and it is overall obtained, with a slightly higher correlation for modified Shepherd (R 2 = 0.985 in Fig. 1). Both models underestimate EOD with an average error of 0.3 mV, which decreases as a function of cycle number. Thus, both models are characterized by a higher accuracy at EOL ( Fig. 1(a)), which indicates that they are both suitable in describing a capacity fade higher than 20%. Therefore, the proposed model can successfully replace modified Shepherd for prediction purposes, without significantly compromising the accuracy of the fit, or the computational time (below 10 s for both models).

Load testing at different temperatures
The datasets at 0°C (Fig. 2), 25°C (Fig. 3), and 45°C (Fig. 4) are fitted with both the proposed model and modified Shepherd. The proposed model displays a slightly higher correlation at room temperature (R 2 = 0.976 in Fig. 3), a comparable one to modified Shepherd at 0°C (Fig. 2), and slightly inferior at 45°C (R 2 = 0.876 in Fig. 4). Even if an error is clearly introduced, especially at 45°C, its average (0.05 V) is equal to ∼1% of the experimental voltage, which is a rather accurate result in the framework of empirical models [21]. The high correlation coefficients (R 2 ≥ 0.9) and the fast run time (< 10 s) make both the proposed model and modified Shepherd good predictive tools at different temperatures. However, the temperature dependence is rather simplified in the proposed model; therefore, to increase its level of description towards this variable, two Arrhenius-based temperature corrections are proposed and discussed in the Supporting Material. Both modifications come at the price of increased input or fitting parameters, thus negatively impacting the complexity of the model, and its computational efficiency. In addition, neither of the two modifications improve the predictability of the proposed model (due to the overall small impact of the Arrhenius terms on the output voltage); for such reasons, they were not considered further in the present work.

Cycle testing
The capacity lost in Fig. 1(a) is displayed in Fig. 5. The cells are repeatedly cycled to high DOD (SOC < 20%, in Fig. 1), which can negatively impact LIBs lifetime [44]: as a result, 30% capacity fade is reached in less than 200 cycles. Based on the capacity fade rate, Fig. 5 can be divided in two regimes of mild (I) and more pronounced degradation (II). A rapid change in the capacity fade rate indicates the coexistence of several degradation phenomena [45]. In addition, both regimes are characterized by a nearly linear trend, indicating the contribution of both electrodes to the overall degradation (in case of dominating Solid Electrolyte Interface (SEI) growth, a square-root-oftime dependency is observed instead) [44].
In such situations, decoupling the ageing processes, and identifying the dominating one(s), is particularly challenging.

Interaction parameter J.
The J term in Eq. (3) captures the sloped voltage curves of LiMO 2 -containing LIBs [33]. In Fig. 6(a), J displays a different sign between charge and discharge (since the slope varies depending on the direction of the cycling), and approximately the same absolute value; slight differences are caused by the cycleinduced hysteresis present even at very low C-rates [46]. The slope of the voltage curves is SOC-dependent, due to a different Gibbs free energy (Δ r G) for Li + intercalation in LiMO 2 [47]. However, the SOC range is constant across the cycling in Fig. 1, and as a consequence, J is expected to remain unchanged. This applies only as far as regime I is concerned, while regime II displays a clear J increase (in absolute value) (Fig. 6(a)). As previously mentioned, a net change in J indicates a change in the lattice of the electrodes, although the proposed model cannot indicate in which of the two. Fig. 6(b), K kin displays the same value for charge and discharge, in agreement with the findings in [48]. On the other hand, K diff is characterized by different values between charge and discharge, with charge being the highest (Fig. 6(c)); this is caused by the inaccurate EOD prediction in Fig. 1(b), since in the SOC considered the highest diffusive limitation should occur in discharge, due to the high re-lithiation hindrance at high DOD [29].

Kinetic and diffusive resistances (K kin and K diff ). In
At both electrodes, the following steps need to take place for Li + intercalation to occur: (a) Li + (de)solvation from the electrolyte; (b) Li + transport through the electrode/electrolyte interface; (c) Li + ⇌ Li (one-electron charge transfer at the electrode/electrolyte interface) [30]. Depending on the operating conditions, one step of Li + intercalation mechanism might become rate-determining [30]; the results in      Fig. 6(b) and (c) indicates it to be Li + solid-state diffusion, since K kin is characterized by a constant value across the two regimes, while a progressive increase in K diff is observed (particularly in regime II). The possible causes behind the increasing hindrance towards Li + solid-state diffusion are here discussed. In first approximation, layers formation on the electrode's surface cannot be considered as the main ageing mechanism; otherwise, a change in K kin would be observed [30]. The hindrance towards diffusion has to be connected to changes in the bulk of the electrodes instead, as the change in J seems to indicate. At low SOC, graphite suffers of mechanical stress, induced by one of the phase transitions associated with its staging [44]. The occurrence of mechanical stress increases the activation energy barrier for diffusion processes in LIBs materials [49], which is compatible with the trends of J and K diff in Fig. 6, and the rapid capacity fade in Fig. 5.

Load testing at different temperature
Power failure in LIBs is generally induced by impedance rise, which increases the total overpotential; as a consequence, the battery reaches the low-voltage cut-off prematurely. In Fig. 7, the power fade of the cells tested under DST protocol is displayed at 0, 25 and 45°C. While the voltage profiles are overall similar, an increasing polarization can be observed which decreasing temperature, particularly in the discharge phase of the cycling. Consequently, the cell at 0°C, characterized by the highest polarization, reaches the low voltage cut-off more rapidly than at 25 and 45°C.
The fitting results at room temperature are first considered, and a comparison at different temperatures follows.
All fitting parameters display a current-dependent spread in Fig. 8; K kin and K diff are directly dependent towards C-rate [29], while J is indirectly dependent of I, since a higher current (applied for a fixed amount of time) leads to a larger SOC range covered during cycling. Overall, the SOC interval in Fig. 3 (∼ 80%) is narrower than in Fig. 1 (0  -80%), and as a consequence, slightly different values of J are obtained (Figs. 6(a) and 8 (a)). J is characterized by constant values until the final cycles, where a rise is observed. K diff behaves similarly (Fig. 8(c)), while K kin steadily increases across the entire cycling ( Fig. 8(b)): consequently, charge transfer seems the limiting process in the power fade considered. For this reason, and for the presence of higher currents in the DST protocol (compared to Fig. 1), K kin values are also significantly higher than in Fig. 6(b). A major cause for impedance rise is the build-up of layers at the electrodes' interfaces [50]. Although charge transfer is generally considered a fast process, it can become rate-limiting depending on the electrode materials and electrolytes involved [30], and a large increase in charge transfer resistance has been observed in LiMO 2 -containing LIBs, due to layers formation at the positive electrode interface [51]. These layers are both inorganic (structural reconstruction at LiMO 2 interface) and organic (electrolyte decomposition induced by O release from LiMO 2 lattice) [51]. This phenomenon could also explain the  Fig. 1(a) with the proposed model.   Fig. 3(a) with the proposed model. increase in J (due to local phase changes at the electrode surface [51]) and K diff (since the reconstructed layer is unsuitable for diffusion [50]) in the final stage of cycling.
Also temperature can affect the rate-determining step in Li + intercalation mechanism [30]. For such reason, the fitting parameters as a function of temperature are presented in Fig. 9.

Interaction parameter J.
Since J is associated with the thermodynamic properties of LIBs, no significant influence of temperature (at least, in the experimental temperature range considered) is expected. This is due to the fact that the steady-state voltage E 0 can be defined as: where Δ r G, Δ r H and Δ r S are the Gibbs free energy, the entropy and enthalpy of the considered reaction (in this case, Li + intercalation). Due to the small value of Δ r S, the effect of temperature on E 0 is mostly negligible for Li + -intercalating electrodes [47]. This is matched by the results in Fig. 9(a), except those at 45°C, affected by a higher error in the fitting (Fig. 4).

Kinetic and diffusive resistances (K kin and K diff ).
The dependency of all polarizations towards temperature is ruled by Arrhenius equation [52], which is re-written as a function of the polarization resistances in [30]. According to this law, higher temperatures are associated with smaller hindrances; even considering the error associated with the fitting in Fig. 4(a)), both parameters follow the expected trend ( Fig. 9(b) and (c)), with the highest hindrance at 0°C, and the lowest at 45°C. This indicates that the simplified temperature dependence in the proposed model can describe the evolution of polarization resistances as a function of temperature to a satisfactory extent. At all temperatures, the trend in K kin indicates that charge transfer remains the limiting process. This seemingly underlines a common ageing phenomenon at all temperatures, exacerbated at low temperature by the slow rate of all undergoing kinetic processes. This would explain the similarity of the voltage profiles in Fig. 7, and the temperaturedependent duration of the cycling.

Conclusions
In the present work, a semi-empirical model is suggested for predicting LIBs terminal voltage under different operational conditions. The major findings of the study are the following: • A semi-empirical model was formulated, entirely based on (simplified) equations describing physical phenomena in LIBs; in this way, all its fitting parameters hold a physical meaning. This makes it a better descriptive model than other representatives of the empirical class.
• The proposed model matched the performance of the up-to-date modified Shepherd chosen, proving to be as good a predictive model. The comparison was performed considering common cycling conditions for LIBs, specifically cycle testing and load testing. During cycle testing, the proposed model suitably described a capacity fade to 30% of the nominal value, with a decreasing instantaneous error indicating higher accuracy at EOL. During load testing (DST protocol was chosen as representative of this category), the proposed model captured power failure dynamics at 0, 25, and 45°C, with slightly worsened performances at the highest temperature.
• The evolution of the fitting parameters under different testing conditions was analyzed, and matched with literature reports. In addition, the variation of the fitting parameters across the different cycling protocols was justified by considering their theoretical dependency towards SOC, current and temperature. This not only confirms their physical meaningfulness, but it is also a promising result towards generalizing the fitted values to other LIBs and cycling protocols.
• Although simplified, the degree of description of the proposed model allowed for decoupling the main ageing processes, and indicated the limiting ones both for a complex cycling protocol at different temperatures (as DST), and for a rather complex capacity fade dynamic. The capacity fade during cycle testing was characterized by a complex ageing rate profile, and a rapid degradation; the proposed model indicated Li + diffusion as the limiting process, correlated with thermodynamic changes. This matched the occurrence of a phase change in graphite at low SOC, responsible of mechanical stress and diffusive hindrance. Regarding load testing, the model indicated charge transfer as the limiting process at all temperatures, possibly caused by layer formation on the surface of the electrodes; the literature indicated the positive electrode as the main responsible in such regard.
In conclusion, a semi-empirical model was delivered, with input parameter and computational cost comparable to up-to-date mathematical models: this makes it a suitable candidate for real-time predictions. The suggested model is thus suitable for both predicting and describing LIBs performances in the framework of a larger system (like EVs or BEES), with the added value of providing insights into the battery dynamics and degradation. This information can be used for quick failure detection, battery sizing and design, system design and operational strategy optimization.

Declaration of interests
The authors have no competing interests to declare.

A.1 Derivation of Eq. (4)
η kin represents the overpotential associated with one-electron charge transfer occurring in LIBs: where Θ represents a site in the electrodes' lattice. Butler-Volmer equation is commonly used in LIBs for expressing η kin [15]: kin kin (8) where i is the kinetic current density (in A m −2 ) associated with Li + intercalation, i 0 the exchange current density, and α the barrier factor. 'c' and 'a' denotes cathodic and anodic reaction respectively. i 0 depends on the concentration of the species involved in Reaction (7), and for one-electron charge transfer, assuming that the thermodynamic equilibrium can be described by Nernst equation [53]: where k c represents the kinetic rate constant for the cathodic reaction, k a for the anodic reaction, c a and c c the concentrations of oxidized and reduced species in Reaction (7). It is important to underline that the assumption of Nernstian equilibrium in LIBs electrodes is a simplification, as most Li + intercalating materials do not fully obey Nernst law [33]. Before combining Eqs. (8 and 9), the barrier factor α is assumed to be symmetric (α = 0.5), and both equations are expressed in terms of I = i × A (where A is the area), since current densities are difficult to quantify [54]. After combining the equations, the following is obtained: Note that in Eq. (10) the following formalism is adopted: Eq. (10) can be re-written as: Note the lumping 2k kin = K kin . At BOL, the positive electrode is fully lithiated ((1 − x) = 1), and the negative electrode fully de-lithiated ((1 − x) = 0). Assuming a charging process, and following the evolution of the positive electrode (which is the source of cyclable Li + in LIBs [29]), a complete Li + depletion of the positive electrode delivers the maximum theoretical capacity (Q max ), according to Faraday's first law of electrolysis [53]: (with M being the molecular weight). In practice, the complete depletion of the positive electrode is avoided due to lattice instability; therefore, Q max in Eq. (13) is replaced by the practical maximum capacity available, which is assumed to be the nominal one (Q nom ). During cycling (incomplete Li + extraction), Eq. (13) has the form: where 0 < (1 − x) < 1 is the fraction of sites occupied by Li + . Knowing that SOC at the time t is defined as [14]: by substituting Eqs. (13 and 14) in Eq. (15): Substituting Eq. (16) in Eq. (12), the following is obtained: which is the non-linearized form of Eq. (4).

A.2 Validity of linear Butler-Volmer approximation
Eq. (12) In order for the linear approximation to be valid, the experimental values of η kin (calculated according to Eq. (4)) should be included in this interval. The separate contributions of each term in Eq. (2) to the total calculated voltage E are displayed in Fig. 10. The dataset chosen for this analysis is the one in Fig. 3(b), due to its broad current range across cycling. Even at the highest current applied (−4 A, 2 C), η kin assumes small, positive values (1e−4), thus fulfilling the requirements of Eq. (18) for linearization.

Appendix B. Sensitivity analysis
A sensitivity analysis is performed, in order to evaluate the impact of a variation in the fitting parameters on the calculated terminal voltage. The method chosen is a stochastic approach based on Monte Carlo Analysis, with 1e + 6 number of samples. The parameters are assumed to have a uniform distribution, and the interval for their variation corresponds to ± 50% of the fitted value for the cycle considered. The cycle testing in Fig. 1 has been chosen as representative, due to the broader SOC range covered.
A randomly chosen charge/discharge curve from Fig. 1(a) is displayed in Fig. 11, located approximately halfway the cycling. 11 (equally spaced) points are selected, and for each of them the maximum and minimum voltage variation is calculated, as well as the correlation between the voltage and the fitting parameters (quantified by means of R 2 in Table 1).  Fig. 11. Calculated voltage for a charge/discharge cycle of Fig. 1(a) (randomly chosen), and associated voltage limits, quantified with Monte Carlo analysis. The specific points of the curve used for the sensitivity analysis are marked with numbers.
Very low values of R 2 in Table 1 indicate the overall low impact of K kin on the calculated voltage, at all stages of cycling. Therefore, the variation of either J or K diff (or both) is expected to impact the voltage the most. At the beginning of cycle, J displays a very high correlation coefficient (R 2 ∼ 1) for both charge and discharge. However, this value steadily decreases as the cycle proceeds, until R 2 = 0.32 is reached at end-of-charge, and R 2 = 0.17 at EOD. R 2 for K diff perfectly mirrors the trend for J: while low correlation is observed at early stages of cycling, R 2 = 0.68 is obtained at the end of charge, and R 2 = 0.83 at EOD. J appears to be the most influential parameter on the calculated voltage at early stages of the cycle, possibly due to the steep voltage variation in the experimental data. Both J and K diff contribute to the voltage variation in the majority of the SOC range covered; finally, at the end of cycling, K diff contribution on the voltage variation dominates, particularly at high DOD (SOC < 20%).

Table 1
Correlation coefficient R 2 for the fitting parameters and the calculated voltage, quantified for the marked points in Fig. 11