Parameter Identification and Sensitivity Analysis for Zero-Dimensional Physics-Based Lithium-Sulfur Battery Models

As the performance of Lithium-ion (Li-ion) batteries approaches its practical limit (up to 500 Wh/kg) [1], researchers are actively seeking alternative, high energy density solutions. The Lithium-Sulfur (Li-S) chemistry is one such solution. It has attracted significant attention from the scientific community for at least four reasons. First, the Li-S chemistry offers very attractive theoretical limits for specific capacity (1672 Ah/kg) and specific energy (2600 Wh/kg) [1]. Prototype Li-S cells have already achieved specific energies around 700 Wh/kg [2]. Second, Li-S batteries have the potential to offer a range of operating temperatures from −40 to 60 °C [3]. This is especially beneficial for low temperature applications. Third, phenomena such as precipitation and the internal charge shuttle effect provide Li-S batteries with some degree of intrinsic overcharge protection [4,5], thereby making them less vulnerable to catastrophic failure compared to their Li-ion counterparts. Fourth, the use of sulfur, instead of rare earth materials, makes Li-S batteries quite appealing in terms of production cost and environmental footprint [6]. In spite of these benefits, the Li-S chemistry suffers from low cycle life and relatively high self-discharge rates [7,8]. Efforts are ongoing to improve the electrochemical performance of Li-S batteries [9]. This paper is part of a complementary effort to improve Li-S battery performance through modeling, estimation, and, ultimately, optimal control. The paper focuses specifically on the problem of parameterizing Li-S battery models from experimental data: a critical first step towards model-based control.

The literature already presents a number of Li-S battery models that can be used for parameter estimation and model-based control. These models fall on a spectrum of fidelity and complexity levels. Equivalent circuit models (ECMs) of Li-S batteries sit at one end of this spectrum. They typically combine empirical open-circuit voltage maps with (often charge- and temperature-dependent) resistor-capacitor circuits [10–13]. Typical uses of such simple models include online state of charge estimation [14–16]. At the other end of the spectrum, one finds physics-based models of the reduction-oxidation, diffusion, precipitation, and self-discharge dynamics in Li-S batteries. These are typically partial differential algebraic equation (PDAE) models of coupled diffusion-reaction dynamics [17–22]. They have the advantage of providing higher-fidelity representations of the underlying battery physics, at the expense of high complexity. Zero-dimensional (0D) models provide an important middle ground between these two modeling extremes. They neglect the diffusion of species inside Li-S batteries, but use the laws of electrochemistry to model the underlying battery reactions. This approach results in lower-order models with fewer parameters to identify and reduced computational complexity. Simplified two-reaction and multi-reaction 0D models introduced in Refs. [23,24] are able to capture the main features of battery voltage response during discharge. The computational complexity of a 0D Li-S battery model depends on the specific choice of which reactions to model, and which to ignore. One potential approach for making this choice is to compare the fidelity with which different 0D models can fit a given experimental data set. Such an exercise is appealing because it can furnish both: (i) a systematic comparison between different battery models using the same data set, and (ii) an experimental parameterization of different 0D Li-S battery models. To the best of the authors’ knowledge, both of these contributions remain relatively unexplored in the literature. The overarching goal of this work is to address this research gap.

Starting from previous work on 0D Li-S models, this work goes one step further to analyze different 0D models capturing different underlying redox reactions. The goal is to answer the following questions: (i) How does each reaction influence the discharge voltage curve in a 0D model? (ii) Do all the model parameters have a significant impact on the discharge voltage curve? Is it necessary to identify them all? (iii) How accurately can the identified models fit the experimental voltage behavior? Is there a model with a specific reaction chain that performs best? The ultimate goal is to provide guidance on the choice of which reactions to model for future applications, such as state of charge estimation.

The remainder of this paper addresses these questions, and is organized as follows. Section 2 describes the 0D model structure adapted from Refs. [23,24], combined with the reduction reactions considered in each model. Section 3 provides a simulation study illustrating the effects of different reactions on the discharge voltage performance. Section 4 performs a parameter sensitivity study to provide guidance on which parameters are critical to identify from experimental data. Section 5 presents the parameter identification results using experimental data. Finally, Sec. 6 summarizes this paper’s conclusions. These conclusions are based on an experimental study at a discharge rate of 0.3C. At this discharge rate, the discharge voltage curve of an Li-S battery exhibits three main features, namely: a high plateau, a low plateau, and a dip point separating these two plateaus. In the high plateau region, the active material, S₈, in the cathode accepts electrons to produce the polysulfide, $Sx2−$ (x can be 8, 6, 4). Further polysulfide reduction takes place in the lower voltage plateau region [25]. Figure 1 shows two discharge characteristic curves, obtained by cycling Li-S cells fabricated as part of this study at 0.25C and 0.125C. Both of these curves exhibit the three main features discussed above, namely: the high plateau, low plateau, and dip point. Moreover, the small differences in discharge voltage between these curves can be partially attributed to simple polarization effect, such as Ohmic polarization. Thus, the paper’s focus on a discharge C-rate of 0.3C can help shed a light on behaviors observed over a range of C-rates, especially C-rates below 0.3C.

