王子赟，張 帥，王 艷，劉子幸，紀(jì)志成

(1.江南大學(xué)物聯(lián)網(wǎng)應(yīng)用技術(shù)教育部工程研究中心，江蘇無錫 214122;2.江南大學(xué)輕工過程先進(jìn)控制教育部重點(diǎn)實(shí)驗(yàn)室，江蘇無錫 214122)

1 引言

系統(tǒng)辨識是控制領(lǐng)域?qū)W者們持續(xù)關(guān)注的一個(gè)熱點(diǎn)問題[1-3].隨著科技的發(fā)展和系統(tǒng)應(yīng)用領(lǐng)域的延伸，系統(tǒng)的復(fù)雜性越來越高，導(dǎo)致噪聲的不確定性已無法采用傳統(tǒng)辨識方法中對噪聲的既有概率假設(shè)進(jìn)行描述，其本質(zhì)在于實(shí)際系統(tǒng)中的噪聲往往是不可測量的，無法直接用假定的概率分布函數(shù)來準(zhǔn)確表述，諸如最小二乘法[4]、極大似然法[5]等經(jīng)典算法應(yīng)用在實(shí)際系統(tǒng)中，往往會出現(xiàn)無法辨識或辨識效果不佳等現(xiàn)象.更廣泛的情況是，系統(tǒng)的噪聲滿足有界性，僅能從集合與空間圖形的角度分析噪聲對系統(tǒng)帶來的影響，噪聲的概率分布情況未知[6].

集員辨識方法假定參數(shù)模型和附加有界誤差，進(jìn)而計(jì)算參數(shù)的可行解集，確保獲得的測量集和所考慮的系統(tǒng)不確定性兼容.在研究系統(tǒng)建模問題時(shí)，系統(tǒng)參數(shù)的可行解集可以被視為一個(gè)多面體，其中多面體的頂點(diǎn)可以由模型結(jié)構(gòu)、測量值和有界噪聲準(zhǔn)確求解，多面體可由其頂點(diǎn)表述[7-9].同時(shí)，在多面體方法的基礎(chǔ)上也可以用多面體錐方法描述參數(shù)可行解集[10].因?yàn)樵诙嗝骟w中所有參數(shù)都可以使得系統(tǒng)正常運(yùn)行，因此這類方法也被稱為參數(shù)可行解集的精確描述法.然而，雖然精確描述方法可以準(zhǔn)確表述參數(shù)可行解集，但是其計(jì)算量過大，在計(jì)算過程中，多面體頂點(diǎn)未知且頂點(diǎn)個(gè)數(shù)會無規(guī)律增加或減少，在計(jì)算多面體頂點(diǎn)時(shí)，也會隨著參數(shù)維數(shù)的增加而變得相當(dāng)復(fù)雜.為了減少計(jì)算負(fù)擔(dān)，需要用近似參數(shù)可行解集包裹參數(shù)可行解集，其中近似參數(shù)可行解集的集合有橢球[11-13]、全對稱多胞體[14-16]和正多胞體[17-18]等，這類方法也被稱為近似描述方法.相比于精確描述，近似描述計(jì)算很簡單，卻也帶來參數(shù)可行解集寬松的缺點(diǎn).在近似描述中，橢球空間描述法因其與卡爾曼濾波的相似性且具有良好性能，近年來受到學(xué)者們的廣泛關(guān)注[19].然而，橢球空間描述法一直存在精度不高、保守性強(qiáng)的缺陷.相比而言，采用正多胞體進(jìn)行可行集描述，迭代過程計(jì)算簡單且圖形規(guī)則.雖然在集員辨識框架中，辨識時(shí)不變系統(tǒng)參數(shù)的方法和理論已日趨成熟，但跟蹤有界誤差設(shè)置中的時(shí)變參數(shù)問題仍存在很多不足，例如全對稱多胞體法中無法確保包裹k時(shí)刻相交區(qū)域的全對稱多胞體為最小全對稱多胞體[20]，正多胞體法中擴(kuò)展系數(shù)選擇參數(shù)變化最大值，無法達(dá)到擴(kuò)展系數(shù)選取最優(yōu)化[21]，文獻(xiàn)[22]中更新步驟繁瑣和文獻(xiàn)[23]中H矩陣維數(shù)增加使得計(jì)算量增加，所以辨識時(shí)變參數(shù)主要問題在于如何擴(kuò)展參數(shù)可行解集和降低計(jì)算量，確保在每一步的變化過程中參數(shù)可行解集包含變化后的全部可行參數(shù)解.與時(shí)不變參數(shù)集員辨識過程相似，但此類遞歸算法仍然沒有最優(yōu)擴(kuò)展參數(shù)可行解集.

本文提出一種基于正多胞體空間擴(kuò)展濾波的有界時(shí)變參數(shù)辨識方法.相比于可行解集對應(yīng)的其他空間結(jié)構(gòu)類型，正多胞體的形狀更加規(guī)則，且各邊都垂直或平行于坐標(biāo)軸，可以清晰展示在跟蹤參數(shù)變化時(shí)參數(shù)可行集的移動狀況.相比于文獻(xiàn)[24]，本文優(yōu)化了擴(kuò)展系數(shù)的選取，不是通過求取參數(shù)變化的最大值簡單確定擴(kuò)展系數(shù)，而是通過構(gòu)造擴(kuò)展系數(shù)方程，求解前k步擴(kuò)展系數(shù)最優(yōu)解得到，同時(shí)在辨識過程中，動態(tài)更新正多胞體約束條件.相比于文獻(xiàn)[25]，本文提出的方法在求解正多胞體時(shí)僅需要求解少量參數(shù)的線性規(guī)劃問題，在擴(kuò)展系數(shù)確定后，再確定每一步正多胞體約束條件，使得擴(kuò)展后的正多胞體包含變化后的參數(shù)，減少了計(jì)算多面體頂點(diǎn)的繁瑣步驟，從而降低了運(yùn)行過程的計(jì)算量.

2 問題描述

線性參數(shù)模型可以寫成以下回歸向量形式:

在本文中，滿足模型(1)的全部可行參數(shù)采用正多胞體空間結(jié)構(gòu)包裹，其約束條件可表示為

其中:A0∈Rn×nθ，b0∈Rn.式(2)定義了時(shí)變參數(shù)θ的變化上界取決于γ.按照γ的取值，系統(tǒng)(1)存在3種分類情況:

1) 時(shí)不變參數(shù)情況，γ=0;

2) 有界時(shí)變參數(shù)情況，0 ＜γ ＜∞;

3) 無界時(shí)變參數(shù)情況，γ=∞.

在本文中，主要研究有界時(shí)變參數(shù)系統(tǒng)中的參數(shù)辨識問題，即0 ＜γ ＜∞.通過在時(shí)不變參數(shù)系統(tǒng)約束條件下構(gòu)造擴(kuò)展系數(shù)方程，將參數(shù)可行解集在擴(kuò)展之后包含參數(shù)變化的最優(yōu)解γ作為最終擴(kuò)展系數(shù)，而不是簡單地選取θ變化的最大值作為擴(kuò)展系數(shù)，從而優(yōu)化了擴(kuò)展系數(shù)的選取過程，解決了當(dāng)前正多胞體辨識時(shí)變參數(shù)保守性高和多面體辨識時(shí)變參數(shù)計(jì)算復(fù)雜度高的缺陷.

3 正多胞體空間描述

定義1描述近似參數(shù)可行集的正多胞體O:

為正多胞體O的一對結(jié)構(gòu)面.

可以看出，正多胞體由約束條件Θ決定，通過求解2nθ個(gè)線性規(guī)劃方程可以得到包裹參數(shù)可行解集的最小正多胞體O*(Θ).

其中:Θ(k)為2nθ個(gè)約束條件可行解集的交集;ei表示nθ維單位矩陣的第i列，i=1，2，…，nθ.通過求解線性規(guī)劃方程可以得到第k步每個(gè)參數(shù)的上下界，從而得到最緊致正多胞體其中:

