CN104714188A - 量测噪声方差阵相匹配的电池电荷状态估计方法与系统 - Google Patents

Abstract

本发明为量测噪声方差阵相匹配的电池电荷状态估计方法与系统，本法步骤为：建立电池模型、电压电流采样，最小二乘法参数辨识模型参数；自适应无迹卡尔曼滤波估计SOC；调整滤波增益的方差阵相匹配的鲁棒自适应无迹卡尔曼滤波，计算系统量测噪声实际残差方差阵值和估计的理论残差方差阵，得到ε_k。ε_k大于1自适应调整因子否则调整滤波增益，求得SOC估计值本系统动力电池所接电压、电流传感器经模数转换模块连接微控制器。微控制器含模型参数辨识模块、调整因子计算模块、滤波增益自适应调整、鲁棒自适应无迹卡尔曼滤波SOC估计模块。SOC直接显示或传送到设备CAN网络。本发明动态调整滤波增益，提高了SOC估计的鲁棒性。

Description

量测噪声方差阵相匹配的电池电荷状态估计方法与系统

技术领域

本发明涉及电动汽车动力电池电荷状态估计领域，具体为量测噪声方差阵相匹配的电池电荷状态估计方法与系统。

背景技术

电池是电动汽车整车系统中最重要组件之一，为了确保电动汽车的安全性，合理管理电池至关重要。近年来，锂离子电池以其能量密度高，循环寿命长，无记忆效应，环境友好等优点成为动力电池应用领域研究的热点，并逐渐成为电动车辆动力电池的最佳选择。

电荷状态(State of Charge，SOC)是电池管理系统(Battery ManagementSystem，BMS)中的重要参数，然而它不能直接通过传感器测量所得，只能依据所建立的模型，运用相应的算法，间接估计得到。由于电池运行时受充、放电电流、温度、自放电和寿命等因素的影响，并且在使用过程中表现出的复杂非线性，对精确估计SOC造成一定难度。

目前SOC估计方法主要有：存在累积误差问题的安时法、需长时间静置的开路电压法、卡尔曼滤波法、需大量数据进行学习的模糊逻辑法和神经网络法等。

对于电池内部的化学特性表现出的复杂非线性特点，将经典卡尔曼滤波的状态估计扩展到非线性随机系统的状态估计中，出现了扩展卡尔曼滤波和无迹卡尔曼滤波(Unscented Kalman Filter，UKF)。UKF克服了扩展卡尔曼滤波中需要计算雅克比矩阵，要求状态函数和量测函数必须是连续可微的缺点，但应用扩展卡尔曼滤波是以已知噪声统计特性为先决条件。当噪声统计特性随实际工况条件剧烈变化时，无迹卡尔曼滤波器不具有自适应能力，易出现滤波的状态估计精度下降甚至发散。为了克服这个缺点，又出现了自适应无迹卡尔曼滤波(Adaptive Unscented Kalman Filter，AUKF)。

AUKF虽然同时对量测噪声方差阵R和过程噪声方差阵Q进行调节时，但R与Q的调整未能保证实际残差方差阵与估计得到的理论残差方差阵完全相匹配，并未区别对待系统过程噪声引起的误差和量测噪声引起的误差，其鲁棒性较差。

发明内容

本发明的目的是设计一种量测噪声方差阵相匹配的电池电荷状态估计方法，R本身自适应调整的方法结合实际残差方差阵与估计的理论残差方差阵相匹配的方法，以方差阵的完全匹配为目标，反向推算由虚拟的量测噪声方差阵与量测噪声方差阵R之比构成的调整因子ε_k,实现滤波增益的自适应调整，即方差阵相匹配的鲁棒自适应无迹卡尔曼滤波算法。对滤波增益进行自适应调节，解决过程噪声引起的误差、量测噪声引起的误差所造成的鲁棒性较差的问题，准确估计电池电荷状态SOC，即使在不同电池电荷状态SOC初值下，也能较快地收敛到真值。

本发明的另一目的是设计一种量测噪声方差阵相匹配的电池电荷状态估计系统，其包括微处理器，其含有实现量测噪声方差阵相匹配的电池电荷状态估计方法的对应的软件模块，可嵌入到使用动力电池的装备中，实现动力电池的实时电荷状态在线估计与显示。

本发明设计的量测噪声方差阵相匹配的电池电荷状态估计方法包括如下步骤：

Ⅰ、建立电池模型及参数辨识

Ⅰ-1、建立电池模型

目前，电动汽车电池模型主要有：电化学模型、神经网络模型、阻抗模型、等效电路模型等。本发明根据电池的工作原理，采用联合模型，相对其他的模型，联合模型结构简单、辨识容易，仿真精度较好，具有通用性。其数学表达式：