This section describes the 0D reaction models used for subsequent analysis. The models build on two existing bodies of work. First, research by Zhang et al. [24] develops a lumped model describing the dynamics of dissolved species concentrations and cathode material porosity in Li-S batteries. Moreover, this model includes the full reduction reaction chain and the effective electrolyte resistance. Second, research by Marinescu et al. [23] introduces another 0D model taking into account the dynamic evolution of different species’ masses as well as the shuttle effect. This model simplifies the chain of sulfur reduction reactions into two steps. Our work combines the underlying structures of these two models and expands them to capture the dynamics of different sets of reduction reactions. The end result is a set of four 0D models ranging from low to high complexity. The simplest model covers a minimum set of reduction reactions modeled in the literature, whereas the most complex model covers the full set of sulfur reduction reactions in Li-S batteries. The reactions considered in each model are listed in Table 1. In this work, the authors adopt the following modeling assumptions from the literature: (i) there is an unlimited amount of lithium supply in the cell with a negligible overpotential on the anode side, as in Ref. [18]; (ii) the shuttle effect of polysulfides is not included due to this paper’s focus on the voltage performance instead of capacity fade [24]; and (iii) only the lowest polysulfide’s precipitation reaction (⁠$2Li++S2−⇌Li2S↓$⁠) is modeled [23,24].

Since the model is zero-dimensional, there exists no mass transport due to diffusion/migration. The only dynamics changing the masses of the dissolved species are the electrochemical reaction and precipitation. These dynamics lead to the time-evolution relations for the masses of various sulfur species and the porosity of the cathode material, which act as the state variables in the model. Below is the resulting differential algebraic equation (DAE) model (assuming discharge), including both state equations and algebraic constraints. Descriptions of the model’s variables and parameters are included in the paper’s Nomenclature section.

The rates of mass change for higher-order dissolved sulfide species are governed by Eq. (1). For species S²⁻, one needs to consider its mass generation from the last reduction reaction and its mass loss due to precipitation, as shown in Eq. (2). The nucleation and growth phenomenon is described by Eq. (3) [23]. The rate of change of precipitate mass is partly driven by the precipitate mass itself. This reflects the fact that the existing precipitate serves as a nucleus for further precipitation/growth, as long as the mass in the electrolyte is above a given saturation mass.

Relative porosity is one of the state variables in the model with the dynamics governed by Eq. (4). This variable equals the current porosity of the cathode material divided by initial porosity and has a direct effect on the active reaction area, as shown in Eq. (7) [17]. When the porosity decreases to zero, all the reactants are blocked by the precipitate from the cathode material surface. This can be one of the indicators of the cell’s full discharge, another full discharge scenario being the reduction of all polysulfides to S²⁻. As the cell approaches complete coverage of the cathode (ɛ → 0), larger overpotentials (η_j) are needed for the same current. This provides one explanation of the low measured voltage at the end of Li-S battery discharge.

The reduction potential of each reaction is given by the Nerst equation (Eq. (5)), assuming that $Ej0$ is the equilibrium potential when the concentration of the dissolved species in reaction j is 1 mol/L [24]. The current generated by the corresponding reduction reaction is described by the Butler–Volmer equation (Eq. (6)). All these currents sum to the only input variable, namely, the external discharge current I, as shown in Eq. (8). Finally, Eq. (9) states the relation between the voltage measured across the battery, the overpotential, and the reduction potential of each reaction.

In the literature, there exist at least two simplifications of the reaction chain on the cathode side, namely, the two-step reduction $S8→S42−→S2−$ [26] and the four-step reduction $S8→S82−→S62−→S42−→S2−$ [20]. Different Li-S battery models, capturing different representations of the overall reaction chain, are built to suit the focus areas of the corresponding research studies. In this section, the authors compare the simulated voltage outputs of the models in Table 1 under a 0.3C discharge and analyze each reaction’s influence on the different features of the discharge voltage curve. openmodelica [27] is used to build and simulate the models with the key parameter values listed in Table 2. The simulated discharge voltage-capacity curves are compared in Fig. 2.

This analysis provides the characteristics of the voltage curve considering different reactions. The bottom line is that simplified representations of the reduction reactions are reasonable in terms of one’s ability to capture theoretical specific capacity, but the influence of different simplifications on the voltage features is significant. One needs to choose the reactions carefully based on the purpose of the given application. For example, studies focused on the capacity fade may choose simpler models, while state of charge estimation may require more complicated reaction models due to the accuracy required to capture the characteristics of output voltage.

One can broadly categorize the parameters of the above models into two groups. The first group includes parameters that are known a priori, such as Faraday’s constant, gas constant, temperature, molar mass, and number of electrons transferred in each reaction. The second group contains unknowns such as standard potentials, exchange current densities, initial masses of the dissolved species, and precipitation-relevant parameters (i.e., the morphology parameters γ, ω, and precipitation rate). Increasing the number of reactions and species considered in the model increases the number of unknown parameters. Hence, the difficulty of parameterization increases. In this section, the authors analyze the sensitivity of the 0D model’s output to the underlying parameters in order to better understand which parameters are critical to identify.

Instead of using a theoretical parameter sensitivity approach [28], this work applies a similar analysis method to the work presented in Refs. [29–31] for a 1D Li-S battery model. Specifically, the authors choose one parameter at a time and investigate the voltage performance changes due to variations in only this parameter. All other parameters are held constant. This parameter sensitivity analysis is performed on all models listed in Table 1. Since Model #4 describes the full reaction chain containing the largest number of parameters (total number of 14), the results of Model #4 are chosen for the following sensitivity discussions. Cell voltage sensitivities to different parameters are plotted for two perturbations from the nominal values.

Figure 3 presents the results of analyzing the cell voltage sensitivity of Model #4 to parameter $Ej0$⁠. One can see the relative dominance of different reactions by observing the corresponding changes in the voltage curve. Reaction $S8→S82−$ (j = 1) has a dominant influence on the high plateau voltage level. Reaction $S82−→S42−$ (j = 2, 3) has a strong effect on the “ending slope” of the high plateau and the dip point between the two plateaus, while slightly changing the low plateau voltage level. The last two reactions (j = 4, 5) occur during the low plateau region, and therefore both of the corresponding standard potentials influence the voltage curve dramatically.

The sensitivity analysis for $ij0$ is shown in Fig. 4. One can observe minor changes in the voltage curve within the chosen range of $ij0$⁠. The exchange current densities of the middle reactions have negligible impacts on voltage, as the voltage curves show a negligible difference. Compared to $Ej0$⁠, which influences the voltage V through the variable E_j, the parameter $ij0$ affects the overpotential η_j though an exponential function. On the one hand, this shows that $Ej0$ and $ij0$ have similar effects on voltage and highlights the challenge of identifying both $ij0$ and $Ej0$ simultaneously from experimental data. On the other hand, bigger variations in $ij0$ produce smaller changes in the voltage curve compared to $Ej0$⁠, which suggests the exchange current density has less of an impact on the model’s voltage performance.

Four parameters essentially govern the precipitation phenomenon, namely, the saturation mass of dissolved S²⁻, the precipitation rate k_p, and the morphology parameters of the relative porosity of the cathode, γ and ω. The variations of each parameter produce changes in the low plateau region without influencing high plateau voltages, as shown in Fig. 5. When the saturation mass increases, more S²⁻ is dissolved in the electrolyte, reducing the potential of the reaction $Sx2−→S2−$⁠, hence providing a lower cell voltage. The precipitation rate k_p and morphology parameter γ determine the flatness of both the beginning and end regions of the low plateau. Larger values of k_p and smaller values of γ provide a flatter low plateau. The parameter ω is the only parameter dramatically affecting the capacity of S₈. With an increase in the rate of change of relative porosity, the discharge process can terminate earlier because of the resulting paucity of active reaction area inside the cathode material.

The authors assume that all the active sulfur is dissolved into the electrolyte at the beginning of discharge. The mass of the dissolved S₈ is the main factor determining the total amount of active material inside the battery, and the initial masses of other species at the beginning of discharge are comparatively very small. In this case, sensitivity analysis is only performed with respect to the initial mass of dissolved S₈. The results are presented in Fig. 6 with negligible impacts on the voltage curve and specific capacity. However, one must note that the voltage curve difference emerges when plotting the voltage with respect to time instead of specific sulfur capacity. Hence, the dip point between the two plateaus changes with respect to time.

From this analysis, the following conclusions are made to guide parameter identification: (i) Although both $Ej0$ and $ij0$ affect the voltage curve, V is more sensitive to $Ej0$⁠. (ii) In the set of parameters relevant to precipitation, correctly capturing the values ω and γ is more essential for capturing the shape of the low plateau. (iii) The initial mass of dissolved S₈ influences the timing of the dip point between the high and low plateaus. Hence, in the parameter identification study, the authors focus on parameters $Ej0$⁠, γ, ω, and $mS80$⁠, with all the other parameters treated as known constants.

