王露霄, 段建東, 張凱, 趙克, 孫力
(1.哈爾濱工業大學 電氣工程及自動化學院,黑龍江 哈爾濱 150001; 2.國網山西省電力公司電力科學研究院,山西 太原 030000)
超級電容是一種功率型儲能器件,與電池等比能量型的器件相比,其等效串聯電阻低、循環效率高、循環壽命長、且其功率密度大,適用于較寬的溫度范圍。超級電容能實現快速充放電,能夠消納尖峰脈沖,因此而廣泛應用于動態變化的儲能系統中[1-3]。為使儲能系統可靠運行,可以采用性能優良的超級電容管理系統對超級電容實行實時監測,防止其過充和過放,同時對其老化狀態進行評估,同時對超級電容的荷電狀態(state of charge,SOC)和健康狀態(state of health,SOH)進行實時在線估計,對儲存能量進行實時監測并在異常情況下及時進行報警保護動作[4-9]。
目前,SOC估計的一種基本方法是安時計量法,是一種開環估計方法,無法及時對原始誤差的積累進行修正,估計精度不高,需要利用開環電壓(open circuit voltage,OCV)-SOC關系來定期糾正電流積分數值計算的累積誤差[10]。超級電容SOC估計的觀測器方法主要有Luenberger觀測器[11-12]、PI觀測器[13]、擴展卡爾曼濾波器[14](extended Kalman Filter,EKF)、無跡卡爾曼濾波器[15](unscented Kalman filter,UKF)、滑模觀測器[16](sliding mode observer,SMO)、人工智能算法等。其中,Luenberger觀測器是經典線性系統觀測器,需要將狀態方程化為線性模型,并考慮狀態觀測器收斂性受非線性導致的模型誤差影響的大小,在實際應用中很難處理;PI觀測器也是一種線性狀態觀測器,與龍伯格觀測器類似,存在控制器參數對于模型誤差的魯棒性問題,同樣較難應用于強非線性系統;卡爾曼濾波器能實現狀態變量的最小均方誤差估計,也能在迭代過程中很好的抑制傳感器噪聲,擴展卡爾曼濾波器能夠很好地處理模型的非線性,是一種次優估計,而無跡卡爾曼濾波器對等效模型準確度敏感且容易受到不確定噪聲干擾,在充滿干擾的超級電容應用場合受限;由于滑模觀測器能夠很好地抑制模型誤差以及輸入擾動,在進行非線性模型的SOC的估計時有較好的效果,但滑模觀測器存在抖振問題,會使其狀態估計的收斂速度變慢。
現有超級電容SOC估計的文獻中,為了使超級電容等效電路模型能夠精確地描述超級電容實際端口特性,常采用三支路RC模型,其每條支路描述動態特性的時間常數都有很大跨度,從幾秒到幾千秒[17-18]。但由于采用三支路RC等效電路進行參數辨識時,存在多組局部最優解,使得超級電容快速參數辨識變得困難[19-21]。
另外,以上方法并沒有考慮對模型參數進行實時觀測。為此利用三支路RC模型簡化得到的二支路RC模型進行狀態估計,在保證模型精度的情況下采用EKF實現超級電容的非線性多參數辨識,利用參數辨識的結果來實現超級電容SOC和SOH的實時估計。
超級電容三支路RC模型可以很好地描述其非線性特性,例如超級電容的電荷再分配現象和自放電現象,即充電后靜置,超級電容的端電壓會下降,放電后靜置,超級電容端電壓會上升的現象。其3個支路因為時間常數的差異而代表了超級電容充放電過程的3個不同階段[22]。因此如圖1所示的三支路RC等效電路模型以其精確性和實用性被廣泛應用于超級電容SOC估計。

圖1 超級電容三支路RC等效電路模型Fig.1 Three branch RC model of supercapacitor
但由于將三支路RC模型用于狀態估計時,待確定的變量數為8個,分別為R1、C10、C11、R2、C2、R3、C3、Rleak,為了降低狀態估計的難度,將代表長時間充電的R3、C3和Rleak支路去除。又由于在大電流充電時,超級電容的充電速度很快,其特性主要與第一、第二支路有關,于是將三支路RC模型降階為二支路RC模型有理論可能。其簡化電路模型如圖2所示。

圖2 超級電容二支路RC等效電路模型Fig.2 Two branch RC model of supercapacitor
如圖1所示,超級電容三支路RC模型輸入為u=i,輸出為y=ut。根據基爾霍夫定律可列寫電壓電流方程為:
(1)
(2)
其中:τ1=(C10+C11u1)R1;τ2=R2C2;τ3=R3C3;R‖=R1‖R2‖R3‖Rleak,“‖”表示并聯。電壓電流方程經化簡可得三支路RC模型狀態方程為:

(3)
其中:
(4)
二支路RC模型為三支路RC模型去除代表大時間尺度的R3、C3和Rleak后的一種特殊情形,同理可得其狀態方程為:
(5)
其中:τ1=(C10+C11u1)R1;τ2=R2C2;Rp=R1‖R2。
卡爾曼濾波算法能很好地處理傳感器噪聲,為了能夠對非線性系統的參數進行準確的識別,可以采用EKF構造觀測器。EKF算法可以用于多種場合,對狀態方程中的恒定參數進行實時辨識便是其應用之一[23]。
系統的離散方程可表達如下:
(6)
假設過程噪聲wk和觀測噪聲vk為相互獨立且均值為0的高斯白噪聲,則具體的EKF迭代方程為:
(7)
其中Fk和Hk為雅可比矩陣,且Fk和Hk為對稱陣,有:
(8)

待定參數也可以作為狀態加入到狀態方程中,從而采用EKF對待定參數進行實時觀測[24]。模型參數以及模型狀態,即各支路電容電壓的辨識都可通過EKF實現,因此可以構造EKF算法同時對簡化模型的模型參數與各支路電容電壓進行在線辨識。

(9)
采用前向歐拉法將式(9)化為離散形式,時間間隔為Ts。為了簡便,省去最后結果的所有指標k,為了同時對C10、C11、R2、C2這4個參數進行在線辨識,也需將其加入到狀態方程中,且帶有時間指標k,而R1當作已知,有
(10)
將4個常值使用前向歐拉法離散可得:
(11)
其中w1(k)、w2(k)、w3(k)、w4(k)為人為加入的過程噪聲,其噪聲方差分別為S1k、S2k、S3k、S4k。

R1i(k)+v(k)。
(12)
由式(10)可得到相應的雅可比矩陣,雅可比矩陣主體如下,其中為了方便,采用轉置,且結果中省略了指標k,有
(13)
雅可比矩陣剩余部分有以下結果,式中I4指4階單位陣,有
(14)
EKF中使用的輸出方程Hk如下。這里直接對其做線性近似,無需再用雅可比矩陣。
(15)
根據以上算法即可實現對超級電容的多個模型參數的統一辨識。
算法有效性可通過卡爾曼濾波的自適應系統辨識理論加以驗證[25]。為直觀分析,將一些量重新定義:
(16)
則可以得到非線性方程為:
(17)

(18)

超級電容SOC的計算有兩種形式,一種是用當前可用電荷與最大可用電荷的比值求SOC,另一種方式是用當前儲存能量與最大儲存能量的比值求SOC,基本計算公式為:
(19)
其中:Q和E分別表示當前超級電容存儲的電荷和能量;Qrate和Erate分別表示超級電容額定的儲存電荷和儲存能量。目前依據能量來計算超級電容SOC的方法被大家廣泛運用,故也采用此種計算方式,計算公式如下:
(20)
其中:C表示超級電容的額定容值;urate為超級電容額定電壓。
超級電容在老化過程中的性能變差主要體現在容值的逐漸減小以及等效串聯電阻(equivalent series resistance,ESR)的逐漸增加,超級電容儲能能力會隨容值減小而降低,其損耗會隨著ESR增大而增大,充放電功率隨著ESR的增大而降低,衡量超級電容老化程度用這兩個實時參數也是公認的方法[25]。SOH的計算方法如下:
(21)
式中:R和C指當前的ESR和容值;下標“E”指超級電容壽命結束時的值;下標“new”指剛出廠時的值。根據Maxwell的常用標準:CE=0.8Cnew,RE=2Rnew。剛出廠時,SOH=100%,但SOH值會隨著超級電容的老化而逐漸降低,當超級電容壽命結束時SOH=0[26]。
在估計超級電容的SOC時,考慮到其電荷再分布現象和估計精確性,一般采用三支路RC模型,但對于SOH的估計則大多采用單支路RC模型,將SOH的估計轉化為對ESR和容值的測量或識別[27]。超級電容的等效串聯電阻R1很小,以表1中50 F超級電容為例,R1僅為16 mΩ,并且隨著超級電容容值的增加,其等效串聯電阻會減小,能減小到1 mΩ以下。使用EKF在線識別算法也未必會比采用ΔU/ΔI方法測量的精度高,且EKF中加入多余的變量會加大參數辨識的難度,所以僅需用EKF對超級電容容值進行識別,使用ΔU/ΔI識別R1。因此對超級電容的SOH進行估計時,僅需使用辨識出的參數C10、C11即可。
為了驗證將三支路RC模型簡化為二支路RC模型進行超級電容狀態估計其估計精度仍然滿足要求。使用Simulink中的Parameter Estimation功能來分別識別超級電容三支路RC模型與二支路RC模型參數。為確定最優參數,給定一個輸入電流,通過迭代算法來計算模型參數值,當實際超級電容輸出響應與仿真模型輸出響應均方誤差最小時,對應模型參數即為最優參數。
采用Parameter Estimation功能對圖3中所示激勵電流和該電流激勵得到的端口電壓數據進行擬合,得到如表1所示的模型參數。由于Rleak跟超級電容充電后的漏電流有關,代表長時間的過程,在參數擬合時其影響可忽略不計,故可把Rleak設為常數以代表已知參數。其余的7個模型參數的擬合結果均列于表1。

圖3 激勵電流Fig.3 Exciting current

表1 三支路RC模型的原始參數與使用Simulink的參數識別結果
在圖3所示的激勵電流激勵下三支路RC等效電路模型的響應得到的擬合參數及其與原始參數的誤差如圖4所示。
從圖4可知,原始參數曲線與擬合參數曲線幾乎重合在一起,在整個擬合過程中兩者誤差不超過±1 mV,擬合效果很好,但是從表1中數據看,只有第一支路的參數識別效果比較理想。由于無法對所有的電壓電流數據進行識別,表1中所呈現的只是在特定電壓電流數據下的辨識結果,在其它電流激勵下電壓響應曲線擬合效果不得而知。故在實際測試中,對于不確定的激勵,很難判斷其擬合時間和擬合效果。離線識別算法難以解決多參數識別問題,在線識別多參數則更為困難,故為了實現參數的精確識別,應適當降低模型參數。

圖4 原始參數與擬合參數對應超級電容模型的端口電壓及二者之差Fig.4 Port voltage of the supercapacitor model corresponding to the original parameters and the fitting parameters and difference between them
因此,將三支路RC模型進行簡化,對簡化得到的二支路RC模型施以圖3所示的電流激勵,通過其電壓響應和電流數據來進行參數擬合,從而對二支路RC模型的5個參數R1、C10、C11、R2、C2進行識別,識別結果如圖5和表2所示。從圖5可知,原始模型和二支路RC模型之間的電壓響應之間的誤差在整個擬合時間段內也一直未超過±1 mV,達到了擬合精度的要求,且有效減少了待定參數數目,降低了等效模型的階數。表2中列出了二支路模型的模型參數迭代初值與最優擬合結果。

表2 二支路模型的參數迭代初值與最終的最優參數

圖5 二支路模型的擬合電壓響應及誤差Fig.5 Fitting voltage response and error of two branch model
盡管迭代計算初值與實際參數不同,但只用了6次迭代就得到了參數的最優擬合結果。可見由于等效模型階數的降低,參數辨識問題的難度也會降低,參數也更容易收斂到最優參數組。由此可見,盡管降低模型階數在一定程度上降低了參數識別的精度,但其擬合誤差依然很小,能夠滿足參數識別的精度要求,且可以降低參數辨識的難度,因此二支路RC模型可以用來作為建模對象。
在仿真電路中對二支路RC模型采用恒流充電,識別其模型參數,具體參數為R1=10 mΩ、C10=100 F、C11=50 F/V,識別算法認為R1恒定,電壓電流采樣時間設置為Ts=0.1 s。在電壓數據中人為疊加高斯白噪聲序列以模擬傳感器噪聲,其均值為0、方差為0.01,式(8)中Rk=0.01。仿真中使用
(22)
仿真時采用連續的幅值為1 A的充放電方波電流為超級電容充電,則可得到圖6所示的仿真波形。

圖6 同時識別C10和C11時的仿真結果Fig.6 Simulation results of simultaneous recognition of C10 and C11
從圖中可以看出,當模型參數初值分別為C10=20 F和C11=20 F/V時,經過EKF迭代,模型參數可收斂到其真值。從放大后的參數辨識波形可以看出C10和C11的辨識誤差均在±5 F和±5 F/V之內。仿真中將R1設為常數,代表已知量,其值在恒流充電中由ΔU/ΔI獲取。仿真也顯示對于C10和C11兩個參數的識別并不需要準確的R1,由于R1的值太小,ESR在充放電過程中引起的壓降很小,可以忽略。
將R1也作為待辨識參數加入到狀態方程中,其過程噪聲與C10和C11的過程噪聲相互獨立,探究其對參數辨識的影響。為此,采用與圖6一樣的激勵電流作為輸入,對支路參數進行在線辨識,其中C10和C11的初值與上一仿真相同,R1的初值取值為1 Ω。辨識結果如圖7所示,由圖可知模型參數最終都能收斂到真值。由放大后的結果可知,增加識別R1后,C10和C11的穩態相對誤差比未增加識別R1時的誤差更小,而R1辨識的穩態相對誤差則達到了50%,其辨識結果不準確,幾乎不可用。

圖7 同時識別C10、C11及R1的仿真結果Fig.7 Simulation results of simultaneous recognition of C10、C11 and R1
兩個電壓狀態與模型參數的真值與迭代初值如表3所示。

表3 二支路RC模型參數與狀態的真實值與迭代計算初值
仿真中采樣時間間隔同樣也為Ts=0.1 s,使用EKF實現二支路模型的狀態估計與模型參數的聯合辨識的結果如圖8所示。

圖8 使用EKF對二支路模型的模型參數辨識與估計SOC的仿真結果Fig.8 Simulation results of using EKF algorithm to recognize model parameters and SOC of the two RC branch model of supercapacitor
圖8(b)為幅值為±0.5 A、周期接近300 s的正負交變方波激勵電流。圖8(a)顯示第一支路電容參數C10和C11都能很快收斂到各自的真實值,第二支路中C2的辨識結果也能收斂到其真值,而R2的辨識結果最終仍存在20%的誤差。根據圖8(b)顯示的結果,EKF能很好地追蹤兩個支路電壓u1和u2,從而能夠準確地對超級電容SOC進行實時估計。在超級電容充放電電流快速變化的情況下,EKF也能夠準確地進行狀態觀測和參數辨識。仿真結果表明,對于非線性模型的模型參數辨識和狀態估計,EKF能夠快速響應并實現對相關參數的準確辨識,當實際模型容易出現參數變化時,可以采用EKF來辨識模型參數,進而對運行結果進行一定的矯正。
為了驗證將超級電容三支路RC等效電路簡化為二支路RC等效電路進行模型參數辨識的準確性和快速性,采用表3所示的模型參數計算初值進行迭代,并采用恒定電流源為超級電容充電,超級電容放電時放電端接0.4 Ω的放電電阻進行放電,整個實驗平臺如圖9所示。

圖9 超級電容多參數在線辨識實驗平臺Fig.9 Experimental platform for multiparameter online identification of supercapacitor
實驗中,采用2 A的充電電流為超級電容進行充電,當超級電容的電壓達到1.5 V時將超級電容切換到0.4 Ω的放電電阻上進行放電,當超級電容的電壓降到0.3 V左右時,又重復上述過程,以驗證所提方法的準確性。實驗波形如圖10所示。
圖10中,i代表充電電流,ut為超級電容實際電壓,u1、u2分別為超級電容二支路RC模型的支路電壓,ue為實際超級電容電壓與第一支路電壓估計值之間的誤差。由圖10(a)波形可知,簡化超級電容模型第一支路能很迅速地追蹤超級電容電壓變化,且二者之間的誤差不超過0.1 V。在圖10(b)中,由于二支路模型的第二條支路時間常數較大,主要用來描述超級電容較長時間充電的充電特性,在恒流充電實驗中,如果能給予足夠長的充電時間,其最終可充電至最大電壓,但在快速充放電實驗中,沒有足夠時間使其達到最大充電電壓,于是第二條支路無法像第一支路一樣快速跟蹤實際電壓變化,只能跟蹤實際電壓的變化趨勢。

圖10 2 A充電電流下超級電容實際電壓與二支路電壓波形Fig.10 Actual voltage and two branch voltage waveforms of supercapacitor under 2 A charging current
在實現二支路RC模型準確估計超級電容實際電壓的同時,也實現了對超級電容模型參數的辨識,辨識結果如圖11所示。
由圖11可知,只有C10從初值20 F迅速收斂到了其真實值42.7 F左右;其余3個參數C11、C2、R2均收斂到接近于實際值的某一數值并維持不變,而均未收斂到表3所列真實值。究其原因,在于表3所列數值為基于三支路RC等效電路模型所測參數,當直接簡化為二支路模型時,忽略了第三支路和和泄漏電流支路分流的影響,其結果必然與三支路模型所測參數之間存在誤差。總之,EKF是一種很好的參數識別方法,在實際應用中可用來修正狀態方程以改善狀態估計結果。

圖11 二支路RC模型的參數辨識結果Fig.11 Identification results of two branch RC model parameters
本文采用EKF對超級電容二支路RC模型的模型參數進行在線識別,并依據參數辨識結果對超級電容的SOC與SOH進行實時計算,得出了如下結論:
以二支路RC等效電路模型為基礎進行超級電容狀態估計和參數辨識,在保證辨識精度的情況下,由于模型階數的降低,實現了較三支路RC模型更快的估計和辨識速度,同時也有效避免了高階系統存在多組局部最優解的情況。
采用文中所述的EKF算法,可識別超級電容等效電路的模型參數,對模型參數實時監控。對于受模型參數影響較大的算法來說,可以利用觀測到的模型參數將其狀態方程實時更新,從而提高辨識精度。
在估計超級電容的SOH時,由于其ESR很小,采用EKF算法也未必能對ESR進行有效估計,故其計算采用階躍響應的ΔU/ΔI即可,其余兩個參數的辨識對ESR的依賴程度很低,故在參數辨識時可將ESR設為常數以降低辨識難度。