**Leveraging Cell Expansion Sensing in State of Charge Estimation: Practical Considerations**

#### **Miriam A. Figueroa-Santos \*, Jason B. Siegel and Anna G. Stefanopoulou**

Department of Mechanical Engineering, University of Michigan, Ann Arbor, MI 48109, USA; siegeljb@umich.edu (J.B.S.); annastef@umich.edu (A.G.S.)

**\*** Correspondence: miriamaf@umich.edu

Received: 11 April 2020; Accepted: 14 May 2020; Published: 22 May 2020

**Abstract:** Measurements such as current and terminal voltage that are typically used to determine the battery's state of charge (SOC) are augmented with measured force associated with electrode expansion as the lithium intercalates in its structure. The combination of the sensed behavior is shown to improve SOC estimation even for the lithium ion iron phosphate (LFP) chemistry, where the voltage–SOC relation is flat (low slope) making SOC estimation using measured voltage difficult. For the LFP cells, the measured force has a non-monotonic F–SOC relationship. This presents a challenge for estimation as multiple force values can correspond to the same SOC. The traditional linear quadratic estimator can be driven to an incorrect SOC value. To address these difficulties, a novel switching estimation gain is used based on determining the operating region that corresponds to the actual SOC. Moreover, a drift in the measured force associated with a shift of the cell SOC–expansion behavior over time is addressed with a bias estimator for the force signal. The performance of Voltage-based (V) and Voltage and Force-based (V&F) SOC estimation algorithms are then compared and evaluated against a desired ±5% absolute error bound of the SOC using a dynamic stress test current protocol that tests the proposed estimation scheme across wide range of SOC and current rates.

**Keywords:** state-of-charge estimation (SOC); linear quadratic estimator; lithium ion battery; iron phosphate; cell expansion; force

#### **1. Introduction**

The primary function of the battery management system (BMS) is to provide an accurate state of charge (SOC) estimation. The SOC represents the amount of charge in ampere-hours (Ah) remaining in a cell divided by its total capacity [1,2]. The BMS traditionally uses current, voltage, and sometimes temperature measurements to estimate the SOC to plan future actions and to prevent over-charging or discharging of cells. Generally, manufacturers provide conservative estimates of remaining energy, since an overestimation of SOC can leave the vehicle stranded. In the case of unmanned air vehicles (UAV), overestimation of SOC might prevent the vehicle from safe landing, since landing maneuvers require very high power, which typically cannot be achieved at very low SOC levels [3]. Underestimating SOC, on the other hand, wastes valuable resources and adds cost and weight to the vehicles, which is critical for robotic platforms.

From the lithium ion batteries, lithium ion iron phosphate (LFP) has been considered for UAV, hybrid electric vehicles (HEV), and electric vehicles (EV) due to their capacity for fast charging, high power capability, and long cycle life [4]. LFP batteries consists of graphite as the negative electrode and lithium iron phosphate (metal oxide) as the positive electrode [5]. Due to the relative flat half-cell potential the positive electrode exerts (also known as the phase-separating cathode active material) [6], the open circuit voltage (OCV) has a relatively flat slope through most of its operating

SOC range (10–95%), as shown in Figure 1. The phase transitions in the graphite material correspond to the different voltage plateaus with respect to SOC, as shown in Figure 1a. This makes SOC estimation difficult under noisy environments [7,8] and inexpensive sensing, such as robotic and automotive applications. Previous work has suggested strain or stress (pressure or force) measurements to augment terminal voltage for SOC estimation [9]. Specifically, the graphite in the negative electrode expands when the lithium ions are intercalating into it, and the positive electrode contracts as the lithium ions leave it causing a change of thickness in the components of the battery [10]. This change in thickness causes the battery to swell. Therefore, the overall observed cell swelling is the summation of the swelling from the positive and negative electrodes. When the cells are constrained to a fixed displacement, as typical in automotive packs, the battery swelling results in an increased force on the fixture. This swelling force can be measured using a load cell [9].

**Figure 1.** Measured voltage and force for the 20 Ah A123 lithium ion iron phosphate battery cycled under low current rate *C*/20. (**a**) Measured discharge (blue), charge (red), and average of discharge and charge (black) voltage for a *C*/20 cycle. For military robots such as the Packbot shown above, SOC estimation is critical to avoid the robot getting stranded. Lithium ion iron phosphate (LFP) batteries are typically used for operation of this robot due to the high power required. (**b**) Measured discharge (blue), charge (red), and average of discharge and charge (black) force for the under current rate *C*/20. The experimental battery fixture consists of: (1) load cell (force sensor); (2) movable plate; (3) lithium ion battery; (4) two aluminum end-plates; (5) temperature sensor. For generalization of the results to the other cell sizes, the force measurements can be converted to pressure by diving the force by the surface area of our battery (*Ab*) , which is the width (*wb* = 0.161 m) times the length (*lb* <sup>=</sup> 0.227 m) of the battery, *Ab* <sup>=</sup> *wb* <sup>×</sup> *lb* m2.

The structural changes of this expansion have been studied from the electrode mechanics point of view with respect to strain/OCV coupling [10,11]. The overall volume change of this expansion results in a monotonic function of SOC for a nickel manganese cobalt (NMC) graphite cell [12]. By measuring the force produced by the expansion and including it in the estimation algorithm, improved SOC accuracy can be achieved as compared to voltage based methods [9,13]. The greatest benefits were observed in the 30–50% SOC region, where the voltage slope was relatively flat, but also at the low SOC level, where voltage drops very fast and is challenging to have an accurate model as the battery ages. The change in measured force vs. SOC over many cycles can also inform better estimates of the battery state of health (SOH) [14]. For the LFP-graphite cell studied here, the anode (negative electrode which is graphite) expands during charging but the simultaneous higher contraction rate of the cathode (positive electrode which is LFP) results in a combined cell contraction in the middle SOC range. The overall result is a non-monotonic F–SOC behavior, as shown in Figure 1b. Due to the F–SOC non-monotonicity, SOC estimation based on measured force is challenging since multiple

SOC equilibrium points can correspond to the same measured force. Additionally, all mechanical measurements (force, pressure, or displacement) even after calibration exhibit drift associated with small changes in material and a shift in the SOC–expansion behavior as the battery ages. Loss of cyclable lithium causes a shift in capacity and a change in the SOC–expansion behavior of the cell. This capacity loss shifts the stoichiometric ratio associated with lithium concentration in each electrode and hence changes the electrode expansion [15] as a function of lithium intercalation (Coulombs stored) and the measured force/pressure versus SOC as the unknown drift addressed here. Moreover the LFP pouch cells used here are supported by poron sheets [16] between the cells instead of the spacers used in [13]. Thermal expansion of the battery and fixture and viscoelastic response of the compliant poron pad introduces an additional aging and drift factor. Predicting and modeling this aging behavior requires extensive resources and is currently by-passed by estimating this unknown bias as proposed in this article. As can be appreciated, the drift is a general problem of measuring the mechanical behavior (force, pressure, or displacement) of all batteries and it is not just a battery chemistry (LFP)-related issue.

Battery cell balancing is critical to extend the range of battery powered vehicles, the pack operating lifetime, and charge and discharge power limits [17]. To achieve fast and accurate cell balancing, accurate SOC estimation is needed [18]. Voltage-based SOC estimation for the LFP chemistry is particularly challenging since the voltage is relatively flat with respect to SOC. The proposed method improves the SOC estimation using leveraging information about the cell expansion and contraction during charging and discharging. Cell-to-cell variability due to aging, as well as the resulting changes in the measured force, was not investigated, but should be the focus of future research. The influence of cell balancing or imbalanced cells in a module whose force is measured should also be further investigated following the initial work by [19]. Since there would typically be only one force (or strain) measurement in a pack of series connected cell, cell balancing techniques [20] will be of high importance.