It is difficult to systematically compare multiple models of Li-S battery dynamics without fitting all of them to the underlying experimental data first. This section presents a parameterization study focusing on the proposed four reaction models. The authors obtain the test data through discharging laboratory-made prototype coin cells.

The authors fabricated 2016-type coin cells using the recipe in Ref. [32]. The cathode material consisted of Carbon/Sulfur composite with 70 wt% sulfur, Ketjenblack EC60JD to form a conductive framework, and polyvinylidene fluoride (PVDF) dissolved in N-Methyl-2-pyrrolidone (NMP) as a binder. A lithium chip of 0.6 mm thickness served as the anode. The electrolyte used here was 1M LiTFSI and 4 wt% LiNO₃ in the dioxolane/dimethoxylethane mixture (DOL/DME2 = 1:1, V/V). The sulfur loading of the cell was calculated at 0.521 mg. A constant discharge current of 0.03 mA was applied in the experiment to obtain the discharge voltage curve. Testing was conducted at room temperature, using ARBIN battery test equipment with special attention given to the accuracy of current control. In particular, a separate capacitor charge/discharge test was used for estimating the bias in the cycler’s current measurements, and this bias was used to correct measured data.

The identification results are listed in Table 3, with the corresponding estimated voltage shown in Fig. 7. The identified parameters for all the models are able to fit the experimental data in terms of the total discharging time. The objective function values corresponding to the best estimated parameters for the models are shown in Fig. 8(a). Interestingly, Model #3 performs the best fit with a minimum objective function value among all the models. The over-simplification of reactions in Model #1 leads to failure in capturing the slope of the voltage in the high plateau region. In contrast, the full reaction Model #4 does not provide a good fit in the middle region where the precipitation phenomenon emerges. This is potentially due to the simplification of the precipitation dynamics in the model. A more sophisticated physics-based precipitation model may improve the accuracy of the full reaction model.

Another observation pertains to the simulation time needed for one discharge cycle, as shown in Fig. 8(b). Intuitively, the computational complexity of the models increases with the number of included reduction reactions. However, the execution time increases from Model #1 to Model #3 are minor, whereas the increase from Model #3 to Model #4 is quite significant. One can conclude the best fit for further control application is Model #3, as it provides the least discrepancy between the experimental data and simulation results, without excessive computational cost.

This paper performs parameter sensitivity analysis and identification for a variety of zero-dimensional physics-based Li-S battery models. The authors provide answers to the questions listed in Sec. 1: (i) The intermediate reactions between S₈ and $S42−$ add slope to the second half of the high plateau, while the reactions between $S42−$ and S²⁻ prolong the flat portion of this plateau and deepen the inter-plateau voltage dip. (ii) Sensitivity analysis provides guidance on the final choice of to-be-identified parameters as $Ej0$ (the voltage level), γ (the low plateau slope), ω (the low plateau duration), and $mS80$ (the dip location). (iii) Taking into account the intermediate reactions affecting the high plateau, the estimated voltage curve fits the experimental data well except at the beginning of discharge. Finally, taking into account the importance of balancing parameter identification performance versus execution time, the authors recommend Model #3 as the best model for use in future control oriented applications under this model structure.

This work is funded by National Science Foundation Grant 1351146. The authors gratefully acknowledge this support. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the National Science Foundation.

There are no conflicts of interest.

i =

species index

j =

reaction index

q =

number of species

p =

number of reactions

v =

cell volume (L)

F =

Faraday’s constant (C/mol)

I =

input (discharge) current (A)

R =

gas constant (J/(K mol))

V =

output voltage across the cell (V)

T =

temperature (K)

a_v =

active reaction area (m²)

$av0$ =

initial active reaction area (m²)

i_j =

current generated by reaction j (A)

$ij0$ =

exchange current density (A/m²)

k_p =

precipitation rate constant (1/(g s))

m_i =

mass of species i (g)

$mSp$ =

mass of precipitated sulfur (g)

$mi0$ =

initial mass of species i (g)

$nSi$ =

number of sulfur atoms in species i

n_j =

number of electrons exchanged in reaction j

s_i,j =

stoichiometric coefficients

E_j =

reduction potential for reaction j (V)

$Ej0$ =

standard reduction potential (V)

M_S =

molar mass of sulfur (g/mol)

S_sat =

saturation mass of S²⁻ (g)

η_j =

overpotential for reaction j (V)

ɛ =

relative porosity

γ =

power of the relative porosity

ω =

relative porosity change rate constant (1/g)