$x_{k+1}=x_k-\left(\frac{{\mathrm{\η}}_k\mathrm{\Δt}}{C_n}\right)I_k+{\mathrm{\ω}}_k---\left(1\right)$

$U_k=a_0-\frac{a_1}{x_k}-a_2{x}_k+a_3\mathrm{In}\left(x_k\right)+a_4\mathrm{In}(1-x_k)-R_0{I}_k+v_k---\left(2\right)$

其中U_k，x_k和I_k分别是当前k时刻电池端电压，SOC和电流。电池放电过程I_k为正值，充电过程I_k为负值；a₀,a₁,a₂,a₃和a₄是通过最小二乘法得到的相关参数；R₀是欧姆内阻。η_k是k时刻库伦效率，C_n是电池的标称容量，△t表示采样间隔。ω_k是k时刻系统过程噪声,ν_k是k时刻量测噪声。采样得到U_k、I_k。

Ⅰ-2、模型参数辨识

模型参数辨识方法采用最小二乘法得：

c＝[a₀,a₁,a₂,a₃,a₄,R₀]^T＝(M^TM)^-1M^TU    (3)其中,U＝[U₁,U₂,...,U_N]^T是端电压序列，M＝[M₁,M₂,...,M_N]^T是递归矩阵， $M_k=[1,\frac{1}{x_k},x_k,\mathrm{In}\left(x_k\right),\mathrm{In}(1-x_k),I_k].$

通过最小二乘法离线估计得到模型参数a₀,a₁,a₂,a₃，a₄和R₀，以此作为方差阵相匹配的鲁棒自适应无迹卡尔曼滤波(RAUKF)算法的电池模型基础。

Ⅱ、自适应无迹卡尔曼滤波(AUKF)

自适应无迹卡尔曼滤波估计SOC过程如下：

Ⅱ-1、初始状态统计特性为：

${\hat{x}}_0=E\left(x_0\right)---\left(4\right)$

$P_0=\mathrm{Var}\left(x_0\right)=E\left[(x_0-{\hat{x}}_0){(x_0-{\hat{x}}_0)}^T\right]---\left(5\right)$

x₀是SOC初始值，是SOC初始值的平均值，P₀初始误差协方差阵。

Ⅱ-2、SOC的估计值的过程

SOC的估计值

${\hat{y}}_k={\hat{x}}_{k|k-1}+K_k(z_k-{\hat{z}}_{k|k-1})---\left(6\right)$

参数更新

$K_k=P_{{\tilde{x}}_k{\tilde{z}}_k}{\left(P_{{\tilde{z}}_k}\right)}^{-1}---\left(7\right)$

${\hat{z}}_{k|k-1}=\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^m{\mathrm{\χ}}_{i,k|k-1}---\left(8\right)$

是k时刻自协方差阵，是k时刻互协方差阵。K_k是k时刻的滤波增益，z_k是k时刻电池端电压的真实值，是k时刻的SOC值。

其中参数：