4 基于正多胞體空間擴(kuò)展濾波的時(shí)變參數(shù)辨識方法

若γ=0，即系統(tǒng)處于時(shí)不變狀態(tài)，有θ(k)=θ，此時(shí)系統(tǒng)的不確定性僅表現(xiàn)在未知但有界的噪聲.因此，對于每一時(shí)刻k，由測量值和回歸向量φ(k)得到包含參數(shù)空間的一個(gè)超平面帶Fk，

第k步正多胞體由前k-1時(shí)刻超平面帶和初始正多胞體相交得到

其中:Θk為第k步多面體為第k步包裹多面體的正多胞體.與求解時(shí)不變系統(tǒng)參數(shù)不同的是，在求解時(shí)變參數(shù)時(shí)，需要改變約束條件來擴(kuò)大正多胞體，否則正多胞體無法包含參數(shù)變化狀態(tài)下的可行解，導(dǎo)致辨識錯誤或者誤報(bào)為空集.本文通過擴(kuò)展系數(shù)擴(kuò)展約束條件，進(jìn)而放大擴(kuò)展正多胞體，最終達(dá)到擴(kuò)展參數(shù)可行解集的目的，擴(kuò)展后的正多胞體和多面體分別表示為記前k-1步擴(kuò)展，第k步不擴(kuò)展正多胞體和多面體分別為

定義2Fi|k(i≤k)是第k時(shí)刻在第i時(shí)刻的基礎(chǔ)上擴(kuò)展超平面帶，

其中Δi|k是由于從時(shí)刻i到k的參數(shù)變化引起的輸出估計(jì)的偏差界限，且有

在系統(tǒng)參數(shù)時(shí)變狀態(tài)下，考慮到參數(shù)變化上界γ，有

注1由式(10)和式(12)，可以得到Fk|k=Fk，即參數(shù)可行解集相同.

第k步非擴(kuò)展的正多胞體由前k-1時(shí)刻擴(kuò)展超平面帶和初始正多胞體相交得到

推論1假設(shè)k總是辨識的最后一步，已知

分別是在Rn上的多面體和正多胞體，那么C(k)⊕v(k)無意義，則k時(shí)刻參數(shù)可行集C(k)滿足

其中:

證假設(shè)有θ∈C(k-1)⊕v(k-1)，將存在α∈C(k-1)和‖βi‖1≤‖vi(k-1)‖1，其中βi，vi(k-1)分別表示β，v(k-1)的第i個(gè)元素，所以θ=α+β，從而得知

當(dāng)i=k時(shí)，θ只需滿足第k步超平面帶，無需擴(kuò)展;當(dāng)i ＜k時(shí)，|AAA|＜(k-i)‖φ(i)‖1，其中|AAA|表示矩陣AAA的行范數(shù)，φ(i)表示A的第i行向量，矩陣A和b隨步數(shù)增加自動更新，得到式(17).所以可以得到θ∈ˉC.而且因?yàn)棣痢蔆，那么一定存在α∈ooo*(C)，因此

證畢.

當(dāng)γ=0時(shí)，第k步正多胞體約束條件表示為

2018年6月15日，在河南省鄭州市鄭東新區(qū)，位于七里河南路與康平路交叉口的宏光意中大廈配樓同文酒店的門前人潮涌動。所有的來賓都滿懷期待地邁向二樓，在兩側(cè)擺滿祝賀花籃的大廳正中是一面墻，墻上呈現(xiàn)河南著名書法家王澄所題的“國學(xué)博覽館”五個(gè)行楷大字。這一天是宏光集團(tuán)公司歷時(shí)5年，由王澄嘔心瀝血、殫精竭慮擔(dān)綱策劃，開創(chuàng)國學(xué)博覽館開館揭幕的大喜之日。

依據(jù)推論1中式(17)和定義2中式(14)，當(dāng)i=k時(shí)，Δi|k=0，當(dāng)γ ＞0時(shí)，第k時(shí)刻非擴(kuò)展正多胞體約束條件:是從1到k-1時(shí)刻擴(kuò)展情況下Ak的行數(shù).

