# SOC Estimator (Adaptive Kalman Filter, Variable Capacity)

State of charge and terminal resistance estimator with adaptive Kalman filter and variable capacity

*Since R2023b*

**Libraries:**

Simscape /
Battery /
BMS /
Estimators

## Description

This block implements an estimator that calculates the state of charge (SOC) and terminal
resistance of a battery by using the Kalman filter algorithms. The terminal resistance,
*R*_{0}, is a state of the estimator. The cell
capacity of the battery is an input to the block.

The SOC is the ratio of the released capacity
*C*_{releasable} to the rated capacity
*C*_{rated}. Manufacturers provide the value of
the rated capacity of each battery, which represents the maximum amount of charge in the battery:

$$SOC=\frac{{C}_{\text{releasable}}}{{C}_{\text{rated}}}.$$

This block supports single-precision and double-precision floating-point simulation.

**Note**

To enable single-precision floating-point simulation, the data type of all inputs and
parameters, except for the **Sample time (-1 for inherited)**
parameter, must be `single`

.

For continuous-time simulation, set the **Filter type** parameter to
`Extended Kalman-Bucy filter`

or ```
Unscented
Kalman-Bucy filter
```

.

**Note**

Continuous-time implementation of this block works only in a double-precision floating-point simulation. If you provide single-precision floating-point parameters and inputs, this block casts them to double-precision floating-point values to prevent errors.

For discrete-time simulation, set the **Filter type** parameter to
`Extended Kalman filter`

or ```
Unscented Kalman
filter
```

and the **Sample time (-1 for inherited)**
parameter to a positive value or `-1`

.

### Equations

These figures show the equivalent circuit for a battery with one-time-constant dynamics and two time-constant dynamics, respectively:

The equations for the equivalent circuit with the terminal resistance
*R*_{0} as an additional state and two
time-constant dynamics are:

$$\begin{array}{l}\frac{dSOC}{dt}=-\frac{i}{3600AH}\\ \frac{d{V}_{1}}{dt}=\frac{i}{{C}_{1}\left(SOC,T\right)}-\frac{{V}_{1}}{{R}_{1}\left(SOC,T\right){C}_{1}\left(SOC,T\right)}\\ \frac{d{V}_{2}}{dt}=\frac{i}{{C}_{2}\left(SOC,T\right)}-\frac{{V}_{2}}{{R}_{2}\left(SOC,T\right){C}_{2}\left(SOC,T\right)}\\ \frac{d{R}_{0}}{dt}=0\\ {V}_{\text{t}}={V}_{0}\left(SOC,T\right)-i{R}_{0}-{V}_{1}-{V}_{2}\end{array}$$

where

*SOC*is the state of charge.*i*is the current.*V*_{0}is the no-load voltage.*V*_{t}is the terminal voltage.*AH*is the ampere-hour rating.*R*_{1}is the first polarization resistance.*C*_{1}is the first parallel RC capacitance.*R*_{2}is the second polarization resistance.*C*_{2}is the second parallel RC capacitance.*T*is the temperature.*V*_{1}is the polarization voltage over the first RC network.*V*_{2}is the polarization voltage over the second RC network.

A time constant *τ*_{1} for the
first parallel section relates the first polarization resistance
*R*_{1} and the first parallel RC
capacitance *C*_{1} using the relationship $${C}_{1}={\tau}_{1}/{R}_{1}.$$

A time constant *τ*_{2} for the second
parallel section relates the second polarization resistance
*R*_{2} and the second parallel RC
capacitance *C*_{2} using the relationship $${C}_{2}={\tau}_{2}/{R}_{2}.$$

For the Kalman filter algorithms, the block uses this state and these nonlinear process and observation functions:

$$\begin{array}{l}x={\left[\begin{array}{ccc}SOC& {V}_{1}& {R}_{0}\end{array}\right]}^{T}\\ f(x,i)=\left[\begin{array}{c}-\frac{i}{3600AH}\\ \frac{i}{{C}_{1}(SOC,T)}-\frac{{V}_{1}}{{R}_{1}(SOC,T){C}_{1}(SOC,T)}\\ 0\end{array}\right]\\ h(x,i)={V}_{0}(SOC,T)-i{R}_{0}-{V}_{1}\end{array}$$

If you set the **Charge dynamics** parameter to ```
Two
time-constant dynamics
```

, for the Kalman filter algorithms, the
block uses this state and these process and observation functions:

$$\begin{array}{l}x={\left[\begin{array}{cccc}SOC& {V}_{1}& {V}_{2}& {R}_{0}\end{array}\right]}^{T}\\ f(x,i)=\left[\begin{array}{c}-\frac{i}{3600AH}\\ \frac{i}{{C}_{1}(SOC,T)}-\frac{{V}_{1}}{{R}_{1}(SOC,T){C}_{1}(SOC,T)}\\ \frac{i}{{C}_{2}(SOC,T)}-\frac{{V}_{2}}{{R}_{2}(SOC,T){C}_{2}(SOC,T)}\\ 0\end{array}\right]\\ h(x,i)={V}_{0}(SOC,T)-i{R}_{0}-{V}_{1}-{V}_{2}\end{array}$$

**Extended Kalman Filter**

This diagram shows the structure of the extended Kalman filter (EKF):

The EKF technique relies on a linearization at every time step to approximate the nonlinear system. To linearize the system at every time step, the algorithm computes these Jacobians online:

$$\begin{array}{l}F=\frac{\partial f}{\partial x}\\ H=\frac{\partial h}{\partial x}\end{array}$$

The EKF is a discrete-time algorithm. After the discretization, the Jacobians for the SOC estimation of the battery are:

$$\begin{array}{l}{F}_{d}=\left[\begin{array}{ccc}1& 0& 0\\ 0& {e}^{\frac{-{T}_{\text{S}}}{{R}_{1}{C}_{1}}}& 0\\ 0& 0& 1\end{array}\right]\\ {H}_{d}=\left[\begin{array}{ccc}\frac{\partial {V}_{\text{OC}}}{\partial SOC}& -1& -i\end{array}\right]\end{array}$$

where *T*_{S} is the
sample time and *V*_{OC} is the
open-circuit voltage.

The EKF algorithm comprises these phases:

**Initialization**$$\widehat{x}(0|0)$$— State estimate at time step 0 using measurements at time step 0.

$$\widehat{P}(0|0)$$— State estimation error covariance matrix at time step 0 using measurements at time step 0.

**Prediction**Project the states ahead (

*a priori*):$$\widehat{x}\left(k+1|k\right)=f\left(\widehat{x}\left(k|k\right),i\right).$$

Project the error covariance ahead:

$$\widehat{P}\left(k+1|k\right)={F}_{d}\left(k\right)\widehat{P}\left(k|k\right){F}_{d}^{T}\left(k\right)+Q,$$

where

**Q**is the covariance of the process noise.

**Correction**Compute the Kalman gain:

$$K(k+1)=\widehat{P}(k+1|k){H}_{d}^{T}(k){({H}_{d}(k)\widehat{P}(k+1|k){H}_{d}^{T}(k)+R)}^{-1},$$

where

*R*is the covariance of the measurement noise.Update the estimate with the measurement

*y(k)*(*a posteriori*):$$\widehat{x}\left(k+1|k+1\right)=\widehat{x}\left(k+1|k\right)+K\left(k+1\right)\left({V}_{t}\left(k\right)-h\left(\widehat{x}\left(k|k\right),i\right)\right).$$

Update the error covariance:

$$\widehat{P}\left(k+1|k+1\right)=\left(I-K\left(k+1\right){H}_{d}\right)\widehat{P}\left(k+1|k\right).$$

**Extended Kalman-Bucy Filter**

This diagram shows the structure of the extended Kalman-Bucy filter (EKBF):

The EKBF is the continuous-time variant of the Kalman filter. In continuous time, the prediction and correction steps are coupled. The EKBF algorithm comprises these phases:

**Initialization**$$\widehat{x}({t}_{0})$$— State estimate at time

*t*_{0}.$$\widehat{P}({t}_{0})$$— State estimation error covariance matrix at time

*t*_{0}.

**Prediction-Correction EKBF algorithm**$$\begin{array}{l}K\left(t\right)=P\left(t\right){H}^{T}\left(t\right){R}^{-1}\left(t\right)\\ \frac{d\widehat{x}\left(t\right)}{dt}=f\left(\widehat{x}\left(t\right),i\left(t\right)\right)+K\left(t\right)\left({V}_{\text{t}}\left(t\right)-h\left(\widehat{x}\left(t\right),i\left(t\right)\right)\right)\\ \frac{dP\left(t\right)}{dt}=F\left(t\right)P\left(t\right)+P\left(t\right){F}^{T}\left(t\right)+Q\left(t\right)-K\left(t\right)H\left(t\right)P\left(t\right)\end{array}$$

where:

$$\begin{array}{l}F\left(t\right)=\frac{\partial f}{\partial x}\\ H\left(t\right)=\frac{\partial h}{\partial x}\end{array}$$

**Unscented Kalman Filter**

This diagram shows the structure of the unscented Kalman filter (UKF):

The EKF locally approximates nonlinear functions with the linear equations obtained from the Taylor expansion by using only the first term of the expansion. In a highly nonlinear system, these solutions are not very accurate.

The UKF uses nonlinear transformations on a set of sigma points that the
algorithm chooses deterministically. This technique is called
*unscented transformation*. The mean and the covariance
matrix of the transformed points are accurate to the second order of the Taylor
series expansion.

The UKF algorithm follows these steps:

**Initialization**$$\widehat{x}(0|0)$$— State estimate at time step 0 using measurements at time step 0.

$$\widehat{P}(0|0)$$— State estimation error covariance matrix at time step 0 using measurements at time step 0.

Generate sigma points and calculate the mean weight and covariance weight for each point.

Choose the sigma points,

*x*^{(i)}(k|k)$$\begin{array}{l}{x}^{(i)}\left(k|k\right)=\{\begin{array}{cc}\widehat{x}\left(k+1|k\right)& i=1\\ \widehat{x}\left(k+1|k\right)+{\left(\sqrt{\left(n+\lambda \right)P\left(k|k\right)}\right)}_{i}& i=2,\dots ,n+1\\ \widehat{x}\left(k+1|k\right)-{\left(\sqrt{\left(n+\lambda \right)P\left(k|k\right)}\right)}_{i}& i=n+2,\dots ,2n+1\end{array}\\ {W}_{m}^{(i)}=\{\begin{array}{cc}\frac{\lambda}{n+\lambda}& i=1\\ \frac{1}{2\left(n+\lambda \right)}& i\ne 1\end{array}\\ {W}_{c}^{(i)}=\{\begin{array}{cc}\frac{\lambda}{n+\lambda}+\left(1-{\alpha}^{2}+\beta \right)& i=1\\ \frac{1}{2\left(n+\lambda \right)}& i\ne 1\end{array}\end{array}$$

where:

*n*is the dimension of the state vector*x*.$$\lambda ={\alpha}^{2}\left(n+\kappa \right)-n,\text{\hspace{0.17em}}\alpha \in \left[0,1\right]$$ describes the distance between the sigma point and the mean point. In a normal distribution,

*β =*2 and*κ =*0.$${\left(\sqrt{\left(n+\lambda \right)P}\right)}_{i}$$ is the

*i*th row or column of $$\sqrt{cP}$$. The block calculates the matrix square root by using numerically efficient and stable methods such as the Cholesky decomposition.

Perform first estimation of the system state matrix:

$$\begin{array}{l}{\widehat{x}}^{(i)}\left(k+1|k\right)=f\left({\widehat{x}}^{(i)}\left(k|k\right),i\left(k\right)\right)\\ \widehat{x}\left(k+1|k\right)={\displaystyle \sum _{i=1}^{2n+1}{W}_{m}^{\left(i\right)}}{\widehat{x}}^{\left(i\right)}\left(k+1|k\right)\end{array}$$

Perform first estimation of the covariance matrix of the state variables:

$$P\left(k+1|k\right)={\displaystyle \sum _{i=1}^{2n+1}{W}_{c}^{(i)}\left({\widehat{x}}^{(i)}\left(k+1|k\right)-\widehat{x}\left(k+1|k\right)\right)\cdot {\left({\widehat{x}}^{(i)}\left(k+1|k\right)-\widehat{x}\left(k+1|k\right)\right)}^{T}}+Q$$

Estimate the measured variables:

$$\begin{array}{l}{V}_{t}^{(i)}\left(k+1|k\right)=h\left({\widehat{x}}^{(i)}\left(k+1|k\right),i\left(k\right)\right)\\ {\widehat{V}}_{t}\left(k+1|k\right)={\displaystyle \sum _{i=1}^{2n+1}{W}_{m}^{(i)}{\widehat{V}}_{t}^{(i)}\left(k+1|k\right)}\end{array}$$

Estimate the covariance of the measurement (

*P*) and covariance between the measurement and the state (_{y}*P*):_{xy}$$\begin{array}{l}{P}_{y}={\displaystyle \sum _{i=1}^{2n+1}{W}_{c}^{\left(i\right)}\left({\widehat{V}}_{t}^{\left(i\right)}\left(k+1|k\right)-{\widehat{V}}_{t}\left(k+1|k\right)\right)\cdot {\left({\widehat{V}}_{t}^{\left(i\right)}\left(k+1|k\right)-{\widehat{V}}_{t}\left(k+1|k\right)\right)}^{T}+R}\\ {P}_{xy}={\displaystyle \sum _{i=1}^{2n+1}\left({\widehat{x}}^{\left(i\right)}\left(k+1|k\right)-\widehat{x}\left(k+1|k\right)\right)\cdot {\left({\widehat{V}}_{t}^{\left(i\right)}\left(k+1|k\right)-{\widehat{V}}_{t}\left(k+1|k\right)\right)}^{T}}\end{array}$$

Compute the Kalman filter gain:

$$K\left(k+1\right)={P}_{xy}{P}_{y}^{-1}$$

Perform the second update of the state matrix and of the covariance of the state variables:

$$\begin{array}{l}\widehat{x}\left(k+1|k+1\right)=\widehat{x}\left(k+1|k\right)+K\left(k+1\right)\left({V}_{t}\left(k+1\right)-{\widehat{V}}_{t}\left(k+1|k\right)\right)\\ P\left(k+1|k+1\right)=P\left(k+1|k\right)-K\left(k+1\right){P}_{y}{K}^{T}\left(k+1\right)\end{array}$$

**Unscented Kalman-Bucy Filter**

This diagram shows the structure of the unscented Kalman-Bucy filter (UKBF):

The derived continuous-time filtering equations of the UKBF are similar to the EKBF equations.

Because the UKF uses matrix square roots in its sigma points, the algorithm obtains the square-root version of the UKBF by formulating the filter as a differential equation for the sigma points. The equations for the square-root UKBF are:

$$\begin{array}{l}K\left(t\right)=X\left(t\right)W{h}^{T}\left(X(t),t\right){R}^{-1}(t)\\ M\left(t\right)={A}^{-1}\left(t\right)\left[X\left(t\right)W{f}^{T}\left(X\left(t\right),t\right)+f\left(X\left(t\right),t\right)W{X}^{T}\left(t\right)+Q\left(t\right)-K\left(t\right)R\left(t\right){K}^{T}\left(t\right)\right]{A}^{-T}\left(t\right)\\ \frac{d{X}_{i}\left(t\right)}{dt}=f\left(X\left(t\right),t\right){w}_{m}+K\left(t\right)\left({V}_{t}\left(t\right)-h\left(X\left(t\right),t\right){w}_{m}\right)+\sqrt{c}{\left[\begin{array}{ccc}0& A\left(t\right)\Phi \left(M\left(t\right)\right)& -A\left(t\right)\Phi \left(M\left(t\right)\right)\end{array}\right]}_{i}\end{array}$$

where

$$X\left(t\right)=\left[\begin{array}{ccc}m\left(t\right)& \dots & m\left(t\right)\end{array}\right]+\sqrt{c}\left[\begin{array}{ccc}0& A\left(t\right)& -A\left(t\right)\end{array}\right]$$ is the sigma-point matrix.

$${\Phi}_{ij}\left(M\left(t\right)\right)=\{\begin{array}{cc}{M}_{ij}\left(t\right),& i>j\\ 0.5{M}_{ij}\left(t\right),& i=j\\ 0,& i<j\end{array}$$ is a function that returns the lower diagonal part of the argument.

$${w}_{m}={\left[\begin{array}{ccc}{W}_{m}^{\left(1\right)}& \cdots & {W}_{m}^{\left(2n+1\right)}\end{array}\right]}^{T}$$

$$W=\left(1-\left[\begin{array}{ccc}{w}_{m}& \cdots & {w}_{m}\end{array}\right]\right)\text{diag}\left(\begin{array}{ccc}{W}_{c}^{\left(1\right)}& \cdots & {W}_{c}^{\left(2n+1\right)}\end{array}\right){\left(I-\left[\begin{array}{ccc}{w}_{m}& \cdots & {w}_{m}\end{array}\right]\right)}^{T}$$

$$c={\alpha}^{2}\left(n+\kappa \right)$$

$${\left[\begin{array}{ccc}0& A\left(t\right)\Phi \left(M\left(t\right)\right)& -A\left(t\right)\Phi \left(M\left(t\right)\right)\end{array}\right]}_{i}$$ is the

*i*th column of the argument matrix.

## Assumptions and Limitations

Process and sensor noises are independent, zero-mean, Gaussian noises.

## Ports

### Input

### Output

## Parameters

## Extended Capabilities

## Version History

**Introduced in R2023b**