王永坤，孫明瑋，劉忠信，陳增強(qiáng)

(南開大學(xué) 自動化系，300071天津)

存在撓性結(jié)構(gòu)或包含撓性部件的一類系統(tǒng)稱為撓性系統(tǒng).撓性系統(tǒng)中的撓性部件分為很多種，例如衛(wèi)星的太陽能帆板［1-3］、機(jī)器人的機(jī)械臂［4］、大型運(yùn)載火箭或具有一定范圍長細(xì)比的有翼飛行器等.此類系統(tǒng)的一個主要特點(diǎn)是系統(tǒng)中存在弱阻尼的諧振模態(tài)，在控制系統(tǒng)設(shè)計(jì)時需要特殊考慮，以保證閉環(huán)系統(tǒng)穩(wěn)定［5］.撓性系統(tǒng)振動模態(tài)諧振頻率的理想值或者是通過理論方法計(jì)算出來的，或者是通過地面振動試驗(yàn)測量出來的.理論計(jì)算中分布質(zhì)量和結(jié)構(gòu)特性存在的誤差，地面試驗(yàn)中空天的不一致性以及測量誤差等，都會導(dǎo)致諧振頻率誤差.飛行器在飛行過程中，振動模態(tài)的一些特征參數(shù)會隨著溫度和質(zhì)量的變化而變化，因此通過理論計(jì)算或者模態(tài)振動試驗(yàn)得到的諧振頻率具有一定的不確定性，在控制系統(tǒng)的設(shè)計(jì)過程中必須考慮其對于閉環(huán)魯棒穩(wěn)定性的影響.文獻(xiàn)［6-7］指出采用H∞回路成形法［8-10］設(shè)計(jì)撓性控制系統(tǒng)時，由于系統(tǒng)中弱阻尼模態(tài)的存在，會減小諧振頻率的實(shí)際攝動范圍，并不能保證系統(tǒng)具有期望的魯棒性，但是并沒有給出計(jì)算諧振頻率實(shí)際攝動范圍的方法，同時文獻(xiàn)［7］計(jì)算出的系統(tǒng)H∞范數(shù)的逆并不代表系統(tǒng)中諧振頻率攝動值的百分比，它實(shí)際是系統(tǒng)綜合攝動范圍的一個指標(biāo).由于諧振頻率分別出現(xiàn)在一次項(xiàng)和二次項(xiàng)中，導(dǎo)致其非線性的耦合攝動形式，傳統(tǒng)的線性或雙線性分析方法不再適用［11-12］.目前還沒有發(fā)現(xiàn)文獻(xiàn)報道存在系統(tǒng)性的計(jì)算方法可以直接確定撓性系統(tǒng)諧振頻率的攝動范圍.針對這個問題，工程上只能采用攝動試驗(yàn)法，通過試湊的方式摸索得到諧振頻率的穩(wěn)定性邊界，該方法不僅費(fèi)時費(fèi)力，而且精度難以得到保證.因此，本文給出一種基于D分割法［13］的系統(tǒng)魯棒穩(wěn)定性分析方法，給出弱阻尼模態(tài)下諧振頻率攝動范圍的幾何求法，并通過一個算例驗(yàn)證了該方法的有效性，同時通過另一算例說明了本方法可直接推廣應(yīng)用到具有多個模態(tài)的撓性系統(tǒng)分析中.

1 問題描述

本文以具有太陽能帆板的衛(wèi)星為例來分析撓性系統(tǒng)諧振頻率的攝動范圍，但文中方法也可直接推廣到有翼飛行器等撓性對象的分析中.一般來說，撓性結(jié)構(gòu)振動具有無窮多個模態(tài)，但是振動能量主要集中在低階項(xiàng)上，因此撓性系統(tǒng)設(shè)計(jì)時常以一階模態(tài)為主，將高階模態(tài)作為未建模動態(tài)來處理.在控制系統(tǒng)設(shè)計(jì)時要求控制器的增益在一階模態(tài)頻率后迅速衰減.圖1給出了衛(wèi)星示意圖.在圖1(a)中，θ1是衛(wèi)星主體與星體的夾角，θ2是衛(wèi)星姿態(tài)角，代表星體傳感器與儀表設(shè)備的夾角，圖1(b)所示為將傳感器安裝到與θ2相連的圓盤上的衛(wèi)星等效機(jī)械系統(tǒng)圖.對圖1所示模型，假設(shè)用彈簧鏈接的兩個質(zhì)體的扭矩常數(shù)為k，粘性阻尼系數(shù)為b，建立線性模型如下［4］:

其中Tc為作用在主體上的控制力矩，J1、J2分別是圖1(b)中兩個質(zhì)體的轉(zhuǎn)動慣量．由式(1)可得含有一維振動模態(tài)的撓性系統(tǒng)模型為

其中

式中c0、c1是模態(tài)幅值，只考慮一階模態(tài)時c0=-c1，ωn和ξ分別是系統(tǒng)的諧振頻率和阻尼比．在衛(wèi)星飛行過程中，參數(shù)k會隨溫度的波動以及外形結(jié)構(gòu)特性而變化，因此由式(3)可知諧振頻率ωn會發(fā)生攝動．從式(2)可以看出ωn分別出現(xiàn)在一次項(xiàng)和二次項(xiàng)中，導(dǎo)致其攝動形式是非線性的，其對于閉環(huán)控制系統(tǒng)的影響不是簡單的線性或者雙線性形式，傳統(tǒng)的分析方法無法直接應(yīng)用［11-12］．

圖1 衛(wèi)星及其雙體

2 基于D分割法的魯棒穩(wěn)定性分析

2.1 D分割法

考慮如下一類反饋控制系統(tǒng)的閉環(huán)特征多項(xiàng)式:

其中ak(z)和bk(z)分別為參數(shù)向量z的連續(xù)函數(shù)，z=［z1，z2，…，zl］Τ．L為系統(tǒng)的時滯常數(shù)，L∈Η?R，zi∈Q?RL，R為實(shí)數(shù)集．

對于給定l階向量z和L，若特征多項(xiàng)式δ(s，L，z)無非負(fù)實(shí)部零點(diǎn)，則稱其漸近穩(wěn)定．由于δ(s，L，z)的零點(diǎn)為z和L的連續(xù)函數(shù)，空間Η×Q被超曲面分割成若干個區(qū)域，超曲面上的點(diǎn)對應(yīng)的δ(s，L，z)至少有一個純虛根或s=∞，這種分解法稱為Η×Q的D分割法．D分割法的實(shí)質(zhì)是將s平面上的虛軸映射為參數(shù)空間Η×Q的超曲面，根據(jù)實(shí)數(shù)根穿越原點(diǎn)和穿越無窮，復(fù)根穿越虛軸得出D分割邊界［14］．

其中

顯然，隨著z和L的連續(xù)變化，穩(wěn)定系統(tǒng)變化到不穩(wěn)定系統(tǒng)中間必然穿越?D．換句話說，研究魯棒穩(wěn)定邊界只需要通過分析?D就可以實(shí)現(xiàn)．

2.2 一維線性不確定系統(tǒng)的魯棒穩(wěn)定性分析

由于本文所研究問題的時滯特性可以忽略，因此考慮如下特征多項(xiàng)式:

其中ai(i=0，1，…，n)是多項(xiàng)式系數(shù)，假設(shè)ai(i=0，1，…，n)均可由特定參數(shù)d的線性函數(shù)表示，d是式(5)中的不確定參數(shù)，則上式可寫為

其中p(s)和q(s)是兩個互質(zhì)的多項(xiàng)式．假設(shè)當(dāng)式(6)中的d=0時，Pc(s)是穩(wěn)定的．

此時，基于D分割法，魯棒穩(wěn)定邊界?D上的p(s)和q(s)分別是

其中f(ω)、e(ω)、h(ω)、g(ω)均為實(shí)數(shù)多項(xiàng)式．由Pc(jω)=0有

消去上式中d，可得

求解式(8)得到的正實(shí)根ωi，表示系統(tǒng)的截止頻率．將ωi代入式(7)可得

定義

其中d+，d-即為在滿足系統(tǒng)魯棒穩(wěn)定前提下，不確定參數(shù)d的正負(fù)攝動邊界．

2.3 二維線性不確定系統(tǒng)的魯棒穩(wěn)定性分析

當(dāng)式(5)的某些系數(shù)中存在d1，d2兩個不確定參數(shù)，并且這些系數(shù)均可表示為d1，d2的雙線性函數(shù)形式時，就成為一個二維線性不確定問題．假設(shè)當(dāng)d1=d2=0時，Pc(s)是穩(wěn)定的．為了確定d1、d2的魯棒穩(wěn)定區(qū)域，首先，令d2=0，采用2.2節(jié)的方法求得此時參數(shù)d1的魯棒穩(wěn)定區(qū)域?yàn)槿缓?，在區(qū)域內(nèi)均勻采樣N個點(diǎn)，對于其中每一個d1，i(i=1，2，…，N)，求解此時關(guān)于不確定參數(shù)d2的魯棒穩(wěn)定區(qū)域;最后，分別連接點(diǎn)和和其中i=1，2，…，N-1，不確定參數(shù)的魯棒穩(wěn)定區(qū)域如圖2陰影部分所示．對于多維線性不確定系統(tǒng)，其魯棒性分析方法類似．對于存在二維或者多維不確定參數(shù)的系統(tǒng)，當(dāng)各個不確定參數(shù)之間存在一定的關(guān)系，而參數(shù)之間的關(guān)系曲線與上述得到的魯棒穩(wěn)定區(qū)域的交點(diǎn)即為各個參數(shù)的攝動邊界．而攝動點(diǎn)一般在攝動邊界內(nèi)選取是有意義的．

圖2 不確定參數(shù)魯棒穩(wěn)定區(qū)域示意圖

2.4 撓性系統(tǒng)的魯棒穩(wěn)定性分析

考慮式(2)所表示的撓性系統(tǒng)，系統(tǒng)的閉環(huán)特征多項(xiàng)式為

其中C(s)是控制器，諧振頻率ωn是攝動因子，分別出現(xiàn)在一次項(xiàng)和二次項(xiàng)中．

現(xiàn)在考慮一個更廣義的二維不確定性系統(tǒng)．將式(9)中ωn的一次項(xiàng)和二次項(xiàng)分別表示為

其中ωn和ωn2分別表示諧振頻率攝動后的一次項(xiàng)和二次項(xiàng)，n是諧振頻率的標(biāo)稱值，d1和d2分別表示攝動的百分比，此時存在d1和d2兩個不確定參數(shù)，轉(zhuǎn)化為二維不確定性問題．采用2.3節(jié)的方法，假設(shè)在給定n、ξ、c0、c1和C(s)的情況下，繪制出不確定參數(shù)魯棒穩(wěn)定區(qū)域即圖3中的陰影部分．

圖3 魯棒穩(wěn)定參數(shù)區(qū)域

對于包含二階振動模態(tài)的撓性系統(tǒng)，即存在ωn1和ωn2兩個不確定參數(shù)．首先令ωn2固定為標(biāo)稱值，通過上述方法求得ωn1的攝動區(qū)間．然后在ωn1的攝動區(qū)間內(nèi)均勻采樣，對于每一個采樣值，可以求出ωn2的攝動邊界值，最終得到［ωn1，ωn2］的整體攝動區(qū)域．對于包含多階振動模態(tài)的撓性系統(tǒng)，分析方法類似．

3 算 例

例1 考慮具有兩個太陽能帆板的衛(wèi)星模型［7］，其中c0=1.731 9×10-5，c1=3.785 9×10-4，ωn1=1.539 rad/s，ξ1=0.003，此類系統(tǒng)具有明顯的撓性，一階模態(tài)的c1大于剛性模態(tài)的c0．

這里直接引用文獻(xiàn)［7］中的控制器設(shè)計(jì)結(jié)果，通過回路成形法得到H∞控制器．首先給被控對象加入串聯(lián)補(bǔ)償環(huán)節(jié)，使其具有期望的開環(huán)特性，即回路成形．然后針對補(bǔ)償后的廣義對象，設(shè)計(jì)H∞控制器，保證閉環(huán)系統(tǒng)的穩(wěn)定性．其中補(bǔ)償函數(shù)選擇為

由式(2)、(10)和(11)可得廣義的系統(tǒng)傳遞函數(shù)為

當(dāng)d1=d2=0時，針對標(biāo)稱系統(tǒng)設(shè)計(jì)的H∞控制器為

因此可得閉環(huán)特征多項(xiàng)式

令d2=0，根據(jù)2.2節(jié)一維不確定參數(shù)分析方法可以計(jì)算出，d1的允許攝動區(qū)間為［-0.611 4，0.312 3］．由于直線d1=d2只通過第一、三象限，因此針對每一個采樣的d1，i;當(dāng)d1，i≥0時，只需求出此時對應(yīng)的d2攝動的上邊界值當(dāng)d1，i＜0時，只需求出此時對應(yīng)的d2攝動的下邊界值，結(jié)果如圖4所示．直線d1=d2與一、三象限邊界的交點(diǎn)分別為(0.138 6，0.138 6)和(-0.351 4，-0.351 4)，表明諧振頻率ωn1增加13.86%或減小35.14%后，系統(tǒng)處于臨界穩(wěn)定狀態(tài)，ωn1的攝動范圍即為［-35.14%，13.86%］．

圖4 諧振頻率ωn1的攝動范圍

攝動系統(tǒng)的階躍響應(yīng)曲線如圖5所示，其中實(shí)線表示標(biāo)稱系統(tǒng)的響應(yīng)曲線，短虛線和長虛線分別代表諧振頻率攝動-35.14%和13.86%時的響應(yīng)曲線．可以看出當(dāng)ωn1的攝動達(dá)到所求的正負(fù)邊界時，系統(tǒng)出現(xiàn)等幅振蕩，驗(yàn)證了本文提出方法的正確性．

圖5 系統(tǒng)的階躍響應(yīng)曲線

而采用文獻(xiàn)［7］的分析方法，求得系統(tǒng)的H∞范數(shù)γ=2.486，則γ-1≈0.402 3，根據(jù)小增益定理可得只要對象的攝動小于40.23%，就可以保證閉環(huán)魯棒穩(wěn)定．而實(shí)際情況是當(dāng)諧振頻率正攝動14%時就達(dá)到了穩(wěn)定邊界，文獻(xiàn)［7］得出的結(jié)論是:此情況的出現(xiàn)是由于對象存在弱阻尼模態(tài)，壓縮了實(shí)際允許的攝動范圍，因此導(dǎo)致H∞回路成形法設(shè)計(jì)的系統(tǒng)并不一定具有魯棒性．

本文經(jīng)過分析，γ-1≈0.402 3表示當(dāng)系統(tǒng)中的不確定參數(shù)ωn1攝動值達(dá)到最大時，對應(yīng)整體不確定項(xiàng)的H∞范數(shù)為0.402 3．定義

其中G和G0分別表示實(shí)際對象和標(biāo)稱對象，ΔG表示對象的整體不確定項(xiàng)，其H∞范數(shù)與Δωn1的關(guān)系如圖6所示，其中橫坐標(biāo)表示諧振頻率的攝動值，縱坐標(biāo)表示‖ΔG‖∞．從圖中可以看出當(dāng)諧振頻率ωn1攝動達(dá)到正向邊界即13.86%時，‖ΔG‖∞為0.407 5，幾乎等于文獻(xiàn)［7］給出的0.402 3．因此，系統(tǒng)H∞范數(shù)的逆并不表示某一參數(shù)的攝動范圍，而是一種結(jié)構(gòu)參數(shù)攝動引起的非結(jié)構(gòu)不確定性變化范圍，二者本質(zhì)上是不同的．

圖6 整體不確定項(xiàng)的H∞范數(shù)

例2 仍采用算例1的模型以及控制器設(shè)計(jì)結(jié)果，并在此基礎(chǔ)上加入一個二階振動模態(tài)，得到如下包含兩個振動模態(tài)的撓性系統(tǒng):

其中c0、c1、ωn1、ξ1保持例1中的數(shù)值不變，c2=1.02×10-4，ωn2=5.2，ξ2=0.002，諧振頻率ωn1和ωn2是系統(tǒng)中的不確定參數(shù)．

首先令ωn2固定為標(biāo)稱值，對ωn1進(jìn)行攝動分析．通過本文提出的方法可得ωn1的攝動比例為d1∈［-0.36，0.163］．對比例1可見，由于二階模態(tài)的能量比較低，對于一階模態(tài)振動頻率的攝動范圍影響不是很大，因此一般在對撓性系統(tǒng)進(jìn)行分析時，通常只需考慮能量較大的低階振動模態(tài)即可．當(dāng)ωn1的攝動達(dá)到臨界值時，系統(tǒng)的階躍響應(yīng)曲線如圖7所示出現(xiàn)等幅振蕩．

在區(qū)間［-0.36，0.163］內(nèi)均勻采樣，對于每一個對應(yīng)的ωn1值，求得ωn2攝動比例的邊界值，最終可得ωn1、ωn2的攝動魯棒穩(wěn)定區(qū)域，即圖8中曲線與坐標(biāo)軸所圍區(qū)域，曲線表示攝動的邊界．需要說明的是二階模態(tài)諧振頻率ωn2的攝動沒有上界，反映了諧振頻率越高對于整個剛體的影響越小，符合實(shí)際經(jīng)驗(yàn)．在求得區(qū)域選取內(nèi)點(diǎn)(-0.2，-0.5)，階躍響應(yīng)如圖9(a)所示，可以看出此時系統(tǒng)處于收斂的趨勢;在邊界上選取(-0.32，-0.713)和(0.07，-0.793)，階躍響應(yīng)如圖9(b)所示，輸出產(chǎn)生等幅振蕩.

圖7 ωn 攝動時系統(tǒng)的階躍響應(yīng)曲線

圖8 諧振頻率ωn1、ωn2的攝動魯棒穩(wěn)定范圍

圖9 ωn1、ωn2攝動時系統(tǒng)的階躍響應(yīng)曲線

4 結(jié) 語

本文針對撓性系統(tǒng)中不確定諧振頻率的魯棒穩(wěn)定性分析問題，設(shè)計(jì)了一種基于D分割法的幾何分析方法，給出了振動模態(tài)諧振頻率的精確攝動范圍.通過數(shù)值算例驗(yàn)證了該方法的有效性，對于實(shí)際工程具有一定參考價值.該方法由于突破了傳統(tǒng)魯棒穩(wěn)定性分析要求攝動參數(shù)必須是線性或者雙線性的形式限制，也可推廣到其他含參數(shù)非線性攝動形式的場合，通過繪制相應(yīng)的幾何曲線求取特定參數(shù)的魯棒穩(wěn)定性范圍.

［1］MCFARLANE D C，GLOVER K.Robust controller design using normalized coprime factor plant descriptions［M］.New York:Springer Verlag，1990:205-208.

［2］劉暾，楊大明.撓性衛(wèi)星動力學(xué)及姿態(tài)控制模型的建立［J］.哈爾濱工業(yè)大學(xué)學(xué)報，1985，17(A8):1-17.

［3］崔乃剛，劉家夫，榮思遠(yuǎn).柔性太陽帆航天器動力學(xué)建模與姿態(tài)控制［J］.哈爾濱工業(yè)大學(xué)學(xué)報，2011，43(7):1-5.

［4］FRANKLIN G F，POWELL J D，ABBAS E N.Feedback control of dynamic systems［M］.4th Ed.Reading，MA:Addison-Wesley，1986:488-491.

［5］游江，趙國良，羅耀華.H∞混合靈敏度的魯棒性設(shè)計(jì)及在撓性系統(tǒng)控制中的應(yīng)用［J］.儀器儀表學(xué)報，2006，27(6):556-558.

［6］何朕，孟范偉，劉偉，等.H∞回路成形設(shè)計(jì)的魯棒性［J］.自動化學(xué)報，2010，36(6):890-893.

［7］何朕，孟范偉，王廣雄，等.撓性系統(tǒng)的控制設(shè)計(jì)和H∞回路成形法［J］.控制理論與應(yīng)用，2013，30(8):1014-1020.

［8］MCFARLANE D C，GLOVER K.A loop shaping design procedure usingH∞synthesis［J］.IEEE Transactions on Automatic Control，1992，37(6):759-769.

［9］KAITWANIDVILAI S，OLRANTHICHACHAT P.Robust loop shaping—fuzzy gain scheduling control of a servopneumatic system using particle swarm optimization approach［J］.Mechatronics，2011，21(1):11-21.

［10］OSINUGA M，PATRA S，LANZON A.Smooth weight optimization inH∞loop-shaping design［J］.Systems＆Control Letters，2010，59(11):663-670.

［11］WEDELL E，CHUANG C H，WIE B.Stability robustness margin computation for structured realparameter perturbations［J］.Journal of Guidance，Control，and Dynamics，1991，14(3):607-614.

［12］WIE B，LU J B.Two real critical constraints for real parametermargin computation ［J］. Journal of Guidance，Control，and Guidance，1994，17(3):561-569.

［13］GRYAZINA E N，POLYAK B T.Stability regions in the parameterspace:D-decomposition revisited［J］.Automatica，2006，42(1):13-26.

［14］GRYZINA E N.The D-decomposition theory［J］.Automation and Remote Control，2004，65(12):1872-1884.