第k時(shí)刻擴(kuò)展正多胞體約束條件為

其中:

ak;i是矩陣Ak的第i行，Δbk-1;i是Δbk-1的第i個(gè)值.此時(shí)，式(24)可以看作一個(gè)關(guān)于γ的線性規(guī)劃方程，其余參數(shù)可由系統(tǒng)測量值和初始定義正多胞體獲得，因此γ滿足線性規(guī)劃方程:

其中:f=(0，…，0，1)，x=(θTγ)T，長度為nθ+1.當(dāng)前k時(shí)刻所有γk得到后，最終的擴(kuò)展系數(shù)

圖1描述了相同情況下擴(kuò)展正多胞體和非擴(kuò)展正多胞體的變化趨勢，虛線包裹區(qū)域代表擴(kuò)展正多胞體，實(shí)線包裹區(qū)域代表無擴(kuò)展正多胞體，黑線包裹區(qū)域代表初始正多胞體，因?yàn)樵跀U(kuò)展多胞體情況下，約束超平面帶變寬，可以看出擴(kuò)展正多胞體尺寸更大，以便包含時(shí)變參數(shù).

圖1 正多胞體空間擴(kuò)展變化情況Fig.1 The space variation of the extended orthotope

基于正多胞體空間擴(kuò)展濾波的時(shí)變參數(shù)辨識方法算法步驟:

步驟1定義辨識步數(shù)為L，定義初始正多胞體，得到初始約束條件A0θ≤b0，γ0=0;

步驟2根據(jù)輸入數(shù)據(jù)u(k)，測量數(shù)據(jù)y(k)和回歸向量φ(k)構(gòu)造第k步非擴(kuò)展正多胞體約束條件，

步驟4更新γ值，γ=max(γ，γk).當(dāng)k=L時(shí)得到最大擴(kuò)展系數(shù)γ，同時(shí)重置k=1，否則設(shè)置k=k+1，返回第2步;

步驟5由式(25)更新擴(kuò)展正多胞體約束條件

步驟6分別取i=1，2，…，nθ，求解正多胞體頂點(diǎn)，

步驟7構(gòu)造擴(kuò)展正多胞體

其中:

步驟8置k=k+1，返回第5步;當(dāng)k=L時(shí)，算法結(jié)束.

5 仿真分析

例1考慮下列線性參數(shù)模型:

其中:θ1，θ2為待辨識的時(shí)變參數(shù)，其真實(shí)的變化情況分別為

e(t)為未知但有界噪聲，取輸入信號u(t)∈U[-5，5]，系統(tǒng)噪聲e(t)∈U[-0.3，0.3].為比較跟蹤參數(shù)變化的情況，定義歐氏距離Δθ=其中θ為參數(shù)變化時(shí)的真實(shí)值，?θ為時(shí)變參數(shù)估計(jì)值[26].

采用本文方法辨識參數(shù)變化值時(shí)效果如圖2-3所示.

圖2 參數(shù)θ1辨識曲線Fig.2 Identification curve of parameter estimate θ1

圖3 參數(shù)θ2辨識曲線Fig.3 Identification curve of parameter estimate θ2

從圖2-3可以看出，基于正多胞體空間擴(kuò)展濾波的有界時(shí)變參數(shù)算法能夠有效辨識時(shí)變參數(shù)，參數(shù)真值一直位于正多胞體上下界范圍內(nèi)，隨著參數(shù)變化，正多胞體估計(jì)值跟蹤變化，具有很好的效果.

接下來設(shè)置參數(shù)θ1不變，參數(shù)θ2=-0.5×sin(π×k×0.01)，與文獻(xiàn)[24]中選擇參數(shù)變化最大值作為擴(kuò)展系數(shù)方法和文獻(xiàn)[25]中多面體方法比較的仿真結(jié)果如圖4-5所示.

從圖4-5均可看出，正多胞體方法和多面體方法可行參數(shù)解集上下界相同，最大值正多胞體上下界要明顯大于其余兩種方法，且正多胞體的中心值相較多面體的中心值和參數(shù)變化最大值方法更加貼近真實(shí)參數(shù)點(diǎn).與兩種方法比較的對比分析仿真結(jié)果分別如圖6-7所示.

圖4 正多胞體和多面體辨識時(shí)變參數(shù)θ1比較情況Fig.4 Comparison of time-varying parameter estimate θ1 between orthotopic and polyhedral based algorithms

圖5 正多胞體和多面體辨識時(shí)變參數(shù)θ2比較情況Fig.5 Comparison of time-varying parameter estimate θ2 between orthotopic and polyhedral based algorithms

圖6中，本文在k=50～140范圍內(nèi)間隔15步取正多胞體和多面體數(shù)據(jù)，記號“x”表示當(dāng)前時(shí)刻參數(shù)真值，記號“o”表示當(dāng)前時(shí)刻正多胞體中心估計(jì)值，記號“*”表示當(dāng)前時(shí)刻多面體中心估計(jì)值，記號□表示最大值正多胞體中心估計(jì)值.圖6中，不同顏色下標(biāo)準(zhǔn)矩形包裹區(qū)域?qū)嵕€表示正多胞體，點(diǎn)劃線表示最大值正多胞體，矩形內(nèi)圖形包裹區(qū)域表示多面體.可以看出正多胞體中每個(gè)參數(shù)的最大和最小值是根據(jù)準(zhǔn)確描述可行參數(shù)解集的邊界上下限確定，所以基于正多胞體空間擴(kuò)展濾波方法和基于多面體的空間擴(kuò)展濾波方法在辨識時(shí)變參數(shù)時(shí)上下界一致.就辨識精度而言，圖6中的正多胞體中心估計(jì)值更接近參數(shù)真值，反映出基于正多胞體空間擴(kuò)展濾波方法針對時(shí)變參數(shù)變化可以做到更有效的跟蹤.從圖7中可以看出，在以歐式距離為依據(jù)判斷算法優(yōu)劣的過程中，考慮到時(shí)變參數(shù)系統(tǒng)中參數(shù)每一時(shí)刻都在變化，估計(jì)值也會隨著參數(shù)變化而時(shí)刻跳動，所以在誤差分析曲線中，誤差曲線并不是逐漸減小最終趨于不變，但是在誤差分析過程中，基于正多胞體空間擴(kuò)展濾波方法的估計(jì)誤差要比基于多面體空間擴(kuò)展濾波方法更小.

圖6 正多胞體和多面體辨識過程中遞歸演化Fig.6 Recursive evolution by the orthotopic and polyhedral algorithms

圖7 正多胞體和多面體誤差曲線對比Fig.7 Comparison of error curves by the orthotopic and polyhedral algorithms

例2為了進(jìn)一步驗(yàn)證基于正多胞體空間擴(kuò)展濾波方法解決有界時(shí)變參數(shù)系統(tǒng)辨識問題的有效性，下面采用風(fēng)力發(fā)電機(jī)漿距子系統(tǒng)作為仿真示例進(jìn)行分析.

在風(fēng)力發(fā)電機(jī)系統(tǒng)中，漿距子系統(tǒng)是控制槳葉和漿距角的重要組成部分，如圖8所示.其結(jié)構(gòu)模型為[27]

其中:β和βa分別為漿距角和角速度大小，βr為漿距參考值，ζ和ωn分別為阻尼系數(shù)和系統(tǒng)自然頻率，ζ和ωn會隨著液壓的變化而變化，主要是由于主線壓力下降，液壓系統(tǒng)中所考慮的故障可能導(dǎo)致參數(shù)動態(tài)變化.這種動態(tài)變化是標(biāo)稱值0.6 rad/s和0.9 rad/s之間的阻尼系數(shù)的變化以及3.42 rad/s和標(biāo)稱值11.11 rad/s之間的系統(tǒng)自然頻率的變化[28].將槳距子系統(tǒng)數(shù)學(xué)模型(33)近似為二階系統(tǒng)[29]

其中:y=β，u=βr.為了滿足參數(shù)在規(guī)定區(qū)間內(nèi)動態(tài)變化，將阻尼系數(shù)和自然頻率的變化分別設(shè)置為

據(jù)文獻(xiàn)[27]取采樣時(shí)間為Ts=0.01 s，對該閉環(huán)系統(tǒng)離散化為A(z)y(t)=B(z)u(t)+e(t).根據(jù)離散化條件，A(z)，B(z)均為多項(xiàng)式，其中:

為了驗(yàn)證所提算法跟蹤參數(shù)變化的有效性，在參數(shù)變化過程中，分析B(z)多項(xiàng)式對應(yīng)的待辨識參數(shù)并描述參數(shù)辨識效果和正多胞體演化情況，并取輸入u(t)∈U[-25，25]，噪聲e(t)∈U[-0.01，0.01].同時(shí)與文獻(xiàn)[24]中選擇參數(shù)變化最大值作為擴(kuò)展系數(shù)進(jìn)行對比，比較結(jié)果如圖9-10所示.

圖8 風(fēng)力發(fā)電機(jī)系統(tǒng)Fig.8 Wind turbine system

圖9 擴(kuò)展系數(shù)不同時(shí)參數(shù)估計(jì)θn1比較Fig.9 Comparison of parameter estimates θn1 when the expansion coefficients are different

圖10 擴(kuò)展系數(shù)不同時(shí)參數(shù)估計(jì)θn2比較Fig.10 Comparison of parameter estimates θn2 when the expansion coefficients are different

從圖9-10中可以看出，采用本文方法和文獻(xiàn)[24]提出的方法分別求解擴(kuò)展系數(shù)，都能有效地跟蹤時(shí)變參數(shù)變化情況.相比而言，本文提出的基于正多胞體空間擴(kuò)展的濾波方法辨識過程更平滑，受噪聲擾動的影響更小，優(yōu)于文獻(xiàn)[24]提出的選擇參數(shù)變化最大值作為擴(kuò)展系數(shù)方法.

圖11 擴(kuò)展系數(shù)不同下的正多胞體遞歸演化Fig.11 Recursive evolution of orthotopes with different expansion coefficients

此外，圖11中在k=50～150之間間隔25步取不同擴(kuò)展系數(shù)正多胞體數(shù)據(jù)，記號“x”表示當(dāng)前時(shí)刻參數(shù)真值，實(shí)線包裹區(qū)域和記號“o”分別表示當(dāng)前時(shí)刻線性規(guī)劃擴(kuò)展系數(shù)正多胞體和正多胞體中心估計(jì)值，虛線包裹區(qū)域和記號“*”表示當(dāng)前時(shí)刻選擇參數(shù)變化最大值擴(kuò)展系數(shù)正多胞體和正多胞體中心估計(jì)值.可以看出基于線性規(guī)劃的正多胞體要更緊致，因此算法保守性更低，中心估計(jì)值也更靠近參數(shù)真值.

6 結(jié)論

本文基于正多胞體濾波過程，結(jié)合線性規(guī)劃求解參數(shù)擴(kuò)展系數(shù)，提出了一種時(shí)變線性系統(tǒng)參數(shù)辨識算法，在給定的噪聲和參數(shù)變化上下界的條件下，僅依據(jù)測量值和系統(tǒng)觀測值求得擴(kuò)展系數(shù)最優(yōu)解.本文所提出的正多胞體空間擴(kuò)展濾波方法在辨識時(shí)變參數(shù)過程中，利用正多胞體的空間擴(kuò)展和移動描述了辨識過程，同時(shí)在計(jì)算過程中，通過求解有限個(gè)線性規(guī)劃問題得到緊致的正多胞體，計(jì)算步驟簡單.此外，本文通過數(shù)值仿真和風(fēng)力發(fā)電機(jī)的漿距子系統(tǒng)為例分別進(jìn)行仿真，可以直觀看出該方法辨識有界時(shí)變參數(shù)的可行性和有效性.