**3. Battery State of Charge Estimation based on Multi-Innovation and Kalman Filter Framework Algorithm**

#### *3.1. Definition and Characteristics of a Battery's State of Charge*

SOC refers to the ratio of the remaining charge margin in the battery to the rated charge capacity of the battery. The calculation formula is as follows.

$$\text{SOC} = \text{SOC}\_0 - \frac{1}{\overline{Q}\_N} \int\_{\tau = t\_0}^{\tau = t} \mu I d\tau \tag{5}$$

where *SOC*<sup>0</sup> is the initial SOC value, *QN* is the maximum discharge capacity, *μ* represents coulomb efficiency and generally is set to 1, and *I* is the load current, and its discharge is specified as positive.

*SOC* can only be estimated through indirect methods with, e.g., voltage and current as inputs, which is related to resistance capacitance parameters and working conditions, and so the estimation is related to the ohmic internal resistance, polarization internal resistance, the temperature and health state of the battery, etc. The system has highly complex nonlinear characteristics, and so the estimation error is large if only using the measurement method. If the data of the two dimensions—system state estimation and measurement results—can be used at the same time, the state of the system model over time can be used as an a priori estimation value, the measured value is updated as an a posteriori estimate, and the algorithm should have a very efficient performance improvement in theory. A Bayesian filter has such characteristics. A Bayesian filter can be derived from the Bayesian full probability formula without any omission of the nonlinear characteristics of the system, but it is difficult to obtain the prior probability distribution of the system in practical use, and the calculation is also too complex. In order to facilitate calculation, in practical application, it is generally assumed that the probability distribution is normally distributed, and the system is a complete information system—that is, the state of the next time is only related to the current time. In this case, the algorithm is simplified to the KF algorithm.

KF and its extended algorithm adopt the recursive filtering algorithm, which has the characteristics of simple calculation. At the same time, the algorithm uses the data of two relatively independent dimensions—state equation and measurement equation—for verification and updating, which can effectively improve the accuracy of the algorithm. Therefore, it is widely used in SOC estimation.

*3.2. Kalman Filter*

The basic formula of KF is as follows:

$$
\omega\_{k+1} = A\_{k+1}\omega\_k + B\_{k+1}u\_{k+1} + \omega\_{k+1} \tag{6}
$$

$$y\_{k+1} = H\_{k+1} x\_{k+1} + \upsilon\_{k+1} \tag{7}$$

where Equation (6) is called the state equation, and Equation (7) is called the measurement equation. *xk* is the state vector at time *k*, *xk*<sup>+</sup><sup>1</sup> is the state vector at time *k +* 1, *yk*<sup>+</sup><sup>1</sup> is the observation vector at time *k* + 1, *uk*<sup>+</sup><sup>1</sup> is the system input vector at time *k* + 1, *A* is the state transition parameter matrix, *B* is the input control parameter matrix, *H* is the observation parameter matrix, *ωk*+<sup>1</sup> ∼ *N*(0, *Qk*<sup>+</sup>1) is the system noise vector, and *υk*+<sup>1</sup> ∼ *N*(0, *Rk*<sup>+</sup>1) is the measurement noise vector. *ωk*+<sup>1</sup> and *υk*+<sup>1</sup> are independent of each other.

The recursive process of the algorithm is as follows.

(1) Predict the system state at the next time.

$$\mathfrak{A}\_{k+1}^{(-)} = A\_{k+1} \mathfrak{A}\_{k}^{(+)} + B\_{k+1} u\_{k+1} \tag{8}$$

(2) Predict the system covariance at the next time.

$$P\_{k+1}^{(-)} = A\_{k+1} P\_k^{(+)} A\_{k+1} \, ^T + Q\_{k+1} \tag{9}$$

(3) Calculate the Kalman gain of the measurement update.

$$\mathbb{K}\_{k+1} = \frac{P\_{k+1}^{(-)} H\_{k+1}{}^T}{H\_{k+1} P\_{k+1}^{(-)} H\_{k+1}{}^T + R\_{k+1}} \tag{10}$$

(4) Update the system status by the measured values.

$$\mathfrak{X}\_{k+1}^{(+)} = \mathfrak{X}\_{k+1}^{(-)} + K\_{k+1} (y\_{k+1} - H\_{k+1} \mathfrak{X}\_{k+1}^{(-)}) \tag{11}$$

(5) Update system covariance by the measured values.

$$P\_{k+1}^{(+)} = (I - K\_{k+1} H\_{k+1}) P\_{k+1}^{(+)} \tag{12}$$