① $P_{{\tilde{z}}_k}=\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^c({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T+R_k$

${W_0}^c=\frac{\mathrm{\λ}}{(n+\mathrm{\λ})}+(1-{\mathrm{\α}}^2+\mathrm{\β})$

${W_i}^c=\frac{1}{2(n+\mathrm{\λ})},i=\mathrm{1,2}$

${\hat{x}}_{k|k-1}=\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^m{\mathrm{\γ}}_{i,k|k-1}$

${\mathrm{\γ}}_{i,k|k-1}={\mathrm{\ξ}}_{i,k-1}=\frac{{\mathrm{\η}}_{k-1}\mathrm{\Δt}}{C_n}{I}_{k-1},i=\mathrm{0,1,2}$

$R_k=F_k+\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^c({\mathrm{\χ}}_{i,k|k-1}-z_k+u_k){({\mathrm{\χ}}_{i,k|k-1}-z_k+u_k)}^T$

② ${W_0}^m=\frac{\mathrm{\λ}}{(n+\mathrm{\λ})}$

${W_i}^m=\frac{1}{2(n+\mathrm{\λ})},i=\mathrm{1,2}$

③ ${\mathrm{\χ}}_{i,k|k-1}=a_0-\frac{a_1}{{\mathrm{\ξ}}_{i,k|k-1}}-a_2{\mathrm{\ξ}}_{i,k|k-1}+a_3\mathrm{In}\left({\mathrm{\ξ}}_{i,k|k-1}\right)+a_4\mathrm{In}(1-{\mathrm{\ξ}}_{i,k|k-1})+R_0{I}_k$

④ $P_{{\tilde{x}}_k{\tilde{z}}_k}=\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^c({\mathrm{\ξ}}_{i,k|k-1}-{\hat{x}}_{k|k-1}){({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T$

${\mathrm{\ξ}}_{0,k-1}={\hat{x}}_{k-1}$

${\mathrm{\ξ}}_{1,k-1}={\hat{x}}_{k-1}+\sqrt{(n+\mathrm{\λ})P_{k-1}}$

${\mathrm{\ξ}}_{2,k-1}={\hat{x}}_{k-1}=\sqrt{(n+\mathrm{\λ})P_{k-1}}$

其中，λ＝α²(n+κ)-n，α＝1e-3，α为正值的比例缩放因子，通过调整α的取值来调节Sigma点与的距离，同时调整α使高阶项的影响达到最小；α的取值范围为0≤α≤1。当系统为弱非线性时，α取较大的正值，α≥0.8当系统非线性程度严重时，α取一个非常小的正值,α≤1×10^-3，以避免采样点非局域效应的影响。

n为状态空间维数，n＝1，

κ为比例参数，选择κ＝0，

对于高斯函数分布的情况，最佳选择是β＝2，

是k-1时刻SOC值，ξ_0,k-1是k-1时刻比例对称修正采样的中心点，ξ_1,k-1和ξ_2,k-1分别是k-1时刻比例对称修正采样相对于中心点的对称点，P_k-1是k-1时刻的误差方差阵，W₀ ^m是一阶统计特性的中心点权值，W_i ^m是一阶统计特性的对称点权值，W₀ ^c是二阶统计特性的中心点权值，W_i ^c是二阶统计特性的对称点权值。

库伦效率η_k-1＝1。

ξ_i,k-1是k-1时刻中心点或者对称点的SOC值，γ_i,k|k-1是基于前一时刻(k-1时刻)状态的当前k时刻的预测的中心点或者对称点的SOC值，η_k-1是k-1时刻的库伦效率，C_n是电池的标称容量，△t表示采样间隔，I_k-1是k-1时刻电流，是基于前一时刻(k-1时刻)状态下k时刻的预测的SOC值，P_k|k-1是基于前一时刻(k-1时刻)状态下k时刻预测的误差方差阵，Q_k-1是k-1时刻过程噪声方差阵。ξ_0,k|k-1是基于前一时刻(k-1时刻)状态的当前k时刻预测比例对称修正采样的中心点，ξ_1,k|k-1和ξ_2,k|k-1分别是基于前一时刻(k-1时刻)状态的当前k时刻的预测比例对称修正采样相对于中心点的对称点，χ_i,k|k-1是基于前一时刻(k-1时刻)状态的当前k时刻的预测中心或者是对称点的电压值，ξ_i,k|k-1是基于前一时刻(k-1时刻)状态的当前k时刻的预测中心点或者对称点的SOC值，I_k是k时刻电流值，R₀是欧姆电阻，a₀,a₁,a₂,a₃,a₄是模型参数，是一步预测的电压值，R_k是k时刻量测噪声方差阵，是k时刻自协方差阵，是k时刻互协方差阵。

u_k是电池模型在k时刻的电压残差，F_k是k时刻电压残差方差阵的近似。Q_k是k时刻过程噪声方差阵。

Ⅲ、方差阵相匹配的鲁棒自适应无迹卡尔曼滤波(Robust AdaptiveUnscented Kalman Filter，RAUKF)

卡尔曼滤波的状态估计结果受到许多参数的影响，例如量测噪声方差阵R、过程噪声方差阵Q及滤波增益K。主要调整量测噪声方差阵R，过程噪声方差阵Q的无迹卡尔曼滤波属于自适应无迹卡尔曼滤波(Adaptive UnscentedKalman Filter，AUKF)；部分调整滤波增益的无迹自适应卡尔曼滤波属于鲁棒自适应无迹卡尔曼滤波(RAUKF)。

在公式(1)-(8)的自适应卡尔曼滤波算法估计SOC过程基础上，本发明结合鲁棒思想，提出了实际残差方差阵与估计的理论残差方差阵有效匹配方法。

鲁棒卡尔曼滤波中的系统信息

$e_k=z_k-{\hat{z}}_{k|k-1}---\left(9\right)$

定义系统实际残差方差阵C_k:

$C_k=\frac{{\mathrm{\Σ}}_{i=k-L+1}^k{e}_i{e_i}^T}{L}---\left(10\right)$

其中L为所选取的平滑窗口，根据经验法，取L＝16；系统理论残差方差阵为两者的比值为：

$\mathrm{\θ}=\frac{\mathrm{tr}\left(C_k\right)}{\mathrm{tr}\left(P_{{\tilde{z}}_k}\right)}---\left(11\right)$

其中tr(·)表示矩阵的迹。

当描述系统噪声统计特性准确时，即系统量测噪声的实际值与给定值一致，用于描述两者的匹配程度的残差方差实际值的迹与估计的理论值的迹之比值θ为1。

当系统量测噪声的实际值与给定值不一致时，θ偏离1。为了使卡尔曼滤波稳定，基于由虚拟的量测噪声方差阵和量测噪声方差阵R之比构成的调整因子ε_k，结合R_k本身自适应调节，对估计的理论残差方差阵进行调整，从而达到调整增益K的目的。

令θ＝1，即进行反推计算，过程如下：

得到调整因子ε_k

$\begin{array}{c}{\mathrm{\ϵ}}_k=\frac{C_k-\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^c({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T}{R_k}\\ =\frac{\tilde{R}}{R}\end{array}---\left(13\right)$

其中 $\tilde{R}=C_k-\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^c({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\χ}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T---\left(14\right)$

为虚拟的量测噪声方差阵。

以θ＝1为目标使得估计的理论残差方差阵能够反映实际残差方差阵，实现其方差的完全匹配。

Ⅳ、滤波增益自适应调整

AUKF算法中，量测噪声方差阵R是定值，认为主要受传感器偏差影响，量测噪声统计特性可根据传感器统计特性得知，但当传感器存在偏差时，量测噪声方差阵R将会比真实R值要小。AUKF算法认为R是固定，未考虑到R本身的变化，只是通过大于1的ε_k系数，以ε_kR形式，修正R。

本发明以量测噪声方差阵R本身自适应调整为前提，从新的角度，分析量测噪声方差阵R的所具有特殊性。量测噪声方差阵R的具体分析公式R_k的适应算法，

$R_k=F_k+\underset{i=0}{\overset{2}{\mathrm{\Σ}}}{W_i}^c({\mathrm{\χ}}_{i,k|k-1}-z_k+u_k){({\mathrm{\χ}}_{i,k|k-1}-z_k+u_k)}^T$

由此式得知，R不仅与传感器有关，还与状态变量有关。

通过公式(13)可知，虚拟的量测噪声方差阵可能为正值，或者是负值，对应的调整因子ε_k可能出现小于0、小于1或者是大于1的情况.所以当R本身自适应调整时，调整因子ε_k情况复杂，不局限于RUKF算法仅考虑的调整因子ε_k大于1的情况。

ε_k大于1时，把ε_k大于1带入公式(13)得知，系统实际残差方差阵C_k大于系统估计的理论残差方差阵为即θ大于1。也就是说，在R固定情况下，以ε_kR形式，修正对应量测噪声方差阵，只是考虑的θ的一种特殊情况，即对R只进行了部分修正。

θ反映的是实际残差方差阵与估计的理论残差方差阵匹配情况，完全匹配时，比值θ应为1，但实际中，系统实际残差方差阵C_k与系统估计的理论残差方差阵不相等，即比值θ可能出现θ>1或者0<θ<1。

ε_k<0时，算法初始阶段，初始值与真实值偏差大且噪声统计特性不准确，导致θ严重偏离1，为趋于0的小数值，此时ε_k可能为负值，ε_k不能作为调节因子直接调节R_k来对估计的理论残差方差阵进行调整，故使用自适应调整因子代替ε_k，取修正R_k，以避免对滤波增益过度调节所带来的协方差负定性问题，δ为衰减因子，0.3<δ<0.7。

当ε_k>0时，取自适应调整因子修正R_k，对估计的理论残差方差阵进行调整，使估计的理论的量测噪声统计特性动态跟踪实际的量测噪声特性。

本发明以实际的量测噪声统计特性为目标，通过反向推算，求取自适应调整因子并与R_k本身自适应调整相结合，以修正估计的理论计算得到量测噪声统计特性，实现实时跟踪实际的量测噪声特性。具体实现过程如下：

自适应调整因子带入公式(13)调整滤波增益K，得到新的滤波增益NK：

${\mathrm{NK}}_k=P_{{\tilde{x}}_k{\tilde{z}}_k}{(\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T+{\widetilde{\mathrm{\&epsiv;}}}_k{R}_k)}^{-1}---\left(15\right)$

$P_k=P_{k|k-1}-{\mathrm{NK}}_k{P}_{{\tilde{z}}_k}{{\mathrm{NK}}_k}^T---\left(16\right)$

P_k是k时刻的误差协方差，NK_k ^T是NK_k的转置。

求得k时刻SOC估计值

${\hat{x}}_k={\hat{x}}_{k|k-1}+{\mathrm{NK}}_k(z_k-{\hat{z}}_{k|k-1})---\left(17\right)$

当量测噪声水平发生变化时，比值θ为θ>1或者0<θ<1，与只调整比值θ大于1的情况不同，方差阵相匹配的RAUKF算法可以调整比值θ不为1的任何情况，提高了算法对随机量测噪声变化的鲁棒性，实现实际残差方差阵与估计的理论残差方差阵的完全匹配。

本发明量测噪声方差阵相匹配的电池电荷状态估计系统，包括微控制器及其连接的显示器，所述微控制器是嵌入式微控制器，动力电池输出端接有电压传感器和电流传感器。电压传感器和电流传感器经模数转换模块连接嵌入微控制器，嵌入微控制器含有模型参数辨识模块、调整因子计算模块、滤波增益自适应调整、鲁棒自适应无迹卡尔曼滤波SOC估计模块。嵌入微控制器连接显示器，还接有CAN(控制器局域网络Controller Area Network)总线接口和/或RS232接口。本系统与动力电池一起，嵌入在使用动力电池的设备中，在一个采样周期内完成电压、电流的采集、电池模型参数辨识及SOC在线估计，所得SOC结果在显示器上显示或直接传送到该设备的控制器局域网络。

与现有技术相比，本发明量测噪声方差阵相匹配的电池电荷状态估计方法与系统的优点为：提高电池电荷状态估计的鲁棒性，将量测噪声方差阵R本身自适应调整、实际残差方差阵与估计的理论残差方差阵相匹配相结合，与过程噪声方差阵Q相比，观察到量测噪声方差阵R对电池电荷状态估计具有更大的影响，对R本身进行自适应调整，以方差阵的完全匹配--残差方差实际值的迹与估计的理论值的迹相等为目标，分析了基于虚拟的量测噪声方差阵和量测噪声方差阵R之比的调整因子ε_k出现的各种情况，进行分类，动态地调整滤波增益，保证估计的理论残差方差阵与实际残差方差阵的完全匹配。与无迹卡尔曼滤波(AUKF)、现有的鲁棒无迹卡尔曼滤波(RAUKF)比较，本发明的方差阵相匹配的鲁棒无迹卡尔曼滤波(RAUKF)方法能准确估计电池电荷状态，即使在不同电池电荷状态初值下，也能较快地收敛到真值，具有良好的鲁棒性和收敛性，对噪声有较强的抑制作用。本发明方法所得电池电荷状态估计值的精度提高达0.64％

附图说明

图1为本量测噪声方差阵相匹配的电池电荷状态估计方法实施例的流程示意图；

图2为本量测噪声方差阵相匹配的电池电荷状态估计系统实施例的结构示意图；

图3为本量测噪声方差阵相匹配的电池电荷状态估计方法实施例所得SOC曲线、无迹卡尔曼滤波(AUKF)和现有的鲁棒无迹卡尔曼滤波(RAUKF)所得SOC曲线以及实际SOC曲线的比较图

图4为与实际SOC值比较，本量测噪声方差阵相匹配的电池电荷状态估计方法实施例所得SOC的误差曲线、无迹卡尔曼滤波(AUKF)和现有的鲁棒无迹卡尔曼滤波(RAUKF)所得SOC的误差曲线的比较图。

具体实施方式

量测噪声方差阵相匹配的电池电荷状态估计方法实施例

本量测噪声方差阵相匹配的电池电荷状态估计方法实施例的流程如图1所示，包括如下步骤：

Ⅰ、建立电池模型及参数辨识

Ⅰ-1、建立电池模型

本例采用联合模型，其数学表达式：

$x_{k+1}=x_k-\left(\frac{{\mathrm{\&eta;}}_k\mathrm{\&Delta;t}}{C_n}\right)I_k+{\mathrm{\&omega;}}_k---\left(1\right)$

$U_k=a_0-\frac{a_1}{x_k}-a_2{x}_k+a_3\mathrm{In}\left(x_k\right)+a_4\mathrm{In}(1-x_k)-R_0{I}_k+v_k---\left(2\right)$

其中U_k，x_k和I_k分别是当前k时刻电池端电压，SOC和电流。电池放电过程I_k为正值，充电过程I_k为负值；a₀,a₁,a₂,a₃和a₄是通过最小二乘法得到的相关参数；R₀是欧姆内阻。η_k是k时刻库伦效率，C_n是电池的标称容量，△t表示采样间隔。ω_k是k时刻系统过程噪声,ν_k是k时刻量测噪声。采样得到U_k、I_k。

Ⅰ-2、模型参数辨识

模型参数辨识方法采用最小二乘法得：

c＝[a₀,a₁,a₂,a₃,a₄,R₀]^T＝(M^TM)^-1M^TU     (3)其中,U＝[U₁,U₂,...,U_N]^T是端电压序列，M＝[M₁,M₂,...,M_N]^T是递归矩阵， $M_k=[1,\frac{1}{x_k},x_k,\mathrm{In}\left(x_k\right),\mathrm{In}(1-x_k),I_k].$

Ⅱ、自适应无迹卡尔曼滤波(AUKF)

自适应无迹卡尔曼滤波估计SOC过程如下：

Ⅱ-1、初始状态统计特性为：

${\hat{x}}_0=E\left(x_0\right)---\left(4\right)$

$P_0=\mathrm{Var}\left(x_0\right)=E\left[(x_0-{\hat{x}}_0){(x_0-{\hat{x}}_0)}^T\right]---\left(5\right)$

x₀是SOC初始值，是SOC初始值的平均值，P₀初始误差协方差阵。

Ⅱ-2、SOC的估计值的过程

SOC的估计值

${\hat{y}}_k={\hat{x}}_{k|k-1}+K_k(z_k-{\hat{z}}_{k|k-1})---\left(6\right)$

参数更新

$K_k=P_{{\tilde{x}}_k{\tilde{z}}_k}{\left(P_{{\tilde{z}}_k}\right)}^{-1}---\left(7\right)$

${\hat{z}}_{k|k-1}=\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^m{\mathrm{\&chi;}}_{i,k|k-1}---\left(8\right)$

是k时刻自协方差阵，是k时刻互协方差阵。K_k是k时刻的滤波增益，z_k是k时刻电池端电压的真实值，是k时刻的SOC值。

Ⅲ、方差阵相匹配的鲁棒自适应无迹卡尔曼滤波

鲁棒卡尔曼滤波中的系统信息

$e_k=z_k-{\hat{z}}_{k|k-1}---\left(9\right)$

定义系统实际残差方差阵C_k:

$C_k=\frac{{\mathrm{\&Sigma;}}_{i=k-L+1}^k{e}_i{e_i}^T}{L}---\left(10\right)$

其中L为所选取的平滑窗口，根据经验法，取L＝16；系统估计的理论残差方差阵为两者的比值为：

$\mathrm{\&theta;}=\frac{\mathrm{tr}\left(C_k\right)}{\mathrm{tr}\left(P_{{\tilde{z}}_k}\right)}---\left(11\right)$

其中tr(·)表示矩阵的迹。

当描述系统噪声统计特性准确时，即系统量测噪声的实际值与给定值一致，用于描述两者的匹配程度的残差方差实际值的迹与估计的理论值的迹之比值θ为1。

当系统量测噪声的实际值与给定值不一致时，θ偏离1。调整因子ε_k结合R_k本身自适应调节，对估计的理论残差方差阵进行调整。

得到调整因子ε_k

$\begin{array}{c}{\mathrm{\&epsiv;}}_k=\frac{C_k-\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T}{R_k}\\ =\frac{\tilde{R}}{R}\end{array}---\left(13\right)$

其中 $\tilde{R}=C_k-\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T---\left(14\right)$

为虚拟的量测噪声方差阵。

Ⅳ、滤波增益自适应调整

本例量测噪声方差阵R进行本身自适应调整，R_k的适应算法如下，

$R_k=F_k+\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-z_k+u_k){({\mathrm{\&chi;}}_{i,k|k-1}-z_k+u_k)}^T.$

ε_k<0时，取自适应调整因子修正R_k、对估计的理论残差方差阵进行调整以避免对滤波增益K_k过度调节所带来的协方差负定性问题，δ为衰减因子，本例即自适应调整因子

当ε_k>0时，取自适应调整因子修正R_k、对估计的理论残差方差阵进行调整，使估计的理论量测噪声统计特性动态跟踪实际的量测噪声特性。

求取自适应调整因子并与R本身自适应调整相结合，以修正R_k理论计算得到量测噪声统计特性，具体如下：

调整因子带入公式(13)调整滤波增益K，得到新的滤波增益NK_k：

${\mathrm{NK}}_k=P_{{\tilde{x}}_k{\tilde{z}}_k}{(\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T+{\widetilde{\mathrm{\&epsiv;}}}_k{R}_k)}^{-1}---\left(15\right)$

$P_k=P_{k|k-1}-{\mathrm{NK}}_k{P}_{{\tilde{z}}_k}{{\mathrm{NK}}_k}^T---\left(16\right)$

P_k是k时刻的误差协方差，NK_k ^T是NK_k的转置。

求得k时刻SOC估计值

${\hat{x}}_k={\hat{x}}_{k|k-1}+{\mathrm{NK}}_k(z_k-{\hat{z}}_{k|k-1}).---\left(17\right)$

量测噪声方差阵相匹配的电池电荷状态估计系统实施例

本量测噪声方差阵相匹配的电池电荷状态估计系统实施例如图2所示，嵌入式微控制器连接显示器，动力电池输出端接有电压传感器和电流传感器。电压传感器和电流传感器经模数转换模块连接嵌入微控制器，嵌入微控制器含有模型参数辨识模块、调整因子计算模块、滤波增益自适应调整、鲁棒自适应无迹卡尔曼滤波SOC估计模块。嵌入微控制器还接有CAN总线接口和RS232接口。本系统与动力电池一起，嵌入在使用动力电池的设备中，在一个采样周期内完成电压、电流的采集、电池模型参数辨识及SOC在线估计，所得SOC结果在显示器上显示或直接传送到该设备的控制器局域网络。

用上述本量测噪声方差阵相匹配的电池电荷状态估计系统实施例、按本量测噪声方差阵相匹配的电池电荷状态估计方法对某型号动力电池的电荷状态SOC的全量程变化进行十小时以上的实验得到电池的电荷状态估计值，同时用无迹卡尔曼滤波(AUKF)方法和现有的鲁棒无迹卡尔曼滤波(RAUKF)方法对相同型号的动力电池的电荷状态进行估计，传统实验所得的该型号动力电池的SOC作为SOC实际值。不同方法所得SOC的估计值曲线如图3所示，其中，纵坐标为SOC值，横坐标为时间，大图时间单位为秒(s)，局部放大图中时间单位为10^-4秒。图中实线为实际SOC曲线，点虚线为AUKF估计SOC曲线，短横虚线为RAUKF估计SOC曲线，点划线为本实施例估计SOC曲线。图中可清楚看到本实施例所得SOC曲线点划线与实际SOC曲线实线最为接近。当本实施例电荷状态初值与实际值相差较大的情况下，也能较快地收敛到真值。

图4所示为与实际SOC值比较，本量测噪声方差阵相匹配的电池电荷状态估计方法实施例所得SOC的误差曲线、无迹卡尔曼滤波(AUKF)和现有的鲁棒无迹卡尔曼滤波(RAUKF)所得SOC的误差曲线，其中，纵坐标为误差值，横坐标为时间，大图时间单位为秒(s)，局部放大图中时间单位为10^-4秒。点虚线为AUKF估计SOC误差曲线，短横虚线为RAUKF估计SOC误差曲线，点划线为本实施例估计SOC误差曲线。从图中清楚看到本实施例的点划线在实验过程中平稳且基本处于最低位置。

上述实施例，仅为对本发明的目的、技术方案和有益效果进一步详细说明的具体个例，本发明并非限定于此。凡在本发明的公开的范围之内所做的任何修改、等同替换、改进等，均包含在本发明的保护范围之内。

Claims (3)

1.量测噪声方差阵相匹配的电池电荷状态估计方法，其特征在于包括如下步骤：

Ⅰ、建立电池模型及参数辨识

Ⅰ-1、建立电池模型

采用联合模型，其数学表达式：

$x_{k+1}=x_k-\left(\frac{{\mathrm{\&eta;}}_k\mathrm{\&Delta;t}}{C_n}\right)I_k+{\mathrm{\&omega;}}_k---\left(1\right)$

$U_k=a_0-\frac{a_1}{x_k}-a_2{x}_k+a_3\mathrm{In}\left(x_k\right)+a_4\mathrm{In}(1-x_k)-R_0{I}_k+v_k---\left(2\right)$

其中U_k，x_k和I_k分别是当前k时刻电池端电压，SOC和电流；电池放电过程I_k为正值，充电过程I_k为负值；a₀,a₁,a₂,a₃和a₄是通过最小二乘法得到的相关参数；R₀是欧姆内阻；η_k是k时刻库伦效率，C_n是电池的标称容量，△t表示采样间隔；ω_k是k时刻系统过程噪声,ν_k是k时刻量测噪声；采样得到U_k、I_k；

Ⅰ-2、模型参数辨识

模型参数辨识方法采用最小二乘法得：

c＝[a₀,a₁,a₂,a₃,a₄,R₀]^T＝(M^TM)^-1M^TU            (3)

其中,U＝[U₁,U₂,...,U_N]^T是端电压序列，M＝[M₁,M₂,...,M_N]^T是递归矩阵， $M_k=[1,\frac{1}{x_k},x_k,\mathrm{In}\left(x_k\right),\mathrm{In}(1-x_k),I_k];$

Ⅱ、自适应无迹卡尔曼滤波

自适应无迹卡尔曼滤波估计SOC过程如下：

Ⅱ-1、初始状态统计特性为：

${\hat{x}}_0=E\left(x_0\right)---\left(4\right)$

$P_0=\mathrm{Var}\left(x_0\right)=E\left[(x_0-{\hat{x}}_0){(x_0-{\hat{x}}_0)}^T\right]---\left(5\right)$

x₀是SOC初始值，是SOC初始值的平均值，P₀为初始误差协方差阵；

Ⅱ-2、SOC的估计值的过程

SOC的估计值

${\hat{y}}_k={\hat{x}}_{k|k-1}+K_k(z_k-{\hat{z}}_{k|k-1})---\left(6\right)$

参数更新

$K_k=P_{{\tilde{x}}_k{\tilde{z}}_k}{\left(P_{{\tilde{z}}_k}\right)}^{-1}---\left(7\right)$

${\hat{z}}_{k|k-1}=\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^m{\mathrm{\&chi;}}_{i,k|k-1}---\left(8\right)$

是k时刻自协方差阵，是k时刻互协方差阵，K_k是k时刻的滤波增益，z_k是k时刻电池端电压的真实值，是k时刻的SOC值；

Ⅲ、方差阵相匹配的鲁棒自适应无迹卡尔曼滤波

鲁棒卡尔曼滤波中的系统信息

$e_k=z_k-{\hat{z}}_{k|k-1}---\left(9\right)$

定义系统实际残差方差阵C_k:

$C_k=\frac{{\mathrm{\&Sigma;}}_{i=k-L+1}^k{e}_i{e_i}^T}{L}---\left(10\right)$

其中L为所选取的平滑窗口，根据经验法，取L＝16；系统估计的理论残差方差阵为两者的比值为：

$\mathrm{\&theta;}=\frac{\mathrm{tr}\left(C_k\right)}{\mathrm{tr}\left(P_{{\tilde{z}}_k}\right)}---\left(11\right)$

其中tr(·)表示矩阵的迹；

当描述系统噪声统计特性准确时，即系统量测噪声的实际值与给定值一致，用于描述两者的匹配程度的残差方差实际值的迹与估计的理论值的迹之比值θ为1；

当系统量测噪声的实际值与给定值不一致时，θ偏离1；调整因子ε_k结合R_k本身自适应调节，对估计的理论残差方差阵进行调整；

令θ＝1，即C_k代入：

$P_{{\tilde{z}}_k}=\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T+{\mathrm{\&epsiv;}}_k{R}_k---\left(12\right)$

得到调整因子ε_k

$\begin{array}{c}{\mathrm{\&epsiv;}}_k=\frac{C_k-\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T}{R_k}\\ =\frac{\tilde{R}}{R}\end{array}---\left(13\right)$

其中 $\tilde{R}=C_k-\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T---\left(14\right)$

为虚拟的量测噪声方差阵；

Ⅳ、滤波增益自适应调整

本例量测噪声方差阵R进行本身自适应调整，R_k的适应算法如下，

$R_k=F_k+\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-z_k+u_k){({\mathrm{\&chi;}}_{i,k|k-1}-z_k+u_k)}^T;$

ε_k<0时，取自适应调整因子修正R_k，δ为衰减因子；

当ε_k为正值时，取自适应调整因子修正R_k；

求取自适应调整因子并与R本身自适应调整相结合，以修正理论计算得到量测噪声统计特性，具体如下：

自适应调整因子带入公式(13)调整滤波增益，得到新的滤波增益NK_k：

${\mathrm{NK}}_k=P_{{\tilde{x}}_k{\tilde{z}}_k}{(\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}{W_i}^c({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1}){({\mathrm{\&chi;}}_{i,k|k-1}-{\hat{z}}_{k|k-1})}^T+{\widetilde{\mathrm{\&epsiv;}}}_k{R}_k)}^{-1}---\left(15\right)$

P_k是k时刻的误差协方差，NK_k ^T是NK_k的转置；

$P_k=P_{k|k-1}-{\mathrm{NK}}_k{P}_{{\tilde{z}}_k}{{\mathrm{NK}}_k}^T---\left(16\right)$

求得k时刻SOC估计值

${\hat{x}}_k={\hat{x}}_{k|k-1}+{\mathrm{NK}}_k(z_k-{\hat{z}}_{k|k-1})---\left(17\right).$

2.根据权利要求1所述的量测噪声方差阵相匹配的电池电荷状态估计方法，其特征在于：

所述步骤Ⅳ中衰减因子0.3<δ<0.7。

3.根据权利要求1或2所述的量测噪声方差阵相匹配的电池电荷状态估计方法设计的量测噪声方差阵相匹配的电池电荷状态估计系统，包括微控制器及其连接的显示器，其特征在于：

所述微控制器是嵌入式微控制器，动力电池输出端接有电压传感器和电流传感器；电压传感器和电流传感器经模数转换模块连接嵌入微控制器，嵌入微控制器含有模型参数辨识模块、调整因子计算模块、滤波增益自适应调整、鲁棒自适应无迹卡尔曼滤波SOC估计模块，还接有CAN总线接口和/或RS232接口；本系统与动力电池一起，嵌入在使用动力电池的设备中，在一个采样周期内完成电压、电流的采集、电池模型参数辨识及SOC在线估计，所得SOC结果在显示器上显示或直接传送到该设备的控制器局域网络。