In this paper, we demonstrate the improvement in the SOC estimation of LFP batteries by using mechanical in addition to electrical measurements that can be implemented in packs or modules of both hard encased and pouch cells [21]. A novel solution to the multiple equilibrium SOC points is proposed based on the piecewise linear (PWL) F–SOC characteristic approximation that is further used in a switching gain model-based linear quadratic estimator (LQE) design that consists of a combination of force and voltage measurements [22]. Due to the drift that appears in most force measurements, as shown in Figure 2, the SOC may not be accurately estimated. Therefore a bias state is added to the LQE in order to capture the drift in the force measurement due to un-modeled changes in battery swelling or creep of the plastic materials. Experimental validation is also performed on the model-based (LQE) controller design using the combinations of force and voltage measurements during realistic battery electric vehicle usage profiles including the Dynamic Stress Test (DST). The performance of a controller designed with a "perfect" model is compared to one with model mismatch in the OCV and F–SOC PWL fit. The simulated model mismatch captures the typical modeling uncertainties or changes in the cell expansion and open circuit voltage due to aging [23].

#### **2. Experimental Setup and Force Behavior**

The battery considered here for the experimental validation is an A123 20 Ah lithium iron phosphate pouch cell with a voltage range of 3.6–2 V. The fixture, as shown in Figure 1b, consists of an active battery cell (3) (with a temperature sensor on top of it (5)) and a dummy (inactive) cell with a compliant rubber pad in between. The active battery cell and inactive cell consists of a laminated aluminum pouch cell and the rubber pad is a Poron 4701-30 from Rogers, which is 1.14-mm thick. The temperature sensor is a multimeasurand GE sensor that consists of three resistance temperature detectors (RTDs) and one eddy current expansion sensor [21]. The dummy cell and the rubber pad is used as a stress absorber to emulate the conditions a cell might experience in a pack configuration. The purpose of the dummy cell is to simulate the compliance of the whole system. The stiffness of poron is much less than that of the battery, thus the compression of the dummy cell is negligible in

terms of the whole system. Therefore, the dummy cell can be used without complicating the force measurement. The dummy cell is held tight by an aluminum bottom end-plate (4), and a movable plate is placed on top of the active cell. This plate has one degree of freedom in the vertical direction with the force sensor (1) placed on top of it. On top of the force sensor, the fixed top aluminum end-plate is bolted to the bottom end-plate to simulate the behavior of a constrained battery pack with fixed distance between the end plates. The force sensor is an Omega (LC305-500) load cell sensor (strain gauge type). The sensor has a 2225 N full scale range with an accuracy of 4.45 N. The load cell and voltage are digitized and recorded by Data Translation DT-9828, which has a voltage accuracy of 2 millivolts. In this paper, the force sensor is used because it is cheaper than the displacement sensors [15], but also less accurate. For example, the accuracy of the displacement sensor used in [15] is 1 μm, which corresponds to 0.35 N of compressive force on the poron, while for the Omega force sensor used the accuracy is listed as a percentage of full scale range (4.45 N). The force only accuracy is the important reason for integrating force and voltage information along with performing the bias estimation.

**Figure 2.** Drift in the force sensor can be observed by comparing the force vs. SOC for two *<sup>C</sup>* <sup>20</sup> cycles which were conducted two months apart with three cycles done in between them with the same battery. The force measurements have a minimum of 3 N and a maximum of 6 N drift across the entire SOC range. The drift could be caused by thermal expansion of the battery, pad, and fixture [13], capacity change, or a combination of all four throughout the life of the cell and module on which measurements are performed. Adjustments in the fixture and changes in the preload will cause larger changes in measured force but can also be modeled as a sensor drift. For generalization of the results to the other cell sizes, the pressure is shown.

The force plotted against SOC (F–SOC) and the open circuit voltage vs. SOC (OCV–SOC) are measured experimentally using a pulse–relax profile. Specifically, a CCCV charge is applied followed by a pulse–relax discharge profile, with the current rate of *<sup>C</sup>* <sup>20</sup> for 12 min (which results in a 1% SOC change) followed by a 1 h of rest period to eliminate the influence of the internal resistance and electrolyte polarization of the battery. The data points at the end of the rest are used to obtain the discharge F–SOC and OCV–SOC. The test is then repeated using charge pulses of equal duty cycle. The average between the measurements at the same SOC from the charge and discharge datasets is obtained for the F–SOC, as shown in Figure 1b. The average F–SOC shown in Figure 1b is modeled in Section 3 and is the best fit parameters that match the force inflection points. An average fit is also obtained from data for OCV–SOC and shown in Figure 3a. The values of the fit for the OCV–SOC (average of discharge and charge) and F–SOC can be found in Appendix C.

An important consideration is the drift observed in the measured force between repeated tests as shown in Figure 2. Two average force cycles with an applied current of *<sup>C</sup>* <sup>20</sup> at different test dates are shown. These cycles have approximately a minimum 3 N and a maximum 6 N of drift in between the second test (performed on 06/2017) and first test (performed on 04/2017). This measurement was taken with the same cell. The drift could be caused by thermal expansion of the battery, pad, and fixture, capacity change or a combination of all four throughout the life of the cell and module on which measurements are performed. That is the reason we have treated the observed bias as an unknown variable. We assume that this unknown/uncertain-origin drift evolves slower than the force measurement from the charge/discharge changes and we therefore can estimate it as an unknown constant bias in real time. This drift in force can affect the SOC estimation if not compensated and, thus, estimation of the sensor bias is used to improve the practical force based SOC estimation.

**Figure 3.** A piecewise linear approximation of the (**a**) open circuit voltage and (**b**) cell swelling force was fit to the experimental measured values for the A123 Lithium Iron Phosphate Battery Cell. The gray dots correspond to the average of measured data at each SOC point from the charge and discharge cycles at the *<sup>C</sup>* <sup>20</sup> rate (see Figure 1). The orange line represents the PWL fit. The horizontal solid lines represent the inflection points of the cell swelling force. The horizontal dashed black lines represent the inflection points of the cell swelling force used for the developed observers. *RI*, *RI I*, and *RIII* represent the force slope regions determined by our PWL model. For the PWL force model, we assume *bL* = *c*<sup>3</sup> and *bH* = *c*<sup>5</sup> since the changes in voltage and force slopes are due to the intrinsic phase transitions in the material of the battery electrodes [15]. Assessment of the SOC estimator robustness is performed by imposing ˆ *bL* = *bL* and ˆ *bH* = *bH* to capture model mismatch for the force and voltage behavior. For generalization of the results to the other cell sizes, the pressure is shown.

For the SOC estimation development, tuning, and comparison, two models with two different levels of fidelity are used:


In the next section, the simulation model is detailed and its efficacy is highlighted in Figure 4 based on a modified DST cycle [24] that is scaled for a 20 Ah Li-ion battery, as shown in Figure 4a. The DST was chosen because it has a current profile that has the combination of the following C-rates and is representative of usage in an electric vehicle: C/4, C/2, 1C, and 2C. If the utilization of the electrode is relatively uniform, the C-rate should not influence the swelling significantly. The electrode expansion depends on the bulk concentration of the electrode solid phase, as opposed to the terminal voltage which depends on the surface concentration [25] and therefore is less sensitive to C-rate. The theory for determining up to what current density the electrode utilization is uniform can be found based on the porous electrode models by Fuller Doyle and Newman [26] and Newman and Tobias [27]. The largest expected contribution of C-rate dependence (or more precisely root mean square (RMS) current) on the result is through thermal expansion of the cell. The experimental profile consists of a Constant-Current/Constant-Voltage (CCCV) charging protocol at a rate of 1C until the battery is fully charged. After a rest period of 30 min, a 1C rate discharge current is applied until it reaches 61% SOC. After the second rest period of 2 h, the modified DST cycle is applied. The resulting voltage, temperature, force, and SOC are shown in Figure 4b–e, respectively. The SOC is calculated by Coulomb counting using a high resolution current sensor and assumed to be the true SOC.

**Figure 4.** Comparison of open loop simulation model results (denoted by the color orange) and experimental measurement of voltage, temperature, force, and SOC obtained from the 20 Ah battery during the Dynamic Stress Test (DST) inside an environmentally controlled chamber set at 25 ◦C ambient conditions. (**a**) Current profile scaled for the 20 Ah A123 battery. (**b**) Comparison of the open loop model terminal voltage and measurement. (**c**) Measured battery temperature. (**d**) Comparison of the open loop modeled force and measured force. After 30 min of rest, the un-modeled dynamics in force excited by cycling of the cell to decay to zero. For generalization of the results to the other cell sizes, the pressure is shown. (**e**) Comparison of state of charge (SOC) measurement.

#### **3. Simulation Model**

The models described here include the drift present in our force data and higher dynamics such as hysteresis present in our voltage. These models were simulated to validate the robustness of our developed estimators for analysis purposes before the experimental implementation. After this validation, the estimators were used with experimental data and their performance was evaluated. The battery model used for our simulations is presented in the following subsections. The SOC (*z*) is simply modeled as

$$\frac{dz}{dt} = -\frac{\mathcal{I}}{\mathcal{C}\_b},\tag{1}$$

where *Cb* is the cell capacity.

#### *3.1. Cell Swelling Force Model*

The force measured at the load cell is modeled using a static non-linearity *Fsim*(*z*), which is a function of state of charge, with additive bias and noise terms given by

$$F(z) = F\_{\text{sim}}(z) + \upsilon\_{\text{F}} + f\_d \tag{2}$$

where *F*(*z*) (N) is the measured force that relies on *Fsim*(*z*), which is the PWL model of the average F–SOC behavior; *vF* is the measurement noise; and *fd* is a constant drift or bias (assumed to be constant but unknown) value present in our force measurement. A piecewise linear representation of the average F–SOC behavior is given by

$$F\_{\rm sim}(z) = \mathbb{C}(z)z + \mathbb{C}\_0(z) = \begin{cases} a\_{m}z + a\_{m0}, & \text{if } z \le b\_{L} \\\ \beta\_{m}z + \beta\_{m0\prime} & \text{if } b\_{L} < z \le b\_{H} \\\ \gamma\_{m}z + \gamma\_{m0\prime} & \text{otherwise} \end{cases} \tag{3}$$

where *αm*, *βm*, and *γ<sup>m</sup>* are the slope parameters; *αm*<sup>0</sup> is the preload or the force sensed at zero state of charge; and *bL* and *bH* denote the SOC where a change in the sign of the slope in the PWL model occurs as shown by the solid gray vertical lines in Figure 3. These parameter values can be found in Appendix C and can be adjusted to simulate model uncertainty and mismatch. The parameters *βm*<sup>0</sup> and *γm*<sup>0</sup> in the PWL model are uniquely determined from the other parameters via constraints of piecewise continuity

$$
\beta\_{m0} = (\alpha\_m - \beta\_m) b\_L + \alpha\_{m0} \tag{4}
$$

$$
\gamma\_{m0} = (\beta\_m - \gamma\_m) b\_H + (\alpha\_m - \beta\_m) b\_L + \alpha\_{m0}.\tag{5}
$$

The operating regions *R* represent the force slope regions determined by our PWL model. The first region is defined as *RI* : *<sup>z</sup>*<sup>ˆ</sup> <sup>∈</sup> [0, <sup>ˆ</sup> *bL*), the second region is defined as *RI I* : *<sup>z</sup>*<sup>ˆ</sup> <sup>∈</sup> [<sup>ˆ</sup> *bL*, ˆ *bH*], and the third region is defined as *RIII* : *<sup>z</sup>*<sup>ˆ</sup> <sup>∈</sup> (<sup>ˆ</sup> *bH*, 1]. The operating regions are shown in Figure 3.

#### *3.2. Terminal Voltage Model*

The terminal voltage in volts is modeled as

$$V\_T = V\_{0\text{c}}(z) - IR - V\_1 - V\_2 - V\_h + \upsilon\_{V,} \tag{6}$$

where *R*[Ω] is the total equivalent series resistance, *I* is the discharge current applied to the battery, *V*<sup>1</sup> and *V*<sup>2</sup> are the voltages due to the two resistance and capacitance (RC) pairs, and *vV* represents the V measurement noise. The OCV characteristic, *Voc*(*z*), is SOC dependent and modeled using

$$V\_{\rm oc}(z) = V\_0 + d(1 - \exp(-fz)) + h\left(1 - \exp\left(-\frac{-k}{1-z}\right)\right) + \lg z + \sum\_{i=1}^3 a\_{\rm v,i} \arctan\left(-\frac{z - b\_{\rm v,i}}{c\_{\rm v,i}}\right), \quad (7)$$

where *V*0, *g*, *d*, *f* , *h*, *k*, *av*,*i*, *bv*,*i*, and *cv*,*<sup>i</sup>* are tuned parameters found in Appendix C. The electric equivalent circuit (EC) battery model [28,29] is used for the simulation in our observer validation. In this study, voltage hysteresis (*Vh*) [30] is also considered for the simulation model

$$\frac{dV\_1}{dt} = \frac{-V\_1}{R\_1C\_1} + \frac{I}{C\_1} \tag{8}$$

$$\frac{dV\_2}{dt} = \frac{-V\_2}{R\_2C\_2} + \frac{I}{C\_2} \tag{9}$$

$$\frac{dz}{dt} = -\frac{I}{\mathcal{C}\_b} \tag{10}$$

$$\frac{dV\_h}{dt} = -\left|\frac{\gamma\_h I}{\mathbf{C}\_b}\right| V\_h + \left|\frac{\gamma\_h I}{\mathbf{C}\_b}\right| H(z, \text{sgn}(I)) \tag{11}$$

where *V*<sup>1</sup> and *V*<sup>2</sup> are the voltages of the RC equivalent circuits, *R*<sup>1</sup> and *R*<sup>2</sup> are the resistors, *C*<sup>1</sup> and *C*<sup>2</sup> are the capacitors of the RC equivalent circuits, *Vh* is the hysteresis voltage, *γ<sup>h</sup>* is the hysteresis rate constant, and *H*(*z*, sgn(*I*)) is a function of SOC and the sign of current (sgn(*I*)) following [30]. The function *H*(*z*, sgn(*I*)) is taken to be half the difference between the charge and discharge OCV measurements, and the parameter values can be found in Appendix A. Although the EC model parameters depend on the battery's SOC and temperature, in this paper, we do not take this dependency into consideration. The constant parameters of the EC model can be found in Table 1. The dynamic equations developed for charge/discharge as well as the measurements in Equations (2) and (6) is used to simulate the battery behavior, and it is numerically discretized with a time step of *Ts* = 1 s.

**Table 1.** Battery Equivalent Circuit Parameters and its values.


#### **4. Voltage Model Parameterization**

Before using our model for SOC estimation, the average Force–SOC and OCV–SOC needs to be obtained and model parameterization of the equivalent circuit parameters is required using a pulsed current profile such as the Hybrid Pulse Power Characterization (HPPC) [28]. The pulse current was obtained and the nonlinear programming solver fmincon was used to find the parameters provided in Table 1. To verify that the parameters are correct, different initialization values were used. From the different initialization values used, fmincon converged to the parameters provided in Table 1.

#### **5. SOC Estimation Model**

The linear quadratic estimator (LQE) also known as the steady-state Kalman filter is used for state estimation. The goal of the observer is to find a gain *K* that converges the initial state to the true state of the system using linear filter equation with measurement error feedback

$$
\hat{\mathbf{x}}\_{t+1} = A\hat{\mathbf{x}}\_{t} + Bu\_{t} + \mathbf{K}(y\_{t} - \hat{y}\_{t}) \tag{12}
$$

where the estimated output equation is given by

$$\mathfrak{F}\_t = \mathbb{C}(\mathfrak{X}\_t)\mathfrak{X}\_t - \mathbb{C}\_0(\mathfrak{X}\_t) - Du\_t \tag{13}$$

with *C*(*x*ˆ*t*) the slope of the estimated measurement and *C*0(*x*ˆ*t*) the affine parameter based on our model and shown in Figure 3. The matrix *D* is the direct transition (or feedthrough) term. The error dynamics are governed by the eigenvalues of *A* − *KC*, which depends on the chosen gain *K* so that it is stable and achieves fast convergence of the SOC estimation error *e* = *z* − *z*ˆ. We use the Discrete Algebraic Riccati Equation (DARE) to find our gain *K* on all our developed observers. The DARE is defined as

$$P = APA^T - APC^T \left[ CPC^T + R \right]^{-1} CPA^T + Q \tag{14}$$

and solved for the estimation error covariance matrix *P*. In this equation, *Q* is the process noise covariance matrix (size *nxn*) and *<sup>R</sup>* <sup>∈</sup> <sup>R</sup>*ny* is the measurement noise covariance matrix (size *mxm*). The values of the diagonal elements of *R* are chosen based on actual sensor measurement noise variance and *Q* is tuned so that the desired transient is achieved. The solution of the estimation error covariance (*P*) for the Kalman filter converges to the solution of Equation (14) for *t* → ∞ if (*A*, *C*) is detectable. In this case, the asymptotically stable observer gain is then computed [31] as

$$K = \left[ \text{CPC}^T + R \right]^{-1} \text{CPA}.\tag{15}$$

#### *5.1. Voltage Only Observer Design*

In the case of voltage only estimation, we neglect the hysteresis dynamic term. This represents a typical model mismatch in voltage measurements. The states of the observer are *x*ˆ*<sup>t</sup>* = [*V*ˆ 1, *V*ˆ 2, *z*ˆ] *<sup>T</sup>* given by the discretized version of Equations (8)–(10) and the PWL approximation of the nonlinearities in voltage measurement *y*ˆ*<sup>t</sup>* is given by the equations in Appendix B. The values of the parameters in the voltage and force models can be found in Appendix C. The gain *K* for this observer is given as

$$K\_V = \begin{bmatrix} K\_1 & K\_2 & K\_3 \end{bmatrix}^T \tag{16}$$

The values for the gains *Ki* with *i* = 1–3 are obtained by tuning the Q matrix. The gain *KV* is found for the eight regions in which the voltage is divided: *R*1: *z* ∈ [0, *c*0], *R*2: *z* ∈ (*c*0, *c*1], *R*3: *z* ∈ (*c*1, *c*2], *R*4: *z* ∈ (*c*2, *c*3], *R*5: *z* ∈ (*c*3, *c*4], *R*6: *z* ∈ (*c*4, *c*5], *R*7: *z* ∈ (*c*5, *c*6], and *R*8: *z* ∈ (*c*6, 1].

The main challenge in this system is the slow convergence of the estimation error due to the almost zero output gain (in C) especially the C(1,3) that corresponds to the *dV dz* in Equation (18). Increasing the *Kv* gain to compensate for the low state to output gain C governed by *dV dz* will amplify voltage sensor noise. Therefore, another observer is developed that uses voltage and force measurement since SOC and bias in the force signal are unobservable by force measurement only.

#### *5.2. Voltage and Force Observer Design*

In the case of force and voltage estimation, the states of the observer are *x*ˆ*<sup>t</sup>* = [*V*ˆ 1, *V*ˆ 2, *z*ˆ] *T* and *y*ˆ = [*V*ˆ *<sup>T</sup>*, *F*ˆ] *<sup>T</sup>* given by the equations in Appendix B. The gains for the V&F observer have the following format

$$K\_{VF} = \begin{bmatrix} K\_{11} & K\_{21} & K\_{31} \\ K\_{12} & K\_{22} & K\_{32} \end{bmatrix}^{T} \tag{17}$$

and the gain *KVF* is found for the eight regions of the voltage, as explained in the previous subsection. The modeled force used for the state estimator, *F*ˆ, is given by

$$\hat{\mathcal{E}}\_{\rm sim}(\mathfrak{z}) = \hat{\mathcal{C}}\_{\rm dF}(\mathfrak{z})\mathfrak{z} + \hat{\mathcal{C}}\_{0}(\mathfrak{z}) = \begin{cases} a\_{m}\mathfrak{z} + a\_{m0}, & \text{if } \mathfrak{z} \le \mathfrak{z}\_{L} \\\ \mathfrak{z}\_{m}\mathfrak{z} + \mathfrak{z}\_{m0}, & \text{if } \hat{\mathfrak{b}}\_{L} < \mathfrak{z} \le \hat{\mathfrak{b}}\_{H} \\\ \gamma\_{m}\mathfrak{z} + \gamma\_{m0}, & \text{otherwise} \end{cases} \tag{18}$$

where *z*ˆ is the estimated SOC from our observer and ˆ *bL* and ˆ *bH* are the inflection points that are shifted by −10%, as shown by the dashed black lines in Figure 3 to emulate model uncertainty and aging. In the case of a "perfect" model, ˆ *bL* = *bL* and ˆ *bH* = *bH*.

#### Switching Logic in Observer Design

The gain (*K*) of the LQE depends on the relationship between the estimated state and the slope of the measurement. In the case of a monotonic function, the estimated state will converge to the true state value when the estimator gain is chosen so that the error dynamics are stable. In the case of a non-monotonic function, the estimated state and the true state can have different slopes in the force output depending on the operating region. Therefore, the estimator gain would have a wrong sign which would lead to the divergence of the estimated state of charge from the true state. This is due to the traditional LQE using the modeled slope at the region of the estimated state. For example, consider the case where the model initialization occurs in the middle section of the SOC range (force decreases and cell contracts as SOC increases), whereas the actual SOC is in the high SOC range (force increases and cell swells as SOC increases); the traditional approach will lead to divergence of the estimated SOC state from the actual SOC due to the difference in slope. To address this difficulty, an algorithm was developed that uses a window of past measurements force data to identify the slope of the non-monotonic F–SOC. The observer gain will need to switch based on a judicious combination of the information at hand, namely


Therefore, the output error injection gain K is a function of the state estimate and the estimated force derivative, e.g., *<sup>K</sup>*(*z*ˆ, *dF*˜ *dz* ), using two sources of information due to its importance in the convergence of the estimation error. The DFDZ is computed as a line fitting problem based on the moving window of past force measurements and the Coulomb counting based SOC integration (*z*˜) over the moving window. The fitted DFDZ line *F*˜ = *dF*˜ *dz <sup>z</sup>*˜ <sup>+</sup> *<sup>F</sup>*˜ <sup>0</sup> parameters are computed using the least-squares estimation as

$$
\begin{bmatrix} \frac{d\tilde{F}}{dz} \\ \tilde{F}\_0 \end{bmatrix}\_k = (\mathbb{L}\_k^T \mathbb{L}\_k)^{-1} \mathbb{L}\_k^T \mathbb{F}\_k \tag{19}
$$

where *F*˜ <sup>0</sup> is the affine parameter in Equation (18). The moving window of force measurements and the design matrix respectively are defined as follows

$$\mathbf{F}\_k = \begin{bmatrix} F\_{k-n} \\ F\_{k-n+1} \\ \vdots \\ F\_k \end{bmatrix}, \quad \mathbb{L}\_k = \begin{bmatrix} \mathbb{z}\_{k-n} & 1 \\ \mathbb{z}\_{k-n+1} & 1 \\ \vdots & \vdots \\ \mathbb{z}\_k & 1 \end{bmatrix} \tag{20}$$

where *k* represents the discrete-time measurement index and *n* is the number of past samples. The Coulomb counting based SOC integration is computed as

$$\overline{z}\_{k} = \overline{z}\_{k-1} - \frac{I\_{k-1} T\_s}{C\_b}. \tag{21}$$

The integration is initialized with *<sup>z</sup>*˜*k*−*<sup>n</sup>* <sup>=</sup> 0. The moving window provides the *dF*˜ *dz* computation as the average slope in a window of *n* prior values of force which causes a delay *δ* in the switching logic as a function of the window size *n*, as shown in Figure 5d. We can see that for *n* = 150 the *δ* is smaller compared to *n* = 450, but, with smaller *n*, the estimated slope is more susceptible to noise. We chose *n* here to be 150 samples with a discretized time of 1 s. The reasoning behind choosing this window size is explained in Section 5.3. Due to this delay and the fact that during this time the actual and the estimated state can be in different segments of the SOC, which could cause divergence, the gain is set to zero when the estimated slope *dF*˜ *dz* has a different sign from the modeled slope. Therefore, the system will run open loop for both voltage and force when there is slope mismatch. This is done in order to avoid instability issues.

