CN114417686B - 一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法 - Google Patents

一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法 Download PDF

Info

Publication number
CN114417686B
CN114417686B CN202210065461.1A CN202210065461A CN114417686B CN 114417686 B CN114417686 B CN 114417686B CN 202210065461 A CN202210065461 A CN 202210065461A CN 114417686 B CN114417686 B CN 114417686B
Authority
CN
China
Prior art keywords
battery
state
covariance
service life
degradation
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202210065461.1A
Other languages
English (en)
Other versions
CN114417686A (zh
Inventor
张九思
李翔
罗浩
蒋宇辰
霍明夷
田纪伦
尹珅
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Institute of Technology
Original Assignee
Harbin Institute of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology filed Critical Harbin Institute of Technology
Priority to CN202210065461.1A priority Critical patent/CN114417686B/zh
Publication of CN114417686A publication Critical patent/CN114417686A/zh
Application granted granted Critical
Publication of CN114417686B publication Critical patent/CN114417686B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/25—Design optimisation, verification or simulation using particle-based methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/08—Probabilistic or stochastic CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/02—Reliability analysis or reliability optimisation; Failure analysis, e.g. worst case scenario performance, failure mode and effects analysis [FMEA]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/04—Ageing analysis or optimisation against ageing

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Medical Informatics (AREA)
  • Software Systems (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Artificial Intelligence (AREA)
  • Secondary Cells (AREA)
  • Tests Of Electric Status Of Batteries (AREA)

Abstract

一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法，本发明涉及电池的自适应在线剩余使用寿命预测方法。本发明的目的是为了解决现有单个电池工作过程中的剩余使用寿命预测方法需要大量历史数据，在一定程度上会受到强马尔可夫特性的制约，以及电池在退化过程中的随机性难以得到量化的问题。过程为：一：提取电池第k次工作过程中的容量数据；二：构建电池的退化模型；三：构建自适应估计退化模型参数；四：计算电池的剩余使用寿命；判断第k次工作循环下的电池的荷电状态SOC是否低于80％，如果SOC高于80％，则令k＝k+1，执行一至四，否则电池报废；五：评估剩余使用寿命的在线预测效果。本发明用于电池使用寿命预测领域。

Description

一种面向单个锂离子电池的自适应在线剩余使用寿命预测 方法
技术领域
本发明涉及锂离子电池的剩余使用寿命预测和统计分析相结合的学科交叉领域，具体涉及面向单个锂离子电池的自适应在线剩余使用寿命预测方法。
背景技术
锂离子电池由于其电化学特性稳定、能量密度高、充放电速度快等优点被广泛地应用于交通运输、航空航天、电力存储、机械制造等领域。
在锂离子电池的工作过程中，由于其内部物理结构和外部环境条件的影响，电池的容量会逐渐下降，这将对实际的用电系统或者设备的可靠性和安全性产生很大的负面影响，严重时甚至会发生灾难性的事故。因此，在用电设备的电池管理系统(BatteryManagement System，BMS)中，准确地获得锂离子电池的退化信息是至关重要的。剩余使用寿命(Remaining useful Life，RUL)是衡量锂离子电池退化趋势的重要性指标。剩余使用寿命指的是系统首次从当前状态到完全失效时所经历的时间。准确的预测锂离子的剩余使用寿命可以充分发挥电池的动力性能，提高电池使用的安全性，其意义在于及时掌握电池的运行情况，协助操作人员及时对电池进行维护和更换，降低电力系统或者设备的运行成本。
目前锂离子电池的剩余使用寿命预测方法从总体上可以分为基于物理化学模型的方法、数据驱动的方法和基于混合经验模型的方法。由于电池的内部工作机理复杂、退化机理尚不明晰，基于物理化学模型的方法在预测剩余使用寿命时存在着一定的困难。另一方面，数据驱动的方法在本质上是黑箱模型，需要大量的历史数据训练算法模型，并且不能量化电池在退化过程中的不确定性。不仅如此，在线数据还被要求遵循与历史数据相似的分布，否则学习类的方法很难保证预测的准确性。基于混合经验模型的方法则是将先验知识融入到算法模型中，根据实际情况，可以通过构建电池的容量和电流、电压以及循环次数等参数之间的算法模型，然后对模型中的参数进行辨识，进而完成剩余使用寿命的预测任务，具备很强的解释性和可靠性。
值得一提的是，传统的数据驱动方法需要大量的历史数据，如何仅仅使用单个电池工作过程中的数据进行剩余使用寿命预测存在着一定的挑战。不仅如此，现有的一些预测方法在一定程度上会受到强马尔可夫特性的制约，如何以自适应的方式利用电池从初始运行到当前工作循环的数据有待于进一步的研究。最后，电池在退化过程中的随机性往往难以得到量化，如何很好地描述退化过程中的不确定性需要深入的研究。
发明内容
本发明的目的是为了解决现有单个电池工作过程中的剩余使用寿命预测方法需要大量历史数据，在一定程度上会受到强马尔可夫特性的制约，以及电池在退化过程中的随机性难以得到量化等问题，而提出一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法。
一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法具体过程为：
步骤一：提取电池第k次工作过程中的容量数据；
步骤二：基于步骤一构建电池的退化模型；
步骤三：基于步骤二构建自适应估计退化模型参数；
步骤四：基于步骤三计算电池的剩余使用寿命；
获得电池的剩余使用寿命后，判断第k次工作循环下的电池的荷电状态SOC是否低于80％，即计算电池当前的容量Qp与额定容量Qf之间的比值是否低于80％，如果SOC高于80％，则令k＝k+1，执行步骤一至步骤四，否则电池报废；
步骤五：基于步骤四评估剩余使用寿命的在线预测效果。
本发明的有益效果为：
本发明面向单个电池提出了一种基于卡尔曼滤波-最大期望-Rauch-Tung-Striebel(RTS)的自适应方法预测电池的剩余使用寿命，从而提高电池使用的安全性，优化电子设备的工作性能。
本发明旨在解决单个电池工作过程中的剩余使用寿命预测问题。在无需大量历史数据的前提下，仅仅依据电池容量数据构建基于维纳过程的卡尔曼滤波动态退化模型描述电池的退化过程。采用基于最大期望-Rauch-Tung-Striebel的自适应地估计动态退化模型中的隐含参数，并计算剩余使用寿命的概率密度函数，进而预测电池的剩余使用寿命。
本发明针对单个电池提出了一种基于卡尔曼滤波-最大期望-Rauch-Tung-Striebel的自适应方法预测电池的剩余使用寿命。在无需大量历史数据的前提下，仅仅依据电池容量数据构建基于维纳过程的卡尔曼滤波动态退化模型描述电池的退化过程。采用基于最大期望-Rauch-Tung-Striebel的自适应地估计动态退化模型中的隐含参数，并计算剩余使用寿命的概率密度函数，进而预测电池的剩余使用寿命。
1、传统的数据驱动方法需要大量的历史数据。然而，在现实情况中，尤其是对于新安装的或者关键性系统中的电池，获得多个同类电池的历史退化数据是十分困难的。本发明所提出的方法仅仅使用单个电池工作过程中的数据进行剩余使用寿命预测。
2、本发明所提出的方法能够减少强马尔可夫特性的制约，通过最大期望-Rauch-Tung-Striebel的自适应方法利用电池从初始运行到当前工作循环的数据进行模型参数的更新。
3、由于电池在退化过程中的随机性往往难以得到量化，本发明提出的方法能够通过计算剩余使用寿命的概率密度函数来描述电池退化过程中的不确定性。
附图说明
图1是本发明的工作流程图；
图2是NASA Ames预测中心电提供的数据集中4个电池的电池容量衰减曲线示意图；
图3a是本发明中剩余使用寿命预测方法对于B0005电池中退化模型参数电池容量的退化速率初始值
Figure GDA0003859163290000034
的自适应估计结果图，Estimated为估计的；
图3b是本发明中剩余使用寿命预测方法对于B0005电池中退化模型参数电池容量的初始值
Figure GDA0003859163290000032
的自适应估计结果图；
图3c是本发明中剩余使用寿命预测方法对于B0005电池中退化模型参数电池容量的退化速率方差初始值
Figure GDA0003859163290000031
的自适应估计结果图；
图3d是本发明中剩余使用寿命预测方法对于B0005电池中退化模型参数电池容量方差初始值
Figure GDA0003859163290000033
的自适应估计结果图；
图3e是本发明中剩余使用寿命预测方法对于B0005电池中退化模型中电池容量退化速率过程噪声方差
Figure GDA0003859163290000035
的自适应估计结果图；
图3f是本发明中剩余使用寿命预测方法对于B0005电池中退化模型中电池容量过程噪声方差β2的自适应估计结果图；
图3g是本发明中剩余使用寿命预测方法对于B0005电池中退化模型参数测量噪声方差
Figure GDA0003859163290000036
的自适应估计结果图；
图4是本发明中剩余使用寿命预测方法对于B0005电池在不同循环下剩余使用寿命的概率密度函数结果示意图；
图5a是本发明中剩余使用寿命预测方法在NASA电池数据集上B0005电池的剩余使用寿命预测值和真实值的结果图；
图5b是本发明中剩余使用寿命预测方法在NASA电池数据集上B0006电池的剩余使用寿命预测值和真实值的结果图；
图5c是本发明中剩余使用寿命预测方法在NASA电池数据集上B0007电池的剩余使用寿命预测值和真实值的结果图；
图5d是本发明中剩余使用寿命预测方法在NASA电池数据集上B0018电池的剩余使用寿命预测值和真实值的结果图。
具体实施方式
具体实施方式一：本实施方式一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法具体过程为：
本发明旨在解决单个电池工作过程中的剩余使用寿命预测问题。在无需大量历史数据的前提下，仅仅依据电池容量数据构建基于维纳过程的卡尔曼滤波动态退化模型描述电池的退化过程。采用基于最大期望-Rauch-Tung-Striebel(RTS)的自适应地估计动态退化模型中的隐含参数，并计算剩余使用寿命的概率密度函数，进而能够在线实时预测电池的剩余使用寿命。本发明的工作流程图如图1所示。
步骤一：提取电池第k次工作过程中的容量数据，过程为：在电池每次工作循环完成后进行充放电测试，获得相应循环下的电池容量数据，作为电池自适应剩余使用寿命预测方法的输入；
步骤二：基于步骤一构建电池的退化模型；
步骤三：基于步骤二构建自适应估计退化模型参数；
步骤四：基于步骤三计算电池的剩余使用寿命；
获得电池的剩余使用寿命后，判断第k次工作循环下的电池的荷电状态SOC是否低于80％，即计算电池当前的容量Qp与额定容量Qf之间的比值是否低于80％，如果SOC高于80％，则令k＝k+1，执行步骤一至步骤四，否则电池报废；
步骤五：基于步骤四评估剩余使用寿命的在线预测效果。
具体实施方式二：本实施方式与具体实施方式一不同的是，所述步骤二中基于步骤一构建电池的退化模型(公式(4))；具体过程为：
仅仅基于步骤一提取的电池容量数据，构建基于维纳过程的卡尔曼滤波动态退化模型描述电池的退化过程，通过对于状态变量和协方差的估计、卡尔曼滤波增益的计算以及状态变量和协方差的更新的步骤获得隐含的参数变量。维纳过程是基于标准布朗(Brownian)运动的一种具备线性漂移系数的扩散过程。其中，布朗运动最开始用于描述微小粒子的随机运动，其本质为均值是0，方差是与时间有关的高斯过程。
步骤二一、构建基于维纳过程的卡尔曼滤波动态退化模型；具体过程为：
将电池容量作为描述电池退化过程中的退化量，维纳过程{q(t),t≥0}表达成公式(1)所示的形式：
q(t)＝αt+βB(t) (1)
其中，q(t)为t时刻下的电池容量；α为描述电池容量的退化速率，与电池的退化过程密切相关，又可称为漂移系数；β为描述电池容量退化过程随时间的不确定性程度，又可称为扩散系数；B(t)为标准布朗运动，且满足B(t)～N(0,t)，N(0,t)表示服从均值为0，方差为t的正态分布；
公式(1)的退化过程可以具体地表示成公式(2)和(3)所示的形式：
Figure GDA0003859163290000051
qk＝qk-1+αk-1(tk-tk-1)+βζk (3)
其中，
Figure GDA0003859163290000052
为描述电池容量退化速率过程中的不确定性，
Figure GDA0003859163290000053
为表示服从均值为0，方差为
Figure GDA0003859163290000054
的正态分布，
Figure GDA0003859163290000055
为
Figure GDA0003859163290000056
的方差，αk为电池在第k次工作循环下的容量退化速率，αk-1为电池在第k-1次工作循环下的容量退化速率，qk为电池在第k次工作循环下的电池容量，qk-1为电池在第k-1次工作循环下的电池容量，tk为第k次工作循环下的时间，tk-1为第k-1次工作循环下的时间，ζk为扩散系数的噪声项，ζk～N(0,tk-tk-1)，N(0,tk-tk-1)为表示服从均值为0，方差为tk-tk-1的正态分布；
考虑到电池容量退化速率在当前循环依赖于前一循环的值，且在电池不同的健康状态是不一致的，具备时变性，故采用
Figure GDA0003859163290000058
来描述电池容量退化速率过程中的不确定性。而引入βζk的目的在于描述公式(1)中的标准布朗运动。
步骤二二、构建卡尔曼状态空间模型；具体过程为：
为了更好地采用卡尔曼滤波算法进行状态估计，公式(2)和(3)可以采用公式(4)中的动态系统状态空间模型进行重构，以描述电池的退化；
Figure GDA0003859163290000057
其中，xk为第k次工作循环下的状态向量，包含了电池容量退化速率和电池容量两个状态变量，xk∈R2×1，R为实数；xk-1为第k-1次工作循环下的状态向量，Ak为第k次工作循环下的状态空间方程系数，Ak∈R2×2；wk为过程噪声项，wk∈R2×1；yk为电池在第k次工作循环下的实际容量，Ck为第k次工作循环下的输出方程yk＝Ckxk+ek的系数，Ck∈R1×2；ek为测量噪声项；
步骤二三、定义卡尔曼状态空间模型参数；具体过程为：
卡尔曼状态空间模型由：隐含状态向量x1:k＝[x1,x2,…,xk]、可观测的测量向量y1:k＝[y1,y2,…,yk]以及状态空间模型参数Ξ＝[Ak,wk,Ck,ek]所构成；
给出卡尔曼滤波的相关定义，具体如公式(6)—公式(9)所示：
Figure GDA0003859163290000061
其中，
Figure GDA0003859163290000062
为退化状态滤波均值，
Figure GDA0003859163290000063
和
Figure GDA0003859163290000064
分别为电池容量退化速率和电池容量这两个状态变量的退化状态滤波均值；E(αk|y1:k)为αk在已知测量数据y1:k条件下的后验期望，E(xk|y1:k)为状态向量xk在已知测量数据y1:k条件下的后验期望，E(qk|y1:k)为qk在已知测量数据y1:k条件下的后验期望；
Figure GDA0003859163290000065
其中，Sk|k为退化状态滤波协方差，var(αk|y1:k)和var(qk|y1:k)分别为电池容量退化速率和电池容量两个状态变量的退化状态滤波方差，cov(αk,qk|y1:k)为电池容量退化速率和电池容量两个状态变量之间的退化状态滤波协方差，cov(xk|y1:k)为状态向量xk在已知测量数据y1:k条件下的协方差矩阵；
Figure GDA0003859163290000066
其中，
Figure GDA0003859163290000067
为退化状态一步预测均值，
Figure GDA0003859163290000068
和
Figure GDA00038591632900000610
分别为两个状态变量的退化状态一步预测均值，
Figure GDA00038591632900000611
为
Figure GDA00038591632900000612
在已知测量数据y1:k-1条件下的先验期望，E(αk|y1:k-1)为αk在已知测量数据y1:k-1条件下的先验期望，E(qk|y1:k-1)为qk在已知测量数据y1:k-1条件下的先验期望；
Figure GDA0003859163290000071
其中，Sk|k-1为退化状态一步预测协方差，var(αk|y1:k-1)和var(qk|y1:k-1)分别为电池容量退化速率和电池容量两个状态变量的退化状态一步预测方差，cov(αk,qk|y1:k-1)为电池容量退化速率和电池容量两个状态变量之间的退化状态滤波协方差，cov(xk|y1:k-1)为状态向量xk在已知测量数据y1:k-1条件下的协方差矩阵；
基于公式(6)—(9)，当获得了电池容量在第k个工作循环的取值后，即可以采用卡尔曼滤波的方法对于隐含的状态向量xk进行估计；如公式(10)—(14)所示：
步骤二四、基于步骤二三估计状态与协方差；
步骤二五、基于步骤二四计算卡尔曼滤波增益；
步骤二六、基于步骤二五更新均值与协方差。
当分别给定电池在退化过程中公式(4)-(5)卡尔曼状态空间方程的均值和协方差的初始值
Figure GDA0003859163290000075
和S0|0后，即可以根据卡尔曼滤波方法进行系统状态的估计和更新。
其它步骤及参数与具体实施方式一相同。
具体实施方式三：本实施方式与具体实施方式一或二不同的是，所述步骤二二中第k次工作循环下的状态向量xk、第k次工作循环下的状态空间方程系数Ak、过程噪声项wk、第k次工作循环下的输出方程yk＝Ckxk+ek的系数Ck、测量噪声项ek的表达式为：
Figure GDA0003859163290000072
其中，wk的协方差矩阵表示为
Figure GDA0003859163290000073
v为输出方程yk＝Ckxk+ek的噪声项，噪声项
Figure GDA0003859163290000074
为测量噪声v的方差。
其它步骤及参数与具体实施方式一或二相同。
具体实施方式四：本实施方式与具体实施方式一至三之一不同的是，所述步骤二四中基于步骤二三估计状态与协方差；具体过程为：
基于公式(6)—(9)，当获得了电池容量在第k个工作循环的取值后，即可以采用卡尔曼滤波的方法对于隐含的状态向量xk进行估计；如公式(10)—(14)所示：
状态估计：
Figure GDA0003859163290000081
其中，
Figure GDA0003859163290000082
为状态向量的前向一步预测结果；
协方差估计：
Figure GDA0003859163290000083
其中，Sk-1|k-1为状态向量xk在第k-1个工作循环的协方差矩阵。
其它步骤及参数与具体实施方式一至三之一相同。
具体实施方式五：本实施方式与具体实施方式一至四之一不同的是，所述步骤二五中基于步骤二四计算卡尔曼滤波增益；具体过程为：
计算卡尔曼滤波增益Kk：
Figure GDA0003859163290000084
其中，
Figure GDA0003859163290000085
为测量噪声的方差。
其它步骤及参数与具体实施方式一至四之一相同。
具体实施方式六：本实施方式与具体实施方式一至五之一不同的是，所述步骤二六中基于步骤二五更新均值与协方差；具体过程为：
状态更新(对于电池容量退化速率和电池容量这两个状态变量的更新)：
Figure GDA0003859163290000086
协方差更新：
Sk|k＝Sk|k-1-KkCkSk|k-1 (14)。
其它步骤及参数与具体实施方式一至五之一相同。
具体实施方式七：本实施方式与具体实施方式一至六之一不同的是，所述步骤三中基于步骤二自适应估计退化模型参数；具体过程为：
采用基于最大期望-Rauch-Tung-Striebel的自适应地估计动态退化模型中的隐含参数。具体而言，针对基于维纳过程的卡尔曼滤波动态退化模型中所存在的未知参数，采用最大期望-Rauch-Tung-Striebel的方法进行自适应参数估计，利用电池从开始工作到当前工作循环的电池容量数据对退化模型中的状态变量进行递归更新。其中未知的参数向量可以表示为：
Figure GDA0003859163290000091
步骤三一、构建联合对数似然函数；具体过程为：
在已知自初始到第k个工作循环的电池容量数据y1:k＝[y1,y2,…,yk]的前提下，构建联合对数似然函数如公式(15)所示：
Figure GDA0003859163290000092
其中，Θ为公式(4)-(5)中退化方程的未知参数估计向量，
Figure GDA0003859163290000093
Figure GDA0003859163290000094
为均值的初始值，S0|0为协方差的初始值，p(x0:k,y1:k|Θ)为x0:k和y1:k在Θ下的联合概率密度函数，x0:k＝[x0,x1,...,xk]，xk为第k个工作循环的状态向量，p(y1:k|x0:k,Θ)为y1:k在x0:k,Θ下的条件概率密度，y1:k|x0:k,Θ为y1:k和x0:k,Θ之间的概率关系，p(x0:k|Θ)为x0:k的先验概率密度函数，y1:k|x0:k为y1:k和x0:k之间的概率关系，p(x0|Θ)为x0的先验概率密度函数，x0|Θ为x0和Θ之间的概率关系，p(xi|xi-1,Θ)为xi在xi-1,Θ下的条件概率密度，xi|xi-1为xi和xi-1之间的概率关系，p(yi|xi,Θ)为yi在xi,Θ下的条件概率密度，yi|xi为yi和xi之间的概率关系，x0为电池未工作时的状态向量，xi为第i个工作循环的状态向量，xi-1为第i-1个工作循环的状态向量，i＝1，2，…，k-1，k；
基于步骤二中公式(11)的卡尔曼滤波的自适应状态更新以及公式(4)中的关系式，获得公式(16)—(18)：
Figure GDA0003859163290000095
xi|xi-1,Θ～N(Aixi-1,Qi) (17)
Figure GDA0003859163290000101
其中，
Figure GDA0003859163290000102
表示服从均值为
Figure GDA0003859163290000103
方差为S0|0的正态分布，N(Aixi-1,Qi)表示服从均值为Aixi-1，方差为Qi的正态分布，
Figure GDA0003859163290000104
表示服从均值为Cixi，方差为
Figure GDA0003859163290000105
的正态分布，Qi为过程噪声wi的协方差矩阵，表示为
Figure GDA0003859163290000106
将公式(16)—(18)代入到公式(15)中，则联合对数似然函数进一步表示为公式(19)所示的形式：
Figure GDA0003859163290000107
其中，tr表示对于矩阵迹的计算；
在已知自初始到第k个工作循环的电池容量数据y1:k＝[y1,y2,…,yk]的前提下，对于在第j次迭代过程，退化模型公式(4)的未知参数估计向量表示为：
Figure GDA0003859163290000108
根据公式(16)—(19)的推导，以及最大期望算法的基本原理。对于第j+1次迭代过程，分为E步和M步，表示成公式(20)—(21)所示：
E步骤：计算
Figure GDA0003859163290000109
M步骤：计算
Figure GDA00038591632900001010
其中，
Figure GDA00038591632900001011
为Θ在
Figure GDA00038591632900001012
下的条件期望，
Figure GDA00038591632900001013
为Θ与
Figure GDA00038591632900001014
之间的概率关系，
Figure GDA00038591632900001015
表示后验概率密度下的期望算子，
Figure GDA00038591632900001016
为在第j+1次迭代过程，退化模型公式(4)的未知参数估计向量；
步骤三二、判断是否满足收敛条件
Figure GDA00038591632900001017
或者达到最大迭代次数，是，执行步骤四；否则，令j＝j+1，然后执行步骤三三；
通过迭代步骤三三至步骤三六直到满足收敛判据为止，从而获得退化模型的未知参数值；
通过不断进行E步骤或者M步骤的迭代，直到满足
Figure GDA0003859163290000111
或者达到最大迭代次数10次为止停止迭代，由此实现关于退化模型中未知参数的自适应估计；
其中ε为10-10。
步骤三三、构建后验条件下的联合对数似然函数；具体过程为：
公式(20)扩展地写成后验条件下的联合对数似然函数，如公式(22)所示的形式：
Figure GDA0003859163290000112
其中，
Figure GDA0003859163290000113
为第j次迭代过程电池容量退化速率和电池容量的协方差初始值，
Figure GDA0003859163290000114
为第j次迭代过程中过程噪声方差矩阵，
Figure GDA0003859163290000115
为第j次迭代过程中测量噪声方差，
Figure GDA0003859163290000116
为第j次迭代过程电池容量退化速率和电池容量状态变量的初始值；
在公式(22)中，
Figure GDA0003859163290000117
和
Figure GDA0003859163290000118
是基于可观测的测量向量(自初始到第k个循环为止所有的电池容量数据)y1:k＝[y1,y2,…,yk]的隐含状态向量x0:k＝[x0,x1,…,xk]的后验条件期望；
步骤三四、根据后验条件期望构建新的联合对数似然函数；具体过程为：
将公式(30)—(33)中的条件期望代入到公式(22)中，可以写成公式(34)所示的形式：
Figure GDA0003859163290000121
其中，
Figure GDA0003859163290000122
为初始状态向量x0在RTS后向平滑运算中的值，
Figure GDA0003859163290000123
为初始状态变量协方差矩阵S0在RTS后向平滑运算中的值；tr表示对于矩阵迹的计算；
步骤三五、计算一阶偏导数获得自适应参数估计解；具体过程为：
基于公式(22)—(34)中E步骤的计算结果，根据公式(22)计算最大期望算法的M步骤，具体过程为：
针对退化模型公式(4)的未知参数估计向量
Figure GDA0003859163290000124
中每一个参数，令每一个参数的偏导数值为0，求解出公式(35)—(38)的结果：
Figure GDA0003859163290000125
Figure GDA0003859163290000126
Figure GDA0003859163290000127
Figure GDA0003859163290000128
其中，
Figure GDA0003859163290000129
为第j+1次迭代过程电池容量退化速率和电池容量状态变量的初始值，
Figure GDA00038591632900001210
为第j+1次迭代过程电池容量退化速率和电池容量的协方差初始值，
Figure GDA00038591632900001211
为RTS后向平滑运算中第i-1个工作循环的协方差更新值，
Figure GDA00038591632900001212
为第j+1次迭代过程中测量噪声方差，
Figure GDA00038591632900001213
为第j+1次迭代过程中过程噪声方差矩阵，
Figure GDA00038591632900001214
Figure GDA00038591632900001215
为第j+1次迭代过程中电池容量退化速率噪声方差，
Figure GDA00038591632900001216
为第j+1次迭代过程中电池容量噪声方差；通过计算
Figure GDA0003859163290000131
可以获得
Figure GDA0003859163290000132
和
Figure GDA0003859163290000133
从而获得退化模型的未知参数值；执行步骤三二。
其它步骤及参数与具体实施方式一至六之一相同。
具体实施方式八：本实施方式与具体实施方式一至七之一不同的是，所述步骤三三中
Figure GDA0003859163290000134
和
Figure GDA0003859163290000135
的求解过程为：
为了计算公式(22)中的在后验条件下的状态变量期望，本发明采用Rauch-Tung-Striebel(RTS)最优平滑算法进行后向平滑实现。
基于卡尔曼滤波算法的前向迭代，获得状态向量均值和协方差的估计如公式(23)和公式(24)所示：
Figure GDA0003859163290000136
Figure GDA00038591632900001318
其中，
Figure GDA0003859163290000137
为RTS后向平滑状态向量的初始值，
Figure GDA0003859163290000138
为RTS后向平滑协方差初始值；
相应地，状态向量xk-1和xk之间的协方差可以表示成公式(25)所示的形式：
Figure GDA0003859163290000139
其中，
Figure GDA00038591632900001310
为状态向量xk-1和xk之间的协方差，cov(xk-1,xk|y1:k)为状态向量xk-1和xk之间的退化状态滤波协方差，Sk-1|k-1为k-1循环的状态向量xk-1的协方差矩阵；
根据公式(23)—(25)中的
Figure GDA00038591632900001311
和
Figure GDA00038591632900001312
即可以进行公式(26)-(29)中的RTS平滑运算；
RTS平滑增益的表达式如公式(26)所示：
Figure GDA00038591632900001313
其中，Gi为RTS平滑增益，Si|i为第i个循环的状态向量xi的协方差矩阵，Si+1|i为第i个循环的根据公式(11)计算的前向一步协方差矩阵；T为转置；
相应地，后向迭代过程中的状态向量
Figure GDA00038591632900001314
和协方差
Figure GDA00038591632900001315
的更新如公式(27)—(28)所示：
Figure GDA00038591632900001316
Figure GDA00038591632900001317
其中，
Figure GDA0003859163290000141
为RTS后向平滑运算中第i个工作循环的状态向量更新值，
Figure GDA0003859163290000142
为卡尔曼前向滤波运算中第i个工作循环的状态向量值，
Figure GDA0003859163290000143
为RTS后向平滑运算中第i+1个工作循环的状态向量更新值，
Figure GDA0003859163290000144
为卡尔曼前向滤波运算中第i个工作循环的状态向量前向一步预测值，
Figure GDA0003859163290000145
为RTS后向平滑运算中第i个工作循环的协方差更新值，
Figure GDA0003859163290000146
为RTS后向平滑运算中第i+1个工作循环的协方差更新值，Si+1|i为卡尔曼前向滤波运算中第i个工作循环的协方差向一步预测值；
第i-1和i循环之间状态向量xi-1和xi之间的协方差可以表示成公式(29)所示的形式：
Figure GDA0003859163290000147
其中，
Figure GDA0003859163290000148
为第i-1和i循环之间状态向量xi-1和xi之间的协方差，
Figure GDA0003859163290000149
为第i和i+1循环之间状态向量xi和xi+1之间的协方差；
基于公式(23)—(29)，求解后验下的条件期望表达式
Figure GDA00038591632900001410
Figure GDA00038591632900001411
和
Figure GDA00038591632900001412
如式(30)—(33)所示：
Figure GDA00038591632900001413
Figure GDA00038591632900001414
Figure GDA00038591632900001415
Figure GDA00038591632900001416
其中
Figure GDA00038591632900001417
为RTS后向平滑运算中第i-1个工作循环的状态向量更新值，
Figure GDA00038591632900001418
为第i和i+1循环之间状态向量xi-1和xi之间的协方差。
其它步骤及参数与具体实施方式一至七之一相同。
具体实施方式九：本实施方式与具体实施方式一至八之一不同的是，所述步骤四中基于步骤三计算电池的剩余使用寿命；具体过程为：
计算电池的剩余使用寿命：经过步骤三中对于退化模型参数的自适应估计，获得步骤2中状态变量的更新值。通过状态变量的取值构建剩余使用寿命的概率密度函数，进而获得电池的剩余使用寿命。
步骤四一、定义电池的剩余使用寿命表达式；具体过程为：
需要指出的是，锂离子电池的荷电状态(State of Charge，SOC)是衡量电池剩余容量的关键性指标。美国先进电池联合会(USABC)将电池的荷电状态SOC定义为电池当前的容量Qp与额定容量Qf之间的比值，计算公式如式(39)所示：
Figure GDA0003859163290000151
当SOC低于80％时，表明电池的性能难以达到设备或者系统的要求，此时需要进行电池的维护或者更换措施；在此基础上，本发明将电池的剩余使用寿命定义为从当前工作时间直到电池的SOC首次下降到80％为止所用的时长；
对于电池的第k个工作循环，表示为公式(40)的形式：
RULk＝inf{lk:Q(tk+lk)≥0.8Q0|Q0＞0} (40)
其中，RULk为剩余使用寿命lk的下确界，tk为电池从初始状态开始工作到第k个工作循环所经历的时间，Q0为电池在初始状态下的容量；
对于公式(1)中的线性维纳过程{q(t),t≥0}，如果漂移系数α服从正态分布，即
Figure GDA0003859163290000152
则q(t)首达电池失效阈值0.8Q0时间的概率密度函数fRUL(t)可以表示成公式(41)所示的形式：
Figure GDA0003859163290000153
其中，
Figure GDA0003859163290000154
表示服从均值为μα，方差为
Figure GDA0003859163290000155
的正态分布；μα为漂移系数α的均值；
Figure GDA0003859163290000156
为漂移系数α的方差；
基于可观测的测量向量(自初始到第k个循环为止所有的电池容量数据)y1:k＝[y1,y2,…,yk](y1:k＝[y1,y2,…,yk]代表了电池从第1个工作循环到第k个工作循环下的数据，这是可测的)估计出电池容量的漂移系数α的概率密度分布如公式(42)所示：
Figure GDA0003859163290000157
其中，E(αk|y1:k)为αk在动态状态空间方程下的期望，y1:k为自初始到第k个循环为止所有的可观测的测量向量(自初始到第k个循环为止所有的电池容量数据)；var(αk|y1:k)为电池容量退化速率的退化状态滤波方差；
步骤四二、计算剩余使用寿命的概率密度函数；具体过程为：
结合公式(40)—(42)，求解出剩余使用寿命lk的概率密度函数
Figure GDA0003859163290000161
如式(43)所示：
Figure GDA0003859163290000162
其中，
Figure GDA0003859163290000163
和
Figure GDA0003859163290000164
分别为第k个循环下电池容量的退化速率αk的方差和均值，
Figure GDA0003859163290000165
为在自初始到第k个循环为止所有的电池容量数据条件下的剩余使用寿命概率密度函数，
Figure GDA0003859163290000166
为在自初始到第k个循环为止所有的电池容量数据条件下的剩余使用寿命概率密度函数，yk为第k个循环下的电池容量；
步骤四三、计算剩余使用寿命的期望值；具体过程为：
电池的剩余使用寿命为针对公式(43)的数学期望，公式(43)的数学期望表达式如公式(44)所示：
Figure GDA0003859163290000167
其中，
Figure GDA0003859163290000168
为剩余使用寿命lk在αk,y1:k条件下的后验概率密度函数，p(αk|y1:k)为αk在y1:k条件下的后验概率密度函数，y1:k为为自初始到第k个循环为止所有的可观测的测量向量。
在完成步骤四三中计算剩余使用寿命的期望值后，判断第k次工作循环下的电池SOC是否低于80％，即计算电池容量yk和电池初始容量Q0之间的比值是否小于80％。如果SOC低于80％，则令k＝k+1，执行步骤一至步骤四，否则执行步骤五。
其它步骤及参数与具体实施方式一至八之一相同。
具体实施方式十：本实施方式与具体实施方式一至九之一不同的是，所述步骤五中基于步骤四评估剩余使用寿命的在线预测效果；具体过程为：
采用绝对平均误差(Mean Absolute Error，MAE)和均方根误差(Root MeanSquare Error，RMSE)这两个指标评估提出的自适应剩余使用寿命的预测效果；
绝对平均误差和均方根误差的表达式如式(45)—(46)所示：
Figure GDA0003859163290000171
Figure GDA0003859163290000172
其中，N为测试数据的样本数量，i为样本的序号，RULpi和RULti分别为第i个样本剩余使用寿命的预测值和真实值。
绝对平均误差和均方根误差的取值越小，则说明本发明提出的自适应剩余使用寿命的在线预测效果越好。
其它步骤及参数与具体实施方式一至九之一相同。
采用以下实施例验证本发明的有益效果：
实施例一：
本发明采用NASA Ames预测中心电提供的数据集验证所提出的电池剩余使用寿命的自适应预测方法。NASA Ames的预测中心提供了型号为18650的4个额定容量为2Ah的锂离子电池(编号为B0005,、B0006,、B0007和B0018)在室温(24℃)下工作所产生的数据。这4个电池在在充电、放电以及阻抗测试这3种模式下进行工作。对于充电模式，电池以1.5A的恒定电流进行充电，直到电池两端的电压达到4.2V，然后以恒定电压继续充电，直到充电电流下降到20mA。对于放电模式，电池以2A的恒定电流进行放电，直到4个电池的电压分别下降到2.7V、2.5V、2.2V和2.5V。阻抗测试则是对电池以频率从0.1Hz到5kHz的电化学阻抗谱(EIS)进行扫描，从而获得电池的内部参数。NASA电池数据集本身直接提供了每一个充放电循环结束时的电池容量，因此可以直接用于电池剩余使用寿命的动态预测。这4个电池的电池容量衰减曲线如图2所示。本发明对于上述的4个电池分别采用本发明提出的自适应的方法进行剩余使用寿命的预测，其具体步骤如下：
步骤1，提取电池工作过程中的容量数据：对NASA电池数据集中的电池过程中的容量数据进行提取，作为电池自适应剩余使用寿命预测方法的输入。
步骤2，构建电池的退化模型：仅仅依据电池容量数据，构建基于维纳过程的卡尔曼滤波动态退化模型描述电池的退化过程，通过对于状态变量和协方差的估计、卡尔曼滤波增益的计算以及状态变量和协方差的更新的步骤获得隐含的参数变量。
步骤3，自适应估计退化模型参数：采用基于最大期望-RTS的自适应地估计动态退化模型中的隐含参数。具体而言，针对基于维纳过程的卡尔曼滤波动态退化模型中所存在的未知参数，采用最大期望-RTS的方法进行自适应参数估计，利用电池从开始工作到当前工作循环的电池容量数据对退化模型中的状态变量进行递归更新。本发明以编号为B0005的电池为例，展示算法的自适应参数估计结果如图3a、3b、3c、3d、3e、3f、3g所示。从图中可以看出，随着电池容量数据的累积，所有的模型参数在20个工作循环以内即可收敛于特定的值，由此可知本发明的自适应参数估计方法有着很好的估计效果。
步骤4，计算电池的剩余使用寿命：经过步骤3中对于退化模型参数的自适应估计，获得步骤2中状态变量的更新值。通过状态变量的取值构建剩余使用寿命的概率密度函数，进而获得电池的剩余使用寿命。本发明以编号为B0005的电池为例，绘制了不同循环下剩余使用寿命的概率密度函数如图4所示。从图4中可以看出，随着电池容量数据的累积，剩余使用寿命的概率密度函数越来越尖锐，这说明由于越来越多的电池容量数据用于模型参数的估计，剩余使用寿命预测的不确定性得以降低，意味着本发明提出的方法具有很好的适应性。
步骤5，评估剩余使用寿命的在线预测效果：采用绝对平均误差(Mean AbsoluteError，MAE)和均方根误差(Root Mean Square Error，RMSE)这两个指标评估提出的自适应剩余使用寿命的预测效果，结果如表格1所示。图5a、5b、5c、5d展示了提出的方法在NASA电池数据集上4个电池的剩余使用寿命预测值和真实值。可以看出，本发明提出的方法能够有效自适应在线预测电池的剩余使用寿命，从而及时掌握电池的运行情况，协助操作人员及时对电池进行维护和更换，降低电力系统或者设备的运行成本。
表1自适应剩余使用寿命的预测结果
Figure GDA0003859163290000181
本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，本领域技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。

Claims (4)

1.一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法，其特征在于：所述方法具体过程为：
步骤一：提取电池第k次工作过程中的容量数据；
步骤二：基于步骤一构建电池的退化模型；
步骤三：基于步骤二构建自适应估计退化模型参数；
步骤四：基于步骤三计算电池的剩余使用寿命；
获得电池的剩余使用寿命后，判断第k次工作循环下的电池的荷电状态SOC是否低于80％，即计算电池当前的容量Qp与额定容量Qf之间的比值是否低于80％，如果SOC高于80％，则令k＝k+1，执行步骤一至步骤四，否则电池报废；
步骤五：基于步骤四评估剩余使用寿命的在线预测效果；
所述步骤二中基于步骤一构建电池的退化模型；具体过程为：
步骤二一、构建基于维纳过程的卡尔曼滤波动态退化模型；具体过程为：
将电池容量作为描述电池退化过程中的退化量，维纳过程{q(t),t≥0}表达成公式(1)所示的形式：
q(t)＝αt+βB(t) (1)
其中，q(t)为t时刻下的电池容量；α为漂移系数；β为扩散系数；B(t)为标准布朗运动，且满足B(t)～N(0,t)，N(0,t)表示服从均值为0，方差为t的正态分布；
公式(1)的退化过程具体地表示成公式(2)和(3)所示的形式：
Figure FDA0003833082910000011
qk＝qk-1+αk-1(tk-tk-1)+βζk (3)
其中，
Figure FDA0003833082910000012
为描述电池容量退化速率过程中的不确定性，
Figure FDA0003833082910000013
Figure FDA0003833082910000014
为表示服从均值为0，方差为
Figure FDA0003833082910000015
的正态分布，
Figure FDA0003833082910000016
为
Figure FDA0003833082910000017
的方差，αk为电池在第k次工作循环下的容量退化速率，αk-1为电池在第k-1次工作循环下的容量退化速率，qk为电池在第k次工作循环下的电池容量，qk-1为电池在第k-1次工作循环下的电池容量，tk为第k次工作循环下的时间，tk-1为第k-1次工作循环下的时间，ζk为扩散系数的噪声项，ζk～N(0,tk-tk-1)，N(0,tk-tk-1)为表示服从均值为0，方差为tk-tk-1的正态分布；
步骤二二、构建卡尔曼状态空间模型；具体过程为：
公式(2)和(3)可以采用公式(4)以描述电池的退化；
Figure FDA0003833082910000021
其中，xk为第k次工作循环下的状态向量，包含了电池容量退化速率和电池容量两个状态变量；xk-1为第k-1次工作循环下的状态向量，Ak为第k次工作循环下的状态空间方程系数；wk为过程噪声项；yk为电池在第k次工作循环下的实际容量，Ck为第k次工作循环下的输出方程yk＝Ckxk+ek的系数；ek为测量噪声项；
步骤二三、定义卡尔曼状态空间模型参数；具体过程为：
卡尔曼状态空间模型由：隐含状态向量x1:k＝[x1,x2,…,xk]、可观测的测量向量y1:k＝[y1,y2,…,yk]以及状态空间模型参数Ξ＝[Ak,wk,Ck,ek]所构成；
Figure FDA0003833082910000022
其中，
Figure FDA0003833082910000023
为退化状态滤波均值，
Figure FDA0003833082910000024
和
Figure FDA0003833082910000025
分别为电池容量退化速率和电池容量这两个状态变量的退化状态滤波均值；E(αk|y1:k)为αk在已知测量数据y1:k条件下的后验期望，E(xk|y1:k)为状态向量xk在已知测量数据y1:k条件下的后验期望，E(qk|y1:k)为qk在已知测量数据y1:k条件下的后验期望；
Figure FDA0003833082910000026
其中，Sk|k为退化状态滤波协方差，var(αk|y1:k)和var(qk|y1:k)分别为电池容量退化速率和电池容量两个状态变量的退化状态滤波方差，cov(αk,qk|y1:k)为电池容量退化速率和电池容量两个状态变量之间的退化状态滤波协方差，cov(xk|y1:k)为状态向量xk在已知测量数据y1:k条件下的协方差矩阵；
Figure FDA0003833082910000027
其中，
Figure FDA0003833082910000028
为退化状态一步预测均值，
Figure FDA0003833082910000029
和
Figure FDA00038330829100000210
分别为两个状态变量的退化状态一步预测均值，
Figure FDA00038330829100000211
为
Figure FDA00038330829100000212
在已知测量数据y1:k-1条件下的先验期望，E(αk|y1:k-1)为αk在已知测量数据y1:k-1条件下的先验期望，E(qk|y1:k-1)为qk在已知测量数据y1:k-1条件下的先验期望；
Figure FDA0003833082910000031
其中，Sk|k-1为退化状态一步预测协方差，var(αk|y1:k-1)和var(qk|y1:k-1)分别为电池容量退化速率和电池容量两个状态变量的退化状态一步预测方差，cov(αk,qk|y1:k-1)为电池容量退化速率和电池容量两个状态变量之间的退化状态滤波协方差，cov(xk|y1:k-1)为状态向量xk在已知测量数据y1:k-1条件下的协方差矩阵；
步骤二四、基于步骤二三估计状态与协方差；
步骤二五、基于步骤二四计算卡尔曼滤波增益；
步骤二六、基于步骤二五更新均值与协方差；
所述步骤二二中第k次工作循环下的状态向量xk、第k次工作循环下的状态空间方程系数Ak、过程噪声项wk、第k次工作循环下的输出方程yk＝Ckxk+ek的系数Ck、测量噪声项ek的表达式为：
Figure FDA0003833082910000032
Figure FDA0003833082910000033
Figure FDA0003833082910000034
Figure FDA0003833082910000035
ek＝v (5)
其中，wk的协方差矩阵表示为
Figure FDA0003833082910000036
v为输出方程yk＝Ckxk+ek的噪声项，噪声项
Figure FDA0003833082910000037
Figure FDA0003833082910000038
为测量噪声v的方差；
所述步骤二四中基于步骤二三估计状态与协方差；具体过程为：
状态估计：
Figure FDA0003833082910000041
其中，
Figure FDA0003833082910000042
为状态向量的前向一步预测结果；
协方差估计：
Figure FDA0003833082910000043
其中，Sk-1|k-1为状态向量xk在第k-1个工作循环的协方差矩阵；
所述步骤二五中基于步骤二四计算卡尔曼滤波增益；具体过程为：
计算卡尔曼滤波增益Kk：
Figure FDA0003833082910000044
其中，
Figure FDA0003833082910000045
为测量噪声的方差；
所述步骤二六中基于步骤二五更新均值与协方差；具体过程为：
状态更新：
Figure FDA0003833082910000046
协方差更新：
Sk|k＝Sk|k-1-KkCkSk|k-1 (14)
所述步骤三中基于步骤二自适应估计退化模型参数；具体过程为：
步骤三一、构建联合对数似然函数；具体过程为：
在已知自初始到第k个工作循环的电池容量数据y1:k＝[y1,y2,…,yk]的前提下，构建联合对数似然函数如公式(15)所示：
Figure FDA0003833082910000047
其中，Θ为公式(4)-(5)中退化方程的未知参数估计向量，
Figure FDA0003833082910000048
Figure FDA0003833082910000049
为均值的初始值，S0|0为协方差的初始值，p(x0:k,y1:k|Θ)为x0:k和y1:k在Θ下的联合概率密度函数，x0:k＝[x0,x1,...,xk]，xk为第k个工作循环的状态向量，p(y1:k|x0:k,Θ)为y1:k在x0:k,Θ下的条件概率密度，y1:k|x0:k,Θ为y1:k和x0:k,Θ之间的概率关系，p(x0:k|Θ)为x0:k的先验概率密度函数，y1:k|x0:k为y1:k和x0:k之间的概率关系，p(x0|Θ)为x0的先验概率密度函数，x0|Θ为x0和Θ之间的概率关系，p(xi|xi-1,Θ)为xi在xi-1,Θ下的条件概率密度，xi|xi-1为xi和xi-1之间的概率关系，p(yi|xi,Θ)为yi在xi,Θ下的条件概率密度，yi|xi为yi和xi之间的概率关系，x0为电池未工作时的状态向量，xi为第i个工作循环的状态向量，xi-1为第i-1个工作循环的状态向量，i＝1，2，…，k-1，k；
Figure FDA0003833082910000051
xi|xi-1,Θ～N(Aixi-1,Qi)(17)
Figure FDA0003833082910000052
其中，
Figure FDA0003833082910000053
表示服从均值为
Figure FDA0003833082910000054
方差为S0|0的正态分布，N(Aixi-1,Qi)表示服从均值为Aixi-1，方差为Qi的正态分布，
Figure FDA0003833082910000055
表示服从均值为Cixi，方差为
Figure FDA0003833082910000056
的正态分布，Qi为过程噪声wi的协方差矩阵，表示为
Figure FDA0003833082910000057
将公式(16)—(18)代入到公式(15)中，则联合对数似然函数进一步表示为公式(19)所示的形式：
Figure FDA0003833082910000058
其中，tr表示对于矩阵迹的计算；
在已知自初始到第k个工作循环的电池容量数据y1:k＝[y1,y2,…,yk]的前提下，对于在第j次迭代过程，退化模型公式(4)的未知参数估计向量表示为：
Figure FDA0003833082910000061
对于第j+1次迭代过程，分为E步和M步，表示成公式(20)—(21)所示：
E步骤：计算
Figure FDA0003833082910000062
M步骤：计算
Figure FDA0003833082910000063
其中，
Figure FDA0003833082910000064
为Θ在
Figure FDA0003833082910000065
下的条件期望，
Figure FDA0003833082910000066
为Θ与
Figure FDA0003833082910000067
之间的概率关系，
Figure FDA0003833082910000068
表示后验概率密度下的期望算子，
Figure FDA0003833082910000069
为在第j+1次迭代过程，退化模型公式(4)的未知参数估计向量；
步骤三二、判断是否满足收敛条件
Figure FDA00038330829100000610
或者达到最大迭代次数，是，执行步骤四；否则，令j＝j+1，然后执行步骤三三；
步骤三三、构建后验条件下的联合对数似然函数；具体过程为：
公式(20)扩展地写成后验条件下的联合对数似然函数，如公式(22)所示的形式：
Figure FDA00038330829100000611
其中，
Figure FDA00038330829100000612
为第j次迭代过程电池容量退化速率和电池容量的协方差初始值，
Figure FDA00038330829100000613
为第j次迭代过程中过程噪声方差矩阵，
Figure FDA00038330829100000614
为第j次迭代过程中测量噪声方差，
Figure FDA00038330829100000615
为第j次迭代过程电池容量退化速率和电池容量状态变量的初始值；
在公式(22)中，
Figure FDA0003833082910000071
和
Figure FDA0003833082910000072
是基于可观测的测量向量y1:k＝[y1,y2,…,yk]的隐含状态向量x0:k＝[x0,x1,…,xk]的后验条件期望；
步骤三四、根据后验条件期望构建新的联合对数似然函数；具体过程为：
将条件期望代入到公式(22)中，可以写成公式(34)所示的形式：
Figure FDA0003833082910000073
其中，
Figure FDA0003833082910000074
为初始状态向量x0在RTS后向平滑运算中的值，
Figure FDA0003833082910000075
为初始状态变量协方差矩阵S0在RTS后向平滑运算中的值；tr表示对于矩阵迹的计算；
步骤三五、计算一阶偏导数获得自适应参数估计解；具体过程为：
基于E步骤的计算结果，根据公式(22)计算最大期望算法的M步骤，具体过程为：
针对退化模型公式(4)的未知参数估计向量
Figure FDA0003833082910000076
中每一个参数，令每一个参数的偏导数值为0，求解出公式(35)—(38)的结果：
Figure FDA0003833082910000077
Figure FDA0003833082910000078
Figure FDA0003833082910000079
Figure FDA00038330829100000710
其中，
Figure FDA00038330829100000711
为第j+1次迭代过程电池容量退化速率和电池容量状态变量的初始值，
Figure FDA00038330829100000712
为第j+1次迭代过程电池容量退化速率和电池容量的协方差初始值，
Figure FDA00038330829100000713
为RTS后向平滑运算中第i-1个工作循环的协方差更新值，
Figure FDA0003833082910000081
为第j+1次迭代过程中测量噪声方差，
Figure FDA0003833082910000082
为第j+1次迭代过程中过程噪声方差矩阵，
Figure FDA0003833082910000083
Figure FDA0003833082910000084
为第j+1次迭代过程中电池容量退化速率噪声方差，
Figure FDA0003833082910000085
为第j+1次迭代过程中电池容量噪声方差；通过计算
Figure FDA0003833082910000086
可以获得
Figure FDA0003833082910000087
和
Figure FDA0003833082910000088
从而获得退化模型的未知参数值；执行步骤三二。
2.根据权利要求1所述一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法，其特征在于：所述步骤三三中
Figure FDA0003833082910000089
和
Figure FDA00038330829100000810
的求解过程为：
基于卡尔曼滤波算法的前向迭代，获得状态向量均值和协方差的估计如公式(23)和公式(24)所示：
Figure FDA00038330829100000811
Figure FDA00038330829100000812
其中，
Figure FDA00038330829100000813
为RTS后向平滑状态向量的初始值，
Figure FDA00038330829100000814
为RTS后向平滑协方差初始值；
相应地，状态向量xk-1和xk之间的协方差可以表示成公式(25)所示的形式：
Figure FDA00038330829100000815
其中，
Figure FDA00038330829100000816
为状态向量xk-1和xk之间的协方差，cov(xk-1,xk|y1:k)为状态向量xk-1和xk之间的退化状态滤波协方差，Sk-1|k-1为k-1循环的状态向量xk-1的协方差矩阵；
根据公式(23)—(25)中的
Figure FDA00038330829100000817
和
Figure FDA00038330829100000818
即可以进行公式(26)-(29)中的RTS平滑运算；
RTS平滑增益的表达式如公式(26)所示：
Figure FDA00038330829100000819
其中，Gi为RTS平滑增益，Si|i为第i个循环的状态向量xi的协方差矩阵，Si+1|i为第i个循环的根据公式(11)计算的前向一步协方差矩阵；T为转置；
相应地，后向迭代过程中的状态向量
Figure FDA00038330829100000820
和协方差
Figure FDA00038330829100000821
的更新如公式(27)—(28)所示：
Figure FDA0003833082910000091
Figure FDA0003833082910000092
其中，
Figure FDA0003833082910000093
为RTS后向平滑运算中第i个工作循环的状态向量更新值，
Figure FDA0003833082910000094
为卡尔曼前向滤波运算中第i个工作循环的状态向量值，
Figure FDA0003833082910000095
为RTS后向平滑运算中第i+1个工作循环的状态向量更新值，
Figure FDA0003833082910000096
为卡尔曼前向滤波运算中第i个工作循环的状态向量前向一步预测值，
Figure FDA0003833082910000097
为RTS后向平滑运算中第i个工作循环的协方差更新值，
Figure FDA0003833082910000098
为RTS后向平滑运算中第i+1个工作循环的协方差更新值，Si+1|i为卡尔曼前向滤波运算中第i个工作循环的协方差向一步预测值；
第i-1和i循环之间状态向量xi-1和xi之间的协方差可以表示成公式(29)所示的形式：
Figure FDA0003833082910000099
其中，
Figure FDA00038330829100000910
为第i-1和i循环之间状态向量xi-1和xi之间的协方差，
Figure FDA00038330829100000911
为第i和i+1循环之间状态向量xi和xi+1之间的协方差；
基于公式(23)—(29)，求解后验下的条件期望表达式
Figure FDA00038330829100000912
Figure FDA00038330829100000913
和
Figure FDA00038330829100000914
如式(30)—(33)所示：
Figure FDA00038330829100000915
Figure FDA00038330829100000916
Figure FDA00038330829100000917
Figure FDA00038330829100000918
其中
Figure FDA00038330829100000919
为RTS后向平滑运算中第i-1个工作循环的状态向量更新值，
Figure FDA00038330829100000920
为第i和i+1循环之间状态向量xi-1和xi之间的协方差。
3.根据权利要求2所述一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法，其特征在于：所述步骤四中基于步骤三计算电池的剩余使用寿命；具体过程为：
步骤四一、定义电池的剩余使用寿命表达式；具体过程为：
将电池的荷电状态SOC定义为电池当前的容量Qp与额定容量Qf之间的比值，计算公式如式(39)所示：
Figure FDA0003833082910000101
将电池的剩余使用寿命定义为从当前工作时间直到电池的SOC首次下降到80％为止所用的时长；
对于电池的第k个工作循环，表示为公式(40)的形式：
RULk＝inf{lk:Q(tk+lk)＜0.8Q0|Q0＞0} (40)
其中，RULk为剩余使用寿命lk的下确界，tk为电池从初始状态开始工作到第k个工作循环所经历的时间，Q0为电池在初始状态下的容量；
对于公式(1)中的线性维纳过程{q(t),t≥0}，如果漂移系数α服从正态分布，即
Figure FDA0003833082910000102
则q(t)首达电池失效阈值0.8Q0时间的概率密度函数fRUL(t)可以表示成公式(41)所示的形式：
Figure FDA0003833082910000103
其中，
Figure FDA0003833082910000104
表示服从均值为μα，方差为
Figure FDA0003833082910000105
的正态分布；μα为漂移系数α的均值；
Figure FDA0003833082910000106
为漂移系数α的方差；
基于可观测的测量向量y1:k＝[y1,y2,…,yk]估计出电池容量的漂移系数α的概率密度分布如公式(42)所示：
Figure FDA0003833082910000107
其中，E(αk|y1:k)为αk在动态状态空间方程下的期望，y1:k为自初始到第k个循环为止所有的可观测的测量向量；var(αk|y1:k)为电池容量退化速率的退化状态滤波方差；
步骤四二、计算剩余使用寿命的概率密度函数；具体过程为：
结合公式(40)—(42)，求解出剩余使用寿命lk的概率密度函数
Figure FDA0003833082910000108
如式(43)所示：
Figure FDA0003833082910000111
其中，
Figure FDA0003833082910000112
和
Figure FDA0003833082910000113
分别为第k个循环下电池容量的退化速率αk的方差和均值，
Figure FDA0003833082910000114
为在自初始到第k个循环为止所有的电池容量数据条件下的剩余使用寿命概率密度函数，
Figure FDA0003833082910000115
为在自初始到第k个循环为止所有的电池容量数据条件下的剩余使用寿命概率密度函数，yk为第k个循环下的电池容量；
步骤四三、计算剩余使用寿命的期望值；具体过程为：
电池的剩余使用寿命为针对公式(43)的数学期望，公式(43)的数学期望表达式如公式(44)所示：
Figure FDA0003833082910000116
其中，
Figure FDA0003833082910000117
为剩余使用寿命lk在αk,y1:k条件下的后验概率密度函数，p(αk|y1:k)为αk在y1:k条件下的后验概率密度函数，y1:k为为自初始到第k个循环为止所有的可观测的测量向量。
4.根据权利要求3所述一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法，其特征在于：所述步骤五中基于步骤四评估剩余使用寿命的在线预测效果；具体过程为：
采用绝对平均误差和均方根误差评估剩余使用寿命的预测效果；
绝对平均误差和均方根误差的表达式如式(45)—(46)所示：
Figure FDA0003833082910000118
Figure FDA0003833082910000119
其中，N为测试数据的样本数量，i为样本的序号，RULpi和RULti分别为第i个样本剩余使用寿命的预测值和真实值。
CN202210065461.1A 2022-01-20 2022-01-20 一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法 Active CN114417686B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210065461.1A CN114417686B (zh) 2022-01-20 2022-01-20 一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210065461.1A CN114417686B (zh) 2022-01-20 2022-01-20 一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法

Publications (2)

Publication Number Publication Date
CN114417686A CN114417686A (zh) 2022-04-29
CN114417686B true CN114417686B (zh) 2023-02-03

Family

ID=81275678

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210065461.1A Active CN114417686B (zh) 2022-01-20 2022-01-20 一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法

Country Status (1)

Country Link
CN (1) CN114417686B (zh)

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2006057469A1 (en) * 2004-11-29 2006-06-01 Lg Chem, Ltd. Method and system for joint battery stateand parameter estimation
CN103399279A (zh) * 2013-08-01 2013-11-20 哈尔滨工业大学 基于ekf方法和ar模型融合型锂离子电池循环寿命预测方法
CN109932656A (zh) * 2019-04-17 2019-06-25 合肥工业大学 一种基于imm-upf的锂电池寿命估计方法
CN112800616A (zh) * 2021-02-05 2021-05-14 中国人民解放军空军工程大学 基于比例加速退化建模的设备剩余寿命自适应预测方法
CN112906213A (zh) * 2021-02-05 2021-06-04 中国人民解放军空军工程大学 一种机载电子设备剩余寿命自适应预测方法
CN112949026A (zh) * 2021-01-19 2021-06-11 中国人民解放军火箭军工程大学 一种考虑年龄和状态依赖的退化设备剩余寿命预测方法
CN113392507A (zh) * 2021-05-25 2021-09-14 西安科技大学 一种锂离子动力电池剩余寿命预测方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10664562B2 (en) * 2013-02-24 2020-05-26 Fairchild Semiconductor Corporation and University of Connecticut Battery state of charge tracking, equivalent circuit selection and benchmarking
CN112505568B (zh) * 2020-11-11 2022-03-15 电子科技大学 一种多电堆固体氧化物燃料电池系统寿命预测方法
CN112949068B (zh) * 2021-03-05 2022-04-19 北京航空航天大学 基于核密度估计的锂离子电池循环寿命预测方法

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2006057469A1 (en) * 2004-11-29 2006-06-01 Lg Chem, Ltd. Method and system for joint battery stateand parameter estimation
CN103399279A (zh) * 2013-08-01 2013-11-20 哈尔滨工业大学 基于ekf方法和ar模型融合型锂离子电池循环寿命预测方法
CN109932656A (zh) * 2019-04-17 2019-06-25 合肥工业大学 一种基于imm-upf的锂电池寿命估计方法
CN112949026A (zh) * 2021-01-19 2021-06-11 中国人民解放军火箭军工程大学 一种考虑年龄和状态依赖的退化设备剩余寿命预测方法
CN112800616A (zh) * 2021-02-05 2021-05-14 中国人民解放军空军工程大学 基于比例加速退化建模的设备剩余寿命自适应预测方法
CN112906213A (zh) * 2021-02-05 2021-06-04 中国人民解放军空军工程大学 一种机载电子设备剩余寿命自适应预测方法
CN113392507A (zh) * 2021-05-25 2021-09-14 西安科技大学 一种锂离子动力电池剩余寿命预测方法

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
Data-Driven Battery Health Prognosis Using Adaptive Brownian Motion Model;Guangzhong Dong , Member, IEEE, Fangfang Yang , Zhongbao Wei;《IEEE TRANSACTIONS ON INDUSTRIAL INFORMATICS》;20200707;全文 *
State-of-Health Estimation for Lithium-Ion;XIAODONG XU1, CHUANQIANG YU1, SHENGJIN TANG 1, XIAOYAN SUN1;《IEEE TRANSACTIONS ON INDUSTRIAL INFORMATICS》;20190814;全文 *
基于粒子滤波-ARIMA的多层次锂离子电池寿命预测;周亚鹏;《中国优秀博士学位论文全文数据库(信息科技辑)》;20220115;24-35页 *
基于自适应维纳过程的剩余寿命预测方法研究.黄泽毅.《中国优秀硕士学位论文全文数据库(信息科技辑)》.2015,第6页、第13至24页. *
李玥锌 ; 刘淑杰 ; 高斯博 ; 胡娅维 ; 张洪潮.基于维纳过程的锂离子电池剩余寿命预测.《工程科学学报》.2017,全文. *
锂离子电池退化建模与剩余寿命预测方法研究;余斌;《中国优秀博士学位论文全文数据库(信息科技辑)》;20220115;59-35页 *
黄泽毅;基于自适应维纳过程的剩余寿命预测方法研究;《中国优秀硕士学位论文全文数据库(信息科技辑)》;20151215;第6页、第13至24页 *

Also Published As

Publication number Publication date
CN114417686A (zh) 2022-04-29

Similar Documents

Publication Publication Date Title
Li et al. Lithium-ion battery state of health monitoring based on ensemble learning
Jiang et al. State of health estimation for lithium-ion battery using empirical degradation and error compensation models
CN110187290B (zh) 一种基于融合型算法的锂离子电池剩余寿命预测方法
Li et al. A new prognostics method for state of health estimation of lithium-ion batteries based on a mixture of Gaussian process models and particle filter
Fan et al. A novel machine learning method based approach for Li-ion battery prognostic and health management
CN113406521B (zh) 一种基于特征分析的锂电池健康状态在线估计方法
CN109738810A (zh) 一种电池剩余使用寿命的综合估算方法
CN112966432A (zh) 一种锂离子电池剩余有效寿命预测方法及装置
CN109633470B (zh) 基于ekf-gpr和日常片段数据的电池实时全充时间的估算方法
CN113125987A (zh) 一种新型混合的锂离子电池健康状态预测方法
Wang et al. Health diagnosis for lithium-ion battery by combining partial incremental capacity and deep belief network during insufficient discharge profile
CN112068003A (zh) 基于线性维纳过程的镉镍蓄电池寿命预测方法和装置
Sandoval-Chileño et al. Robust State of Charge estimation for Li-ion batteries based on Extended State Observers
CN113392507A (zh) 一种锂离子动力电池剩余寿命预测方法
CN114740388A (zh) 一种基于改进tcn的锂电池剩余寿命状态评估方法
CN112098874A (zh) 一种考虑老化情况下的锂离子电池电量预测方法
CN115201686A (zh) 一种不完备充放电数据下的锂离子电池健康状态评估方法
Audin et al. Auto-encoder LSTM for Li-ion SOH prediction: A comparative study on various benchmark datasets
CN116401770A (zh) 基于电池数字孪生模型与机器学习的快充策略设计方法
Bak et al. Accurate estimation of battery SOH and RUL based on a progressive lstm with a time compensated entropy index
CN114417686B (zh) 一种面向单个锂离子电池的自适应在线剩余使用寿命预测方法
Sahajpal et al. Accurate long-term prognostics of proton exchange membrane fuel cells using recurrent and convolutional neural networks
Zhang et al. Remaining useful life prediction of lithium-ion batteries based on TCN-DCN fusion model combined with IRRS filtering
Tao et al. A support vector regression-based prognostic method for li-ion batteries working in variable operating states
Basia et al. Comparison of data driven algorithms for SoH estimation of Lithium-ion batteries

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant