【滤波估计】双卡尔曼滤波SOC和SOH联合估计【含Matlab源码 2335期】

最新推荐文章于 2024-04-30 10:44:36 发布

Matlab武动乾坤

最新推荐文章于 2024-04-30 10:44:36 发布

阅读量1.1k

点赞数 26

分类专栏： Matlab信号处理（进阶版）文章标签： matlab

本文链接：https://blog.csdn.net/KeepingMatlab/article/details/135243635

版权

Matlab信号处理（进阶版）专栏收录该内容

490 篇文章 39 订阅

订阅专栏

⛄一、双卡尔曼滤波SOC和SOH联合估计

1 引言
为实现节能降耗，降低污染，发展节能环保、不依赖化石燃料的电动汽车取代传统燃油车，已成为当今世界汽车行业的重点发展方向。锂离子动力电池准确可靠的状态估计是电动汽车安全运行的基础[1]，其主要包括荷电状态（SOC）和健康状态（SOH）。SOC直接反映了电池剩余电量的大小，其准确估计直接关系到电动汽车的能量动力分配。

SOH是电池老化程度的一项重要指标，通常表现为电池的能量密度、功率密度、容量的衰减和内部电阻增大[2]，电池状态的准确估计可以使电池得到充分合理的利用，避免电池突发故障造成的危害，对于电动汽车的安全运行具有重要意义[3]。

目前国内应用最多的电池SOC测量方法是安时积分法，它计算简单，但受初始SOC误差的影响较大，且会随时间增长出现较大的累计误差[4]。另外，电池SOC测量方法还有开路电压法、卡尔曼滤波法、神经网络法。开路电压法容易实现但需要通过静置校准OCV值，不利于实际运用[5]。

神经网络法通过训练大量的样本数据进行估算，准确性受训练的方法及训练量大小影响较大，计算量太大。卡尔曼滤波法使用递推迭代的方法对SOC进行估算[6]。最常用的为扩展卡尔曼滤波（EKF）和无迹卡尔曼滤波（UKF):EKF估算过程中容易由于线性化过程中方差矩阵的非正定性导致估计值不收敛；UKF利用无迹变换处理均值和协方差，可以有效提高对非线性分布统计量的估算精度，但是系统噪声的预定变量设置不当通常会导致系统误差和发散[6]。

国内外常用的SOH估计方法有定义法、电化学阻抗谱分析法、容量衰减法等。定义法是以电池的定义为基础，应用电池内部特征量变化前后的比值关系来对电池的进行估计，此方法需要反复充放电实验，在实际中很难实现[7]。电化学阻抗谱分析法对仪器的精度要求高，投入较大，性价比不高。容量衰减法是通过数学建模实验得出容量衰减的规律，属于定性的方法分析SOH[8]。

针对目前电池状态估计精度和实用性不理想的问题，常以表征电池老化的主要参数容量和内阻来定义SOH。相对来说，电池容量与SOH之间的关系更为直接且容易采集[4]，所以这里选择容量作为估计SOH的状态变量，其定义为在一定条件下锂离子动力电池能够充入或放出的最大容量占其额定容量的百分比[9]，如下式所示：
在这里插入图片描述
基于以上分析，目前锂离子动力电池状态参数的研究中，多是将SOC与SOH分开研究，对两者在电池状态参数估算过程中的相互影响考虑较少[10]。

为此，这里提出自适应无迹卡尔曼滤波算法（AUKF），并在其基础上，采用双自适应无迹卡尔曼滤波（DAUKF）对SOC和SOH联合估算，通过两个滤波器的结合，实时更新计算模型中的参数。最后通过实验验证了算法的可行性及估算精度。

2 UKF算法原理
无迹卡尔曼滤波算法（UKF）是通过UT变换，在状态变量附近选取采样点，对于给定的非线性系统，空间状态方程可以由下式描述：
在这里插入图片描述
式中：f—线性状态方程函数；h—非线性观测方程函数。用UKF

算法估计状态变量具体过程如下所示：

(1）算法初始化：
在这里插入图片描述

(3)2n+1个Sigma点集进一步预测，i=1,2，⋯，2n+1:

(4）加权求取均值和协方差：

(5）再次使用UT变换，产生新的Sigma点集：

(6）求新的Sigma点集的观测预测值，利用加权求和求得系统预测的均值和协方差：
在这里插入图片描述
(7）计算Kalman增益矩阵：

(8）最后，计算系统的状态更新和协方差更新：

⛄二、部分源代码

clear;clc;
close all
%% 模型参数
R1=0.0019;
R2=0.0035;
C1=23340;
C2=501270;
load(‘discharge.mat’);%放电数据
load(‘OCV_SOC.mat’);%OCV-SOC关系
Ts=0.1;%采样间隔
%% 矩阵
A=[1-1Ts/R1/C1 0 0;0 1-1Ts/R2/C2 0;0 0 1];%系统矩阵
C=[-1 -1 0];
D=0;
%% 初始值
Q=0.00000001eye(3);%系统误差协方差
R=1;%测量误差协方差
Xekf=[0;0;0.8];%[U1,U2,SOC]初始值
Pa_ekf=[0.0037;30.233600]; %初始值 [R0 最大可用容量]
P0=0.01*[0.1 0 0;0 0.1 0;0 0 1];%状态误差协方差初始值
%% 赋值
tm=discharge(1,:)‘;%时间
Cur=-discharge(2,:)’;%电流
Vot=discharge(3,:)‘;%测量得到的端电压
RSOC=discharge(4,:)’;%SOC真实值-安时法计算得到
T=length™-1;%时间
%% emf-soc关系
x=OCV_SOC(2,:);%SOC
y=OCV_SOC(1,:);%OCV
p=polyfit(x,y,8);%多项式参数值
%% 初始化
L_discharge=length(discharge);
Uoc=zeros(1,L_discharge);%OCV
H=zeros(L_discharge,3); %dUt/dx
Vekf=zeros(1,L_discharge);%估计得到的端电压值
K=zeros(3,L_discharge);%卡尔曼增益
C_1=zeros(L_discharge,2);
C_2=zeros(L_discharge,2);
d_g_PA=zeros(L_discharge,2); %dUt/d参数（内阻，容量）
K_pa=zeros(2,8137); %卡尔曼增益，参数估计中
Pa=zeros(1,L_discharge); %容量估计值
Pa(1)=Pa_ekf(2);
SOH=zeros(1,L_discharge); %SOH
SOH(1)=1; %SOH初始值
Uoc(1)=p(1)*Xekf(3)^{8+p(2)*Xekf(3)}7+p(3)*Xekf(3)^{6+p(4)*Xekf(3)}5+p(5)*Xekf(3)^{4+p(6)*Xekf(3)}3+p(7)*Xekf(3)^2+p(8)Xekf(3)+p(9);%OCV
Vekf(1)=Uoc(1)+CXekf-Cur(1)*Pa_ekf(1);%估计得到的端电压值
counter=0; %计数