**Figure 5.** Simulation (without bias in the force measurement) for the Switch Observer V&F developed in our previous work [22]. This figure shows the impact of high gains corresponding to solution of the linear quadratic estimator using *Q*<sup>1</sup> = diag([2, 0.1, 0.1]) and *R* = diag([1, 1]) with moving windows of length (MW) = 150 and MW = 450. If the gain is lowered by using *Q*<sup>2</sup> = diag([2, 0.1, 1*e* − 2]), as shown by the purple dashed line, the error is decreased by 5% with same 150pt moving window length. (**a**) Comparison of the simulated state of charge (SOC). (**b**) Comparison of the state of charge (SOC) error with the dashed lines representing the target ±5% bound. (**c**) Comparison of the simulated force. For generalization of the results to the other cell sizes, the pressure is shown. (**d**) Comparison of the true slope with the estimated slope *dF*˜ *dz* and observer output. Due to the greater delay (*δ*) before applying the zero gain for the MW = 450 case, a higher error in SOC estimation occurs during the transition from *RI I* to *RI*. (**e**) Comparison of the feedback gain *K*<sup>32</sup> from force to SOC based on the switching logic. When the observer and estimated slopes mismatch, the gain is set to zero.

The switching logic for our gain (*K*(*z*ˆ, *dF*ˆ *dz* )) is given by the following rules:

$$K\_{VF}(3,2) = \begin{cases} \begin{array}{llll} \mathcal{K}\_{\mathcal{R}\_I} & \text{if } \frac{d\mathcal{F}}{d\boldsymbol{z}} > 0 \quad \text{and} \quad \mathfrak{L} < \hat{\mathfrak{b}}\_L \\\ \mathcal{K}\_{\mathcal{R}\_{II}} & \text{if } \frac{d\mathcal{F}}{d\boldsymbol{z}} < 0 \quad \text{and} \quad \hat{\mathfrak{b}}\_L < \mathfrak{z} < \hat{\mathfrak{b}}\_H \\\ \mathcal{K}\_{\mathcal{R}\_{III}} & \text{if } \frac{d\mathcal{F}}{d\boldsymbol{z}} > 0 \quad \text{and} \quad \mathfrak{L} > \hat{\mathfrak{b}}\_H \\\ 0 & \text{otherwise.} \end{array} \tag{22}$$

Note that the gains *KRI* , *KRI I* , and *KRIII* are not constant. The set of gains *KRI* consists of the four gains that correspond to the four regions in voltage in region *RI*, *KRI I* consists of the two gains that correspond to the two regions in voltage in region *RI I*, and *KRIII* consists of the two gains that correspond to the two regions in voltage in region *RIII*, as shown in Figure 3. The gains are therefore a function of *z*ˆ and *dF*˜ *dz* and they switch to the corresponding region depending on the value of the estimated slope and estimated SOC.

#### *5.3. SOC Estimation Error during Switching*

Going open loop, during the time interval when there is a slope mismatch, is not sufficient to avoid divergence of the estimated state of charge with high feedback gain. To verify this, we analyze the Luenberger Observer for the SOC state using the measured force only. The goal of this analysis is to determine the impact of the gain *K* on divergence of the state estimate when the true model and state estimate are operating on opposite regions of the output non-linearity. Using Equations (12) and (13), we can write the error dynamic (*e* = *z* − *z*ˆ*t*) for the observer assuming *x*ˆ is in *RI* and *x* is in *RI I* as:

$$\pounds\_{t+1} = A\pounds\_t - K(y\_t - \hat{y}\_t) = A\pounds\_t - K(\beta\_m \mathbf{x} + \beta\_{m0} - \mathbf{a}\_m \mathbf{\hat{x}} - \mathbf{a}\_{m0} + \mathbf{a}\_m \mathbf{x} - \mathbf{a}\_m \mathbf{x}) \tag{23}$$

$$\mathfrak{e}\_{t+1} = (A - \mathcal{K}\mathfrak{a}\_m)\mathfrak{e}\_t - \mathcal{K}((\beta\_m - \mathfrak{a}\_m)\mathbf{x} + \beta\_{m0} - \mathfrak{a}\_{m0}).\tag{24}$$