where *x*ˆ (−) *<sup>k</sup>*+<sup>1</sup> represents an a priori estimate of x at time *k* + 1, *x*ˆ (+) *<sup>k</sup>*+<sup>1</sup> represents an a posteriori estimate of *x* at time *k* + 1, *x*ˆ (+) *<sup>k</sup>* represents an a posteriori estimate of *x* at time *k*, *<sup>P</sup>*(−) *<sup>k</sup>*+<sup>1</sup> represents an a priori estimate of *<sup>P</sup>* at time *<sup>k</sup>* + 1, *<sup>P</sup>*(+) *<sup>k</sup>*+<sup>1</sup> represents an a posteriori estimate of *<sup>P</sup>* at time *<sup>k</sup>* + 1, and *P*(+) *<sup>k</sup>* represents an a posteriori estimate of *P* at time *k. Kk*<sup>+</sup><sup>1</sup> represents the Kalman gain at time *k* + 1.

KF requires that the state equation and measurement equation of the system are linear, which is difficult to apply to nonlinear tasks such as the battery state of charge. Instead, the EKF and UKF algorithms are generally used in such tasks.

#### *3.3. Extended Kalman Filter Algorithm*

The EKF algorithm is transformed into the KF algorithm by expanding the Taylor series of state and measurement functions and ignoring the second-order and higher-order terms. The basic formula is as follows:

$$x\_{k+1} = f(x\_{k'} u\_{k'} k, \omega\_k) \tag{13}$$

$$y\_{k+1} = h(\mathbf{x}\_{k+1}, k) + \upsilon\_{k+1} \tag{14}$$

where *uk* is the system input vector at time *k*, *ω<sup>k</sup>* ∼ *N*(0, *Qk*) is the system noise vector, and *υ<sup>k</sup>* ∼ *N*(0, *Rk*) is the measurement noise vector. *ω<sup>k</sup>* and *υ<sup>k</sup>* are independent of each other. The recursive process of the algorithm is as follows.

(1) Predict the system state at the next time.

$$\mathfrak{X}\_{k+1}^{(-)} = f(\mathfrak{X}\_{k\prime}u\_k, k) \tag{15}$$

(2) Predict the system covariance at the next time.

Let

$$F\_k = \frac{\partial f}{\partial \mathcal{X}}\Big|\_{\mathbf{x}\_k = \mathbf{t}\_k^{(+)}}\tag{16}$$

Get

$$P\_{k+1}^{(-)} = F\_k P\_k^{(+)} F\_k^{\;T} + Q\_k \tag{17}$$

(3) Calculate the Kalman gain of the measurement update. Let

$$H\_{k+1} = \frac{\partial h}{\partial \mathbf{x}}\Big|\_{\mathbf{x}\_{k+1} = \mathbf{x}\_{k+1}^{(-)}}\tag{18}$$

Get

$$K\_{k+1} = \frac{P\_{k+1}^{(-)} H\_{k+1}{}^T}{H\_{k+1} P\_{k+1}^{(-)} H\_{k+1}{}^T + R\_{k+1}} \tag{19}$$

(4) Update the system status by the measured values.

$$\mathfrak{X}\_{k+1}^{(+)} = \mathfrak{X}\_{k+1}^{(-)} + K\_{k+1} \left\{ y\_{k+1} - h \left[ \mathfrak{X}\_{k+1'}^{(-)}, k+1 \right] \right\} \tag{20}$$

(5) Update system covariance by the measured values.

$$P\_{k+1}^{(+)} = P\_{k+1}^{(-)} - K\_{k+1} \left[ H\_{k+1} P\_{k+1}^{(-)} H\_{k+1}^T + R\_{k+1} \right] \tag{21}$$

where *x*ˆ (−) *<sup>k</sup>*+<sup>1</sup> represents an a priori estimate, and *x*ˆ (+) *<sup>k</sup>*+<sup>1</sup> represents an a posteriori estimate of *x* at time *k+1*; *x*ˆ (+) *<sup>k</sup>* represents an a posteriori estimate of *<sup>x</sup>* at time *<sup>k</sup>*; *<sup>P</sup>*(−) *<sup>k</sup>*+<sup>1</sup> represents an a priori estimate of *P* at time *k+1*; *P*(+) *<sup>k</sup>*+<sup>1</sup> represents an a posteriori estimate of *P* at time *k +* 1, and *P*(+) *<sup>k</sup>* represents an a posteriori estimate of *P* at time *k. Kk*<sup>+</sup><sup>1</sup> represents the Kalman gain at time *k +* 1.

The EKF algorithm omits the influence of second-order and higher-order terms and is only applicable to the case of weak nonlinearity. For the case of heavy nonlinearity, the error is large, and more complex nonlinear filtering algorithms are often required, such as the UKF.

#### *3.4. Unscented Kalman Filter Algorithm*

The UKF takes the KF as the basic framework—the basic formula is consistent with that of the EKF, and unscented transformation is applied to realize nonlinear application scenarios. The recursive algorithm flow is as follows.

(1) Let the a posteriori state estimation and covariance at time *k* be *x* (+) *<sup>k</sup>* and *<sup>P</sup>*(+) *<sup>k</sup>* , respectively. (2) Calculate sampling points.

$$\begin{cases} \begin{array}{c} \chi\_k^0 = \mathfrak{x}\_k^{(+)} \\ \chi\_k^i = \mathfrak{x}\_k^{(+)} + \sqrt{(L+\eta)P\_{\text{xx}}}i = 1,2,\ldots,L \\ \chi\_k^i = \mathfrak{x}\_k^{(+)} - \sqrt{(L+\eta)P\_{\text{xx}}}i = L+1,L+2,\ldots,2L \end{array} \end{cases} \tag{22}$$

where *L* is the length of the state vector, and the weight value is calculated as follows:

$$\begin{cases} \begin{array}{c} \eta = a^2(L+k\_i) - L \\ W\_m^0 = \frac{\eta}{L+\eta}, W\_m^i = \frac{1}{2(L+\eta)}, i = 1, 2 \dots 2L \\ W\_c^0 = \frac{\eta}{L+\eta} + 1 - a^2 + \beta, W\_c^i = \frac{1}{2(L+\eta)}, i = 1, 2 \dots 2L \end{array} \end{cases} \tag{23}$$

where subscript *c* represents the weight of covariance; subscript *m* represents the weight of mean square deviation; *η* represents the scaling ratio, and *α* represents the distribution state of sampling points—when *α* is large, it indicates a greater weight of sigma points at the average value. *β* is a weight, which is used to combine the dynamic differences of higher-order terms in the equation. In this paper, we set *L* = 3, *α* = 0.01, *ki* = 0, *β* = 2.

(1) Update the a priori state value *x* (−) *<sup>k</sup>*+<sup>1</sup> and a priori variance value *<sup>P</sup>*(−) *<sup>k</sup>*+1.

$$\chi\_{k+1}^{i} = f(\chi\_{k'}^{i}, u\_{k'}k) \tag{24}$$

$$\chi\_{k+1}^{(-)} = \sum\_{i=0}^{2L} \mathcal{W}\_m^i \chi\_{k+1}^i \tag{25}$$

$$P\_{k+1}^{(-)} = \sum\_{i=0}^{2L} \left( W\_c^i (\chi\_k^i - \mathfrak{x}\_{k+1}^{(-)}) (\chi\_k^i - \mathfrak{x}\_{k+1}^{(-)})^T \right) + Q\_{k+1} \tag{26}$$

where *Qk* is the system noise covariance matrix; in this paper, we set *Qk*<sup>=</sup> 0.00000001 <sup>∗</sup> <sup>⎡</sup> 100 ⎤


$$y\_{k+1}^i = h(\chi\_{k+1}^i, k+1) \tag{27}$$

$$\mathcal{Y}\_{k+1} = \sum\_{i=0}^{2L} \mathcal{W}\_m^i \mathcal{Y}\_{k+1}^i \tag{28}$$

$$P\_{k+1}^{yy} = \sum\_{i=0}^{2L} \left( \mathcal{W}\_c^i (\mathcal{Y}\_{k+1}^i - \mathcal{Y}\_{k+1}) \left( \mathcal{Y}\_{k+1}^i - \mathcal{Y}\_{k+1} \right)^T + R\_{k+1} \tag{29}$$

$$P\_{k+1}^{xy} = \sum\_{i=0}^{2L} \mathcal{W}\_c^i (\chi\_{k+1}^i - \mathfrak{x}\_{k+1}^{(-)}) (y\_{k+1}^i - \mathfrak{y}\_{k+1})^T \tag{30}$$

where *Rk* is the measurement noise covariance matrix; in this paper, we set *Rk* = 1. (3) Update the a posteriori state value *x* (+) *<sup>k</sup>*+<sup>1</sup> and a posteriori state error covariance *<sup>P</sup>*(+) *k* using the measured value *yk*+1.

$$K\_{k+1} = \frac{P\_{k+1}^{xy}}{P\_{k+1}^{yy}}\tag{31}$$

$$\mathbf{x}\_{k+1}^{(+)} = \mathbf{x}\_{k+1}^{(-)} + K\_{k+1} (\mathbf{y}\_{k+1} - \mathbf{\hat{y}}\_{k+1}) \tag{32}$$

$$P\_{k+1}^{(+)} = P\_{k+1}^{(-)} - \mathcal{K}\_{k+1} P\_{k+1}^{yy} \mathcal{K}\_{k+1} \tag{33}$$

where *Kk*<sup>+</sup><sup>1</sup> represents Kalman gain at time *k* + 1.

Compared with the EKF algorithm, the UKF can adapt to systems with stronger nonlinearity and can achieve third-order approximation accuracy in the case of Gaussian distribution and second-order approximation accuracy in the case of non-Gaussian distribution; the UKF does not need to calculate the Jacobian matrix, but the number of sampling points is 2*n* + 1, and the overall amount of calculation is larger than the EKF algorithm.