From this error dynamic equation, we notice the error converges to a non-zero quantity so that *e* tends toward zero only when *x* = *bL*, and the steady state error *ess* = ((*β<sup>m</sup>* − *αm*)*x* + *βm*<sup>0</sup> − *αm*0)/*α<sup>m</sup>* is independent of the gain K. In the case, where the sign of the model is updated based on the measured slope of the force signal ( *dF*˜ *dz* ), the growth in SOC estimation error is bounded by the number of samples in the filter or the moving window size (*MW*). If the current is bounded (which it is in our case), then the divergence in the state estimate is bounded by the integral of the current and the switching time *δ*. When the estimator model is updated to the correct slope, the observer begins to converge again. However, noise in the measurement of force could still result in divergence of the estimate if the gain *K* is too high. This is shown in Figure 5. We can see that just by changing the *MW* = 150 to *MW* = 450 by using the same gains given by *Q* = diag([2, 0.1, 0.1]) and *R* = diag([1, 1]). This is because the gain is large and the delay in *dF*˜ *dz* crossing zero is 1 min greater than *MW* = 150. Therefore, the divergence in SOC, because of the slope mismatch, will grow for a longer period of time with increasing filter length and a higher error is achieved due to the delay in switching the gain to zero.

Now, if we decrease the gain by using *<sup>Q</sup>* = diag([2,0.1,1 × <sup>10</sup>−2]) and *<sup>R</sup>* = diag([1, 1]) without bias estimation and *n* = 150, we notice that at the same region the error decreases by 5%. Therefore, low gains should be used to avoid the divergence and a window size of 150 is chosen since this is the lowest window that provides sufficient noise rejection. The previous error dynamic analysis can be done for the other mismatch slope areas to determine the minimum SOC state estimation error.

#### *5.4. Bias Influence in SOC Estimation*

In this section, the proposed estimation is tested under biased force measurements. In Figure 6, the simulated force data have a bias or drift of 3 N present in the force measurement. This drift in force affects the SOC estimation, as shown in Figure 6. Application of the switched model observer based on voltage and the force measurements based on previous work [22] without taking into account the

bias estimation results in error larger than 5%. Therefore, we assumed a constant bias state ( ˆ *fd*,*t*) that is given by

$$f\_{d,t+1} = f\_{d,t} \tag{25}$$

**Figure 6.** Simulation assuming accurate modeling using Switch Observer V&F Bias observer with an emulated bias of 3 N and 10 % (0.1) initial SOC error. (**a**) Comparison of the simulated state of charge (SOC). Is clear from the results that Switch Observer V&F without bias estimation diverges from the simulated SOC. Therefore, bias estimation is needed in the developed observer. (**b**) Comparison of the state of charge (SOC) error with the dashed lines representing the target ±5% bound. SOC errors greater than 20% are obtained with the Switch Observer V&F without bias estimation. It is shown that with bias estimation the SOC error is within the 5% estimation error bound (EEB). (**c**) Comparison of the simulated force and observer outputs. The bias state permits deviation in the force output, without compromising SOC estimation. For generalization of the results to the other cell sizes, the pressure is shown. (**d**) Comparison of the simulated terminal voltage. (**e**) Comparison of the true slope with the estimated slope *dF*˜ *dz* and observer output. (**f**) Comparison of the feedback gain *K*<sup>32</sup> from force to SOC based on the switching logic. As shown in region A, denoted by the black dashed horizontal lines, when the observer and estimated slopes mismatch, the gain is set to zero. (**g**) The observer estimates the bias state which converges to the true bias as shown by the black horizontal dashed lines.

where *fd*,*<sup>t</sup>* is the constant bias or drift term that augments the electrical states. The estimator now estimated the electrical states (except hysteresis) and the force drift with *x*ˆ*<sup>t</sup>* = [*V*ˆ 1, *V*ˆ 2, *z*ˆ, ˆ *fd*,*t*] *T*. Therefore, the gains for the V&F Bias observer have the following format

$$K = \begin{bmatrix} K\_{11} & K\_{21} & K\_{31} & K\_{14} \\ K\_{12} & K\_{22} & K\_{32} & K\_{24} \end{bmatrix}^T \tag{26}$$

The values for *K*1*<sup>i</sup>* and *K*2*<sup>i</sup>* with *i* = 1–4 are obtained by tuning the Q and R matrices.

#### **6. Simulation Results without Model Mismatch in F–SOC and OCV–SOC**

The standard Dynamic Stress Test (DST) profile is repeated back to back and modified by adding a constant current to periodically recharge the battery at 1/6 C rate, exercising a wider range of SOC, as shown in Figure 4a. A measurement noise variance of 5 mV for voltage and 0.05 N for force was chosen based on the variance of the experimental data. As for the drift value, we chose a value of 3 N based on the monthly drift observed between repeated characterization experiments. Our observer works if the initial error in bias is within 3 N. The LQE estimator is initialized with an SOC error of ±10% in order to evaluate convergence. In the case of the DST cycle, the true SOC state may be approaching the estimated value or diverging from the true SOC value depending on the initial SOC estimation error. The objective is to stay within the ±5% estimation error bound (EEB) for SOC denoted by the dashed lines in the Figure 6b. We chose to simulate initial conditions around 60–80% and the SOC swing of the whole cycle around 20–80% SOC. This range was chosen due to the challenge associated with flat OCV–SOC profile present for voltage (around 40–60%) and the negative slope in F–SOC (around 35–70%). Therefore, the simulated "measured" data are initialized at 61% SOC. The weights chosen to obtain the gains *K*1*<sup>i</sup>* with *i* = 1–4 and *K*2*<sup>i</sup>* with *i* = 1–4 are shown in Table 2. The gain *K*<sup>32</sup> is shown in Figure 6f to illustrate when the algorithm uses the estimated slope mismatch to zero the observer gain and run open loop.

**Table 2.** Control weights for the different sensors using the simulated data. V, terminal voltage; V&F, fusion of both sensors; V&F Bias, fusion of both sensors with bias state.


The SOC estimation using the V measurement and V&F measurement with and without the bias state estimation are shown in Figure 6a. In all three cases, a 10% initial estimation error is assumed. The inflection points in the force with respect to SOC are denoted by the dotted horizontal lines in Figure 6a. In this case, the correct values for *bL* and *bH* are used for the simulated data and the observer. The switched model applies zero gain (as shown in Figure 6f in Region A) when the model slope (based on *z*ˆ) does not agree with the estimated slope using the least squares on the moving window with Equations (19) and (20). Using the previously developed V&F observer without bias estimation [22], the SOC error does not remain within the ±5% bound with a 3 N bias in the measurement, even though the error in the force signal is small, as shown in Figure 6a,c. For the proposed switching force-and voltage-based LQE with bias (Switch Observer V&F Bias) before *t* = 50 min there is an SOC estimation error of 10%, even though the force estimation is matching our simulated data, due the error in bias state estimation shown in Figure 6g. Between *t* = 25 and *t* = 42 min, there is a slope mismatch in *dF*˜ *dz* due to the delay of our moving window. After *t* = 42 min, the slopes match again, and the correct non-zero feedback gain is applied and the state of charge estimation error convergences within our ±5% EEB, as shown in Figure 6b. Even though both the proposed switching force and voltage based LQE (Switch Observer V&F Bias) and the LQE based on voltage (Observer V) estimate voltage

accurately (as shown in Figure 6d, it can be observed that the Switch Observer V&F Bias has faster convergence than the estimate based on V alone, as shown in Figure 6b. The faster convergence is due to the addition of the force signal that produces a lower output error injection gain from V to SOC in the regions where the voltage signal has flat slope. It can be appreciated from Figure 6g that the bias term is able to estimate the drift value slowly in our force measurement. The slow convergence of the bias estimate to the 3 N value denoted by the dotted black line can be seen in Figure 6g. Therefore, as time progresses, the estimated force converges to the modeled data, as shown in Figure 6c. This is due to the chosen inflection points in F–SOC function being the same as the modeled data (ˆ *bL* = *bL* = *c*<sup>3</sup> and ˆ *bH* = *bH* = *c*5). We can observe that the estimated bias value oscillates around the true bias value. The model error is being attributed to the bias state by the estimation algorithm. The Switch Observer V&F Bias has lower root mean square error (RMSE), faster time convergence to the denoted SOC EEB (referred to as Time to 5% EEB in Table 3), and reduced maximum absolute SOC error after the force measurement has converged to an SOC error (referred to as Max SOC Error in Table 3) compared to Observer V as shown in Table 3. Therefore, the advantage of Switch Observer V&F Bias is the fast convergence in the region of 40–60% SOC while having more accurate SOC estimation due to the low RMSE values.

**Table 3.** Comparison of the RMSE index for different initial estimate error and different sensors using the simulated data. V, terminal voltage; V&F Bias, fusion of both sensors with bias state.


#### **7. Simulation Results with Inflection Point Mismatch in F–SOC and OCV–SOC**

In reality, we do not have a "perfect" model that captures the battery data. Moreover, during aging, the capacity loss shifts the inflection points as compared with the fresh cell. Therefore, to simulate a more realistic application model mismatch is included in the F–SOC and OCV–SOC observer by shifting the inflections points *bL* and *bH* of by <sup>−</sup>10% (ˆ *bL* and ˆ *bH*) to represent capacity loss on the negative electrode [15]. During capacity loss, these inflections points shift for both functions of F–SOC and OCV–SOC because they correspond to the electrochemical and mechanical model having the same phase transitions. The aging effect on the swelling as the battery is cycled is that the inflection points will shift by approximately 10%. The proposed V&F bias SOC method will work if the initial unknown bias is within 3 N. The bias estimator will converge to the true value, when the SOC is outside the middle region (where there is a multiplicity of state of charge). A larger error in the initial bias will be "corrected" by visiting 100% or 0% SOC based on the cell voltage feedback. In the middle region (*RI I*), the feedback of force error is split between the SOC and bias in the middle region without a strong feedback from the terminal voltage and therefore can have a persistent SOC error. Checking the observer performance under inflection points mismatch due to capacity loss is important since it captures the robustness needed as the battery ages. Therefore, we want to check the performance of the developed Switching Observer V&F Bias under this model mismatch and both voltage and force have −10% inflections points shift. We initialize our simulated "measured" data at 61% SOC. The weights chosen to obtain the gains *K*1*<sup>i</sup>* with *i* = 1–4 and *K*2*<sup>i</sup>* with *i* = 1–4 are shown in Table 4. The gain *K*<sup>32</sup> for feedback of the force error to the SOC state is shown in Figure 7f. The results from the "perfect" model observer (Switch V&F Bias), as shown in Figure 6, are compared with the observer with inflection mismatch (Switch V&F Bias Mismatch).

**Figure 7.** Simulation of the impact of model mismatch in the inflection points of the force vs. SOC curve. A 10% error in ˆ *bL* and ˆ *bH* is tested for the observer (Switch Observer V&F Bias) with an emulated bias of 3 N and 10% (0.1) initial SOC error. (**a**) Comparison of the simulated state of charge (SOC). (**b**) Comparison of the state of charge (SOC) error with the dashed lines representing the target ±5% bound. The SOC error converges to within the ±5% error bound. The largest SOC errors are observed near the switching points denoted by the horizontal dashed lines in Figure 7a. (**c**) Comparison of the simulated force and observer output. For generalization of the results to the other cell sizes, the pressure is shown. (**d**) Comparison of the simulated terminal voltage and observer output. (**e**) Comparison of the true slope with the estimated slope *dF*˜ *dz* and observer output. (**f**) Comparison of the feedback gain *K*<sup>32</sup> from force to SOC based on the switching logic. The gain is set to zero when there is mismatch in the estimated slope and that based on the observer SOC. (**g**) The Bias state estimate converges slowly. The impact of model mismatch in the inflection points of the force vs. SOC curve can be seen by comparing the (Switch Observer V&F Bias Mismatch) and (Switch Observer V&F Bias) which uses the correct value. The accurate observer converges to the true value of 3 N, denoted by the dashed horizontal line whereas mismatch leads to a constant over-estimate of about 1 N.

**Table 4.** Estimation Weights for the different sensors using the simulated data with and without model mismatch. V&F Bias fusion of both sensors without mismatch; V&F Bias Mismatch, fusion of both sensors with model mismatch.


The SOC estimation using the V measurement and V&F Bias measurement with 10% initial estimation error are shown in Figure 7a. The inflection points in the force with respect to SOC are denoted by the dotted horizontal lines in Figure 7a. As shown in Table 5, we obtain lower RMSE, faster convergence to the desired SOC estimation error bound (EEB), and smaller absolute SOC error after the force measurement has converged to an SOC error value for the V&F Bias observer compared to V&F Bias Mismatch observer when initialized at 10%. For −10% SOC initialization, the V&F Bias Mismatch observer has lower RMSE and lower time convergence to the desired SOC estimation error bound (EEB) compared to V&F Bias, as shown in Figure 7b. The force bias estimate for V&F Bias Mismatch converges faster to the constant bias value of 3N. Due to this faster convergence, the RMSE and the convergence to the desired SOC estimation error bound (EEB) is lower for the model with mismatch due to the initial condition. In both −10% and 10% SOC initialization, the maximum absolute SOC error after the force measurement has converged to an SOC error lower for V&F Bias than V&F Bias Mismatch. The non-zero state estimation error is due to the mismatch present in the F–SOC and OCV–SOC function (*bL* = *c*<sup>3</sup> and *bH* = *c*5). To understand this, we need to analyze the error dynamics equation. We know from Section 3 the form of our model and the observer form is given in Equation (12). Therefore, denoting our error as *e* = *x* − *x*ˆ and using our model and observer model equations, we obtain the dynamic error equation as

$$
\dot{\varepsilon} = (A - K\mathbb{C})\varepsilon - K\Delta\mathbb{C}\hat{\varepsilon} - K\Delta\mathbb{C}\_0 \tag{27}
$$

where <sup>Δ</sup>*<sup>C</sup>* <sup>=</sup> *<sup>C</sup>* <sup>−</sup> *<sup>C</sup>*<sup>ˆ</sup> and <sup>Δ</sup>*C*<sup>0</sup> <sup>=</sup> *<sup>C</sup>*<sup>0</sup> <sup>−</sup> *<sup>C</sup>*<sup>ˆ</sup> 0. From the dynamic error equation, we notice that the bias error and the SOC error will not converge to 0 due to the terms −*K*Δ*C*0. Therefore, the bias will converge to a value that is not the true value of the drift due to this error, as shown in Figure 7g. The SOC error will converge to a value but it will not converge to 0, as shown in Figure 7b, due to the model mismatch. The magnitude of the estimation error varies depending on the force region we are operating in. The SOC estimation error, due to model mismatch, will grow as the force sensor drift increases due to the terms −*K*Δ*C*0. For the given model tuning and 10% shift in the force curve with respect to SOC, the force-based observer only achieves better SOC estimation than the voltage only case if the uncorrected force sensor drift is less than 3 N initially.

**Table 5.** Comparison of the RMSE index for different initial estimate error and different sensors using the simulated data with and without model mismatch. V, terminal voltage; V&F Bias, fusion of both sensors; V&F Bias Mismatch, fusion of both sensors with model mismatch.


#### **8. Experimental Data Results Using F–SOC and OCV–SOC**

With a 10% initial error in our SOC estimate, we obtain the results shown in Figure 8. The observer was initialized to 51% SOC, whereas the true state was 61% SOC, to highlight the performance in the middle SOC region where voltage based techniques are less effective. The current waveform discussed in Section 2 was applied to the battery and the observer. As in the simulated results, the proposed Switch Observer V&F Bias has faster convergence than the SOC estimation based on V alone as shown in Figure 8b. The V&F Bias observer exhibits lower RMSE, faster time convergence to the denoted SOC estimation error bound (EEB), and reduced maximum absolute SOC error, as shown in Table 6 compared to V observer. The bias term oscillates around the estimated bias value of 6.8 N, as shown in Figure 8g. These larger oscillations in the bias could be due to our force data having additional dynamics besides the bias term. According to [13], the force has a dynamic term that is temperature dependent, and the ambient chamber temperature may oscillate within ±1 ◦C. Therefore, this dynamic term should be added to our force model. There are some large errors in SOC estimation (approximately around 15%) at low SOC. This is due to our piecewise linear (PWL) fit.

**Initial SOC Error Parameters Observer V Switch Observer V&F** Time to 5% EEB [min] 95.5 0.18 +10% Max SOC Error [%] 11.89 11.79 RMSE 0.0741 0.0611 Time to 5% EEB [min] 6.95 6.54

−10% Max SOC Error [%] 11.89 11.79

RMSE 0.0594 0.0574

**Table 6.** Comparison of the RMSE index for different initial estimate error and different sensors for data validation. V, terminal voltage; V&F Bias, fusion of both sensors.

In the area near 20% SOC, the PWL fit is less accurately, as shown in Figure 1, which results in increased SOC estimation error of around 10%. To better capture the non-linearity, the piecewise linear approximation could be further divided into more regions to provide a better fit. The weights chosen to obtain the gains *K*1*<sup>i</sup>* with *i* = 1–4 and *K*2*<sup>i</sup>* with *i* = 1–4 are shown in Table 7. The gain *K*<sup>32</sup> that satisfies sign or is zero is shown in Figure 8f. The bias estimation state is initialized at 10 N.

**Table 7.** Control Weights for the different sensors for data validation: V,terminal voltage; V&F Bias, fusion of both sensors.


**Figure 8.** Experimental validation of the developed observer (Switch Observer V&F Bias) with 10% (0.1) initial SOC error. The offline experimentally measured current, voltage, and force data were fed into the model to assess performance. (**a**) Comparison of high accuracy coulomb counting based state of charge (SOC) with observer estimates. (**b**) Comparison of the state of charge (SOC) error with the dashed lines representing the target ±5% bound. The observer with force bias state estimation demonstrates better performance. The large error at low SOC for both observers is due to model mismatch the piecewise linear OCV–SOC fit at low SOC. (**c**) Comparison of measured and estimated force. (**d**) Comparison of the experimental and observer modeled terminal voltage. (**e**) Comparison of the estimated slope *dF*˜ *dz* and observer model. (**f**) Comparison of the feedback gain *K*<sup>32</sup> from force to SOC based on the switching logic. (**g**) The Bias state estimate fluctuates around the average value (denoted by the horizontal dashed line). This could be due model mismatch, where the force error is exciting the bias state estimate.

#### **9. Conclusions**

In this paper, a new switching estimator design for a battery with the lithium ion LFP chemistry that integrates the non-monotonic F–SOC relation is proposed, verified by simulation, and validated using experimental data with respect to the SOC estimation accuracy. The estimator is based on switching PWL models that are scheduled according to the identified slope of the F–SOC operating point with a bias state in order to capture the drift exerted by the fixture and battery in our force measurement. Two different sensor scenarios, namely V and V&F Bias fusion, are compared, where it is concluded that the V&F Bias sensor improves the rate of SOC estimator convergence. This is

due to the information added from the steeper, hence more informative, F–SOC characteristic than OCV–SOC relation despite large errors in the voltage and force models. The bounds on SOC estimation accuracy depend on the chosen inflection points of the F–SOC function. If they are correct, then the accuracy is better for the observer with both F and V than V only. Our future work will focus on determining if the drift on our force measurement is due to creep exerted by the poron and thermal expansion of the fixture or due to creep exerted by a degraded battery influenced by compressive stresses [32], or a combination of both. The thermal expansion term and swelling dynamic as a function of temperature will be considered. Data from an aged swelling cell will be used with the developed SOC estimator and the State of Health of the battery (SOH) will be estimated through the capacity fade. The inflection point model mismatch on the F–SOC function will be further studied in an aged cell since the inflection points change as the battery degrades or fades. Due to the tight manufacturing tolerances, the variability in thickness should not be a significant contributor to force measurement uncertainty. In terms of aging, the expected variability in the cell expansion is a subject of future studies. Cell-to-cell variability due to aging in the resulting force should also be investigated in future work based on initial findings from [14,15]. Finally, the measured force is the result of all cell's expansion (summation) in a constrained module, therefore different levels of degradation for individual cells when charged in series would result in smoothing (convolution) of the force sensed on the module-level.

#### **10. Disclaimer**

Reference herein to any specific commercial company, product, process, or service by trade name, trademark, manufacturer, or otherwise, does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or the Department of the Army (DA). The opinions of the authors expressed herein do not necessarily state or reflect those of the United States Government or the DA, and shall not be used for advertising or product endorsement purposes.

**Author Contributions:** Methodology, Software Simulation, Experimental validation, Visualization, and Writing—original draft preparation, M.A.F.-S.; Conceptualization, Methodology, Supervision, Funding acquisition, and Writing—review and editing, J.B.S.; and Conceptualization, Supervision, Funding acquisition, and Project administration, A.G.S. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the National Science Foundation (NSF) Graduate Research Fellowship Program and the DOE's Advanced Research Projects Agency-Energy (ARPA-E) AMPED program under award number DE-AR0000269.

**Acknowledgments:** The authors would like to acknowledge funding from the National Science Foundation (NSF), specifically the Graduate Research Fellowship Program, the technical and financial support of the Automotive Research Center (ARC) in accordance with Cooperative Agreement W56HZV-14-2-0001, and the U.S. Army Combat Capabilities Development Command (CCDC) Ground Vehicle Systems Center (GVSC) in Warren, MI.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Abbreviations**

The following abbreviations are used in this manuscript:



#### **Appendix A. Voltage Hysteresis State H Function**

$$H(z, \operatorname{sign}(I)) = \operatorname{sign}(I)(a\_{12}z^{12} + a\_{11}z^{11} + a\_{10}z^{10} + a\_{9}z^{9} + a\_{8}z^{8} + a\gamma z^{7} + a\varsigma z^{6} + a\varsigma z^{5})$$

$$+ a\_{4}z^{4} + a\_{3}z^{3} + a\_{2}z^{2} + a\_{1}z + a\_{0})\tag{A1}$$

where *ai* with *i* = 0–12 are the tuned parameters and their values are found in Table A1.


**Table A1.** Voltage hysteresis state H function parameters and its values.

#### **Appendix B. Discrete Terminal Voltage Model**

The discrete terminal voltage measurement equation is defined as

$$V\_{T,t} = \mathcal{V}\_{\alpha\varepsilon}(z\_t) - IR - V\_{1,t} - V\_{2,t} + \upsilon\_{V,t} \tag{A2}$$

where *vV*,*<sup>t</sup>* is the V measurement noise. The piecewise linear (PWL) approximation of the OCV characteristic is modeled as

$$\mathcal{V}\_{\rm oc}(z\_{t}) = \begin{cases} \begin{array}{ll} \zeta z\_{t} + \zeta\_{0} & \text{if } z\_{t} \leq c\_{0} \\ \eta z\_{t} + \eta\_{0} & \text{if } c\_{0} < z\_{t} \leq c\_{1} \\ \theta z\_{t} + \theta\_{0} & \text{if } c\_{1} < z\_{t} \leq c\_{2} \\ \kappa z\_{t} + \kappa\_{0} & \text{if } c\_{2} < z\_{t} \leq c\_{3} \\ \sigma z\_{t} + \sigma\_{0} & \text{if } c\_{3} < z\_{t} \leq c\_{4} \\ \mu z\_{t} + \mu\_{0} & \text{if } c\_{4} < z\_{t} \leq c\_{5} \\ \varrho z\_{t} + \varphi\_{0} & \text{if } c\_{5} < z\_{t} \geq c\_{6} \\ \lambda z\_{t} + \lambda\_{0} & \text{if } z\_{t} \geq c\_{6} \end{array} \tag{A3}$$

where *ζ*, *η*, *θ*, *κ*, *σ*, *μ*, *ϕ*, and *λ* are the slope parameters; *ζ*<sup>0</sup> is the minimum voltage sensed at fully discharged state; and *c*0, *c*1, *c*2, *c*3, *c*4, *c*5, and *c*<sup>6</sup> are the piecewise point parameters. The parameters

*η*0, *θ*0, *κ*0, *σ*0, *μ*0, and *λ*<sup>0</sup> are uniquely determined from the other parameters via constraints of piecewise continuity.

$$
\eta\_0 = (\mathbb{\zeta} - \eta)c\_0 + \mathbb{\zeta}\_0 \tag{A4}
$$

$$
\theta\_0 = (\eta - \theta)\mathbf{c}\_1 + (\mathbb{J} - \eta)\mathbf{c}\_0 + \mathbb{J}\_0 \tag{A5}
$$

$$\kappa\_0 = (\theta - \kappa)c\_2 + (\eta - \theta)c\_1 + (\zeta - \eta)c\_0 + \zeta\_0 \tag{A6}$$

$$\boldsymbol{\sigma}\_{0} = (\boldsymbol{\kappa} - \boldsymbol{\sigma})\boldsymbol{c}\_{3} + (\boldsymbol{\theta} - \boldsymbol{\kappa})\boldsymbol{c}\_{2} + (\boldsymbol{\eta} - \boldsymbol{\theta})\boldsymbol{c}\_{1} + (\boldsymbol{\zeta} - \boldsymbol{\eta})\boldsymbol{c}\_{0} + \boldsymbol{\zeta}\_{0} \tag{A7}$$

$$\mu\_0 = (\sigma - \mu)c\_4 + (\kappa - \sigma)c\_3 + (\theta - \kappa)c\_2 + (\eta - \theta)c\_1 + (\zeta - \eta)c\_0 + \zeta\_0 \tag{A8}$$

$$\varphi\_0 = (\mu - \varrho)c\mathfrak{z} + (\sigma - \mu)c\mathfrak{z} + (\kappa - \sigma)c\mathfrak{z} + (\theta - \kappa)c\mathfrak{z} + (\eta - \theta)c\mathfrak{z} + (\zeta - \eta)c\mathfrak{z} + \zeta\mathfrak{o} \tag{A9}$$

$$\lambda\_0 = (\varrho - \lambda)c\_6 + (\mu - \varrho)c\_5 + (\sigma - \mu)c\_4 + (\kappa - \sigma)c\_3 + (\theta - \kappa)c\_2 + (\eta - \theta)c\_1 + (\zeta - \eta)c\_0 + \zeta\_0 \tag{A10}$$

The OCV characteristic and its PWL approximation are shown in Figure 3. Note that *c*<sup>3</sup> and *c*<sup>5</sup> also correspond to the inflection points in the PWL approximation of the force. The reason for this is that both functions of F–SOC and OCV–SOC have the same inflection points due to the electrochemical and mechanical model having the same phase transitions.

#### **Appendix C. Force and OCV Function Values that Represent Average Data**

**Table A2.** Parameter values for the F–SOC function that represent average data. This model is also used for PWL F–SOC without model mismatch (*bL* and *bH*) and with model mismatch (ˆ *bL* and ˆ *bH*).



**Table A3.** Parameter values for the OCV–SOC function in Equation (7).

**Table A4.** Parameter values for the PWL OCV–SOC function in Equation (A3) without model mismatch. Same values are used with model mismatch except for *c*<sup>3</sup> = 0.34 and *c*<sup>5</sup> = 0.69.


#### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
