王夢菲,張軍,*
1. 北京控制工程研究所, 北京 100190 2. 空間智能控制技術國家級重點實驗室,北京 100190
由于未建模動態(tài)和建模的不確定性,標稱模型下設計的控制器在實際中性能會出現(xiàn)偏差。在軌時航天器常受到內外干擾作用,特別是帶有大型撓性子結構的現(xiàn)代航天器,其撓性結構易受干擾激發(fā)產(chǎn)生振動,嚴重影響姿態(tài)指向精度和穩(wěn)定度,甚至導致航天器失穩(wěn),如美國的Explore I[1]。因此,大型撓性航天器高精度、高穩(wěn)定度姿態(tài)控制是一個熱點問題也是控制界的難題。
文獻[2-4]研究了經(jīng)典PID控制方法在撓性航天器姿態(tài)控制中的應用;當模型存在不確定性、未建模模態(tài)發(fā)生振動以及耦合振動發(fā)生變化時,其控制性能下降,且PID控制對外部干擾的魯棒性也有待加強。文獻[5]針對結構振動,基于狀態(tài)反饋設計了時間最優(yōu)的LQR控制器,但不適合高頻范圍內未建模模態(tài)的振動抑制,且不能滿足模型不確定性時的魯棒控制。文獻[6]改進了自適應滑??刂品椒ǎ畲笙薅鹊販p少了不確定性、擾動和柔性動態(tài)坐標測量等帶來的困難,并在一定程度上改善了滑??刂乒逃械亩墩瘳F(xiàn)象。文獻[7]證明了基于姿態(tài)角和角速度的反饋對于大范圍攝動的穩(wěn)定性,并針對含參數(shù)不確定性和未建模動態(tài)的姿態(tài)控制系統(tǒng),設計了基于線性矩陣不等式(linear matrix inequality, LMI)的魯棒H∞控制器。文獻[8]針對撓性航天器姿態(tài)機動問題,以姿態(tài)達到目標角度且保持穩(wěn)定的時間為適應度值,將粒子群優(yōu)化(particle swarm optimization, PSO)算法與黃金分割控制等結合,設計了一種離線控制器;相對于常規(guī)的控制方法,該方法的優(yōu)點在于:可根據(jù)系統(tǒng)實際動力學特性和控制要求選擇最優(yōu)機動路徑,實現(xiàn)姿態(tài)快速機動和穩(wěn)定。文獻[9-11]將神經(jīng)網(wǎng)絡、模糊控制理論應用于姿態(tài)控制系統(tǒng)設計;但模糊控制效果往往依賴于操作者經(jīng)驗;為提高可靠性,神經(jīng)網(wǎng)絡一般需在線學習,受星載計算機限制。
上述撓性航天器姿態(tài)控制方法大多是針對單一控制目標而設計的,但在工程實際中,控制系統(tǒng)的設計需考慮魯棒性、快速性、精度和執(zhí)行機構輸出飽和等多目標要求,這就屬于典型的多目標設計問題,關于此類問題的研究較少。
針對上述多目標魯棒控制問題,本文在自適應PSO算法中加入變異操作,提出一種差分粒子群優(yōu)化(differential particle swarm optimization, DPSO)算法,以提高粒子的多樣性和搜索性能;同時將該算法與魯棒D-穩(wěn)定極點約束以及Pareto最優(yōu)結合,涉及的數(shù)學公式較少,設計靈活。相比于傳統(tǒng)帶極點配置的LMI方法,本文方法不受特定指標和約束模式的限制,可避免凸優(yōu)化過程,并減小求解的保守性[12]。相對于文獻[13-14]提到的借助特征結構配置的參數(shù)化方法,本文優(yōu)化過程無需將多目標轉化為一個加權目標函數(shù),可避免加權系數(shù)的選擇困難,特別是多目標存在矛盾關系的情況;而且文獻中方法需將指標轉化為特定的參數(shù)化形式,對參數(shù)攝動敏感度的優(yōu)化無法保證系統(tǒng)魯棒穩(wěn)定。本文方法適用的指標形式更廣,能實現(xiàn)所允許的參數(shù)不確定范圍內的魯棒D-穩(wěn)定;可為大型撓性航天器魯棒控制器設計提供一種思路。
帶有一對太陽翼的航天器姿態(tài)動力學方程可表示為[15]:
(1)
本文僅考慮航天器姿態(tài)定向保持和姿態(tài)鎮(zhèn)定情況,假設姿態(tài)為小角度,此時姿態(tài)運動學為
(2)
Tg=Agθ
(3)
(4)
(5)
考慮模型不確定性,矩陣A、B可表示為[16]:
式中:下標0表示為標稱系統(tǒng)下的矩陣;D、F1和F2為適當維數(shù)的常矩陣;Δ為未知矩陣且滿足ΔTΔ≤I。忽略二階不確定項并只考慮轉動慣量的不確定性可得:
(6)
式中:d∈[0,1),當d=0時系統(tǒng)不存在不確定性;
diag[4(J22+J33),3(J11+J33),J11+J22]
當不考慮撓性模態(tài)時,e12、e13、e21、e31均為零矩陣。
證明:首先推導矩陣相加求逆結果,如對于矩陣M有(M+ΔM)-1=Μ-1+N,其中ΔM為不確定項,N為待求解的矩陣。將等式兩邊分別左乘M+ΔM得:
于是
(M+ΔM)-1=Μ-1-(MΔM-1Μ+Μ)-1
同理可推導:
保留ΔM的一次項,近似可得
(M+ΔM)-1=Μ-1-M-1ΔMM-1
那么忽略二階不確定項,考慮轉動慣量的不確定性有:
記
(7)
同理可得:
(8)
(9)
由式(7)~(9)可得式(6)。得證。
(1)PSO算法
PSO算法是一種受鳥群覓食啟發(fā)的隨機搜索算法[17],具有結構簡單、參數(shù)少和搜索速度快等優(yōu)點,已廣泛應用到自然科學和工程科學多個領域。該算法中速度和位置更新如下。
式中:k(k=1,2,…,G)為當前迭代次數(shù),G為最大迭代次數(shù);vi和xi分別為第i個粒子的速度和位置向量,取值區(qū)間分別為[vmin,vmax]和[xmin,xmax],i∈{1,…,N},N為群體個體總數(shù);pi和g分別為個體及全局極值,pbest(i)、gbest為對應的適應度值;c1、c2為加速常數(shù)即學習率;r1、r2為[0,1]范圍隨機數(shù);w(k)=wmax-k(wmax-wmin)/G為慣性權重,取值區(qū)間為[wmin,wmax]。
(2)差分進化(differential evolution, DE)算法
DE算法是一種魯棒性較強的隨機搜索算法,其基本思想來源于遺傳算法,這里選取差分進化算法DE/rand/1/bin策略,種群參數(shù)定義同上述PSO算法,關鍵步驟如下[18-19]。
1)變異,計算變異向量:
υi(k+1)=xl1(k)+F(k)·[xl2(k)-xl3(k)]
2)交叉,計算試驗變量:
式中:交叉算子CR=0.3×[1+rand(0,1)];randb(j)表示產(chǎn)生[0,1]之間隨機數(shù)發(fā)生器的第j個估計值;rnbr(i)∈(1,2,…,D)表示一個隨機選擇的序列,以確保ui至少從υi中獲得一個參數(shù)。
3)選擇,DE算法按照貪婪準則將試驗向量ui(k+1)與當前種群中的目標向量xi(k)進行比較,較優(yōu)的向量將在下一代種群中出現(xiàn)。
(3)DPSO算法
DPSO算法在PSO算法的基礎上加入DE算法的變異操作,在保證粒子群朝著好的方向發(fā)展的同時增加了粒子的多樣性,提高粒子的全局搜索能力,有效避免了早熟至局部最優(yōu)解的情況[20]。算法流程和部分參數(shù)定義同PSO算法,其速度和位置更新如下。
(10)
F(k)定義同DE算法,r3為[0,1]范圍內的隨機數(shù),xU(k)和xL(k)為當代粒子中的較優(yōu)和較劣解,以保證粒子向有益的方向發(fā)展。
由于撓性模態(tài)難以觀測或直接測量,本節(jié)參考文獻[7],采用姿態(tài)角和角速度反饋設計控制器,下面關于多目標約束的描述建立在此基礎上。
考慮轉動慣量的不確定性和未建模動態(tài)、各種干擾對系統(tǒng)的影響,以及避免控制飽和的現(xiàn)象,本節(jié)在控制器設計時加入魯棒D-穩(wěn)定約束,并對控制能量和內外干擾抑制相關指標進行優(yōu)化。極點約束和兩個典型優(yōu)化指標的顯示表達以及Pareto最優(yōu)理論描述如下。
2.2.1 魯棒D-穩(wěn)定
(1)含撓性模態(tài)
設穩(wěn)定裕度參數(shù)為h(h<0),那么式(5)系統(tǒng)對應如圖1所示,魯棒D-穩(wěn)定實現(xiàn)為:
(11)
式中:K∈Rm×n為相應維數(shù)的反饋增益矩陣,在優(yōu)化的過程中視為已知,此時控制輸入u=BKy;ε>0為標量。優(yōu)化時可將標稱系統(tǒng)下的撓性模態(tài)極點限制在一定變化范圍內,如實部區(qū)間為[a1,a2](a1,a2<0),虛部絕對值區(qū)間為[b1,b2](b1,b2>0)。
圖1 極點配置區(qū)域(1)Fig.1 Pole assignment region(1)
證明:參考文獻[16],考慮不確定性和保證系統(tǒng)一定的魯棒穩(wěn)定裕度有:
P[(A0+ΔA)+(B0+ΔB)KC]+
[(A0+ΔA)+(B0+ΔB)KC]TP-hI<0
式中:P為對稱正定矩陣。記Y=P(A0+B0KC)+(A0+B0KC)TP-hI,則有
Y+PDΔ(F1+F2K)+
(F1+F2K)TΔT(PD)T<0
上述矩陣不等式對所有滿足ΔTΔ≤I的不確定矩陣Δ成立,當且僅當存在ε,使得
Y+εPDDTP+ε-1(F1+F2K)T(F1+F2K)<0
應用矩陣的Shur補性質,可以等價為:
上式左右同乘矩陣diag{P-1,I},并記X1=P-1,便可得到式(11),得證。
(2)不含撓性模態(tài)
針對不確定系統(tǒng),對給定區(qū)域D,若存在對稱矩陣X1、X2和標量ε1、ε2>0同時滿足:
圖2 極點配置區(qū)域(2)Fig.2 Pole assignment region(2)
(12)
(13)
式中:A0、B0、C不含撓性模態(tài)信息且此時C為是單位陣;此處X1可不等于X2,ε1可不等于ε2,以減小保守性。則該不確定系統(tǒng)為魯棒D-穩(wěn)定。證明過程同式(11)。
2.2.2 (內外)干擾抑制能力
式中:Ac=A+BKC。建立從w到y(tǒng)的傳遞函數(shù)矩陣Gyw(s),那么根據(jù)H2控制理論思想,用‖Gyw(s)‖2定義干擾對輸出影響的大小[13]:
2.2.3 控制能量
現(xiàn)代航天器多以動量輪或控制力矩陀螺為姿態(tài)控制執(zhí)行機構,為防止因執(zhí)行機構飽和而導致的系統(tǒng)性能變差或失穩(wěn)現(xiàn)象,需加入控制能量約束。直接將飽和非線性因素加入控制器設計會帶來分析困難,所以參考文獻[13],基于航天器姿態(tài)定向保持時姿態(tài)角偏差和姿態(tài)角速度均為小量的特點,通過優(yōu)化‖K‖2的值來盡可能減少控制能量‖u‖。記
f2(K)=‖K‖2
2.2.4 Pareto最優(yōu)原則
一般優(yōu)化指標之間的關系難以確定,且不乏相互沖突的情況,因此無法找到多目標均為最優(yōu)的唯一解,只能找到在多目標間平衡的次優(yōu)解,此時需應用Pareto最優(yōu)解的概念。假設存在性能指標函數(shù)f1(x),…,fs(x)(s≥2),其中x為搜索空間φ中的向量,與K對應,有如下定義:
x=[K(1,1),…,K(1,n),K(2,1),…,
K(2,n),…,K(m,n)]
Pareto最優(yōu)解:如果不存在同時滿足下式的x∈φ,那么x*∈φ則為Pareto最優(yōu)解或非劣解(此解不唯一)[21]:
(14)
(15)
如果存在滿足上述條件的解x,那么則稱該解Pareto占優(yōu)于解x*。
創(chuàng)建由粒子每代更新產(chǎn)生的非劣解和占優(yōu)解組合而成的外部解集;因為想得到每個指標都相對較優(yōu)而不是某單個指標最優(yōu)的解,所以為了提高算法的快速性,不按擁擠度法[12]而是按距離評價法更新該解集,保留距離較小的解。做歸一化處理,對應x的距離d(x)定義如下:
式中:maxfs和minfs分別表示第s個目標函數(shù)在外部解集中取得的最大和最小值。以s=2為例,如圖3所示。
圖3 距離評價法Fig.3 Distance evaluation
基于DPSO算法的控制器設計步驟如下。
步驟1:初始化粒子群。包括算法參數(shù)和位置、速度向量,其中初始位置向量由PID控制下的增益矩陣和LMI方法下求得的幾組反饋增益矩陣定義,初始速度向量為零向量;初始化粒子i(i=1,2,…,N)的歷史最優(yōu)位置pi和最優(yōu)值p1best(i)和p2best(i)并選取第一個粒子為全局最優(yōu)粒子,對應位置g和最優(yōu)值g1best和g2best。迭代次數(shù)k=1,且xU(k)=xL(k)=x1(k)。
步驟2:更新粒子速度和位置。k=k+1;按式(10)更新速度和位置,并做邊界條件處理;判斷是否滿足魯棒-D穩(wěn)定極點約束,若滿足進入步驟3,否則進入步驟6。
步驟3:更新個體粒子最優(yōu)解和全局最優(yōu)解。計算個體粒子適應度值(性能指標)f1(xi),f2(xi),按Pareto占優(yōu)更新最優(yōu)解;將每個粒子當代更新的非劣解和占優(yōu)解存入外部解集。
步驟4:更新外部解集。計算解集中每個解對應的距離,找到距離最短的解,然后將其他解對應的適應度值與該解比較,保留下非劣解和占優(yōu)解;判斷保留下的解的個數(shù)是否大于集合上限N,若是,進入步驟5,否則剔除距離較大的部分解。
步驟5:計算xU(k)和xL(k)。計算各粒子的適應度值f1(xi),f2(xi),并求取所有粒子適應度值的平均值;將每個粒子對應的適應度值與平均值比較,在Pareto占優(yōu)的粒子中任選一個將其解作為xU(k),在剩余粒子中任選一個將其解作為xL(k);當不存在Pareto占優(yōu)的粒子時,xU(k)=xL(k),取任意粒子位置向量。
步驟6:判斷是否滿足k≤G,如滿足返回步驟2,否則進入步驟7。
步驟7:畫圖,并定義仿真參數(shù)。
以某大型航天器為例,取其一對翼板的前兩階撓性模態(tài)進行控制器設計和仿真。標稱系統(tǒng)轉動慣量、耦合系數(shù)、帆板撓性約束模態(tài)角頻率對角陣和阻尼系數(shù)分別為:
圖2中參數(shù)α=50°,h1=-0.2,h2=-0.02;撓性模態(tài)極點相關參數(shù)h=-0.001,a1=-0.03,a2=-0.001 5,b1=0.25,b2=0.35。優(yōu)化算法參數(shù)取值如表1所示;以控制能量為適應度值為例,PSO、DE和DPSO三種算法的對比優(yōu)化結果見圖4。
從圖4可看出DE算法進化較慢,PSO算法過早收斂,容易陷入局部極值,本文提出的DPSO算法尋優(yōu)能力最好,因此采用DPSO算法進行多目標優(yōu)化。圖5為基于Pareto最優(yōu)原則的多目標優(yōu)化結果,可以看出在距離評價法下,外部解集中的非劣解集合接近Pareto前沿[12]。
表1 三種算法下參數(shù)設置
圖4 三種算法下進化曲線Fig.4 Evolution curve under three algorithms
圖5 外部解集Fig.5 Set of external solutions
為驗證本文控制器的魯棒性能、干擾抑制能力和控制能量優(yōu)化效果,仿真時不妨取轉動慣量(不同于J)和外部干擾為:
圖8所示為三軸姿態(tài)受到的干擾作用,圖9~圖11為閉環(huán)系統(tǒng)仿真結果。由圖9知航天器姿態(tài)穩(wěn)定時間相對于PID控制減小,基本不超過100 s;在300 s后,PID控制姿態(tài)誤差在1.35×10-4(°)以內,而本文方法下姿態(tài)誤差在6.2×10-5(°)以內,減小了約54%,說明對干擾進行了有效抑制;圖10顯示在本文仿真算例下,兩種方法控制力矩相差不大,本文方法控制力矩峰值約為6.99N·m,PID控制峰值約為8.42N·m;由圖11知,本文控制方法下的航天器前兩階撓性模態(tài)振動幅度與PID控制相當,且處于較小的量級,其中與姿態(tài)耦合作用較強的第一階撓性模態(tài)振動衰減更快,在100 s左右接近穩(wěn)定(比PID快約250 s)。綜上,本文方法對航天器高精度高穩(wěn)定度魯棒姿態(tài)控制具有良好的效果。
圖 6 極點配置結果(圖中線條與圖2對應)Fig.6 Result of pole assignment
圖7 開閉環(huán)系統(tǒng)奇異值曲線Fig.7 Singular value of open-loop and closed-loop system
圖8 姿態(tài)干擾Fig.8 Attitude interference
圖9 姿態(tài)角Fig.9 Attitude angle
圖10 控制力矩Fig10 Control torque
圖11 撓性模態(tài)位移Fig.11 Flexible modal displacement
本文針對大型撓性航天器的高精度高穩(wěn)定度姿態(tài)控制問題,將改進智能算法、現(xiàn)代控制理論和多目標Pareto最優(yōu)思想結合,提出一種基于輸出反饋和DPSO算法的多目標魯棒控制器設計方法,并推導了含轉動慣量不確定性的系統(tǒng)模型和魯棒D-穩(wěn)定的LMI表達。主要結論如下:
1)相較于DE算法和PSO算法,本文提出的DPSO算法尋優(yōu)效果更好。在魯棒區(qū)域極點約束下,利用DPSO算法和Pareto最優(yōu)原則,對控制能量和魯棒性能進行優(yōu)化得到的控制器實現(xiàn)了閉環(huán)系統(tǒng)的魯棒D-穩(wěn)定,且起到了一定的振動抑制作用。
2)仿真驗證了本文方法的有效性。在本算例中,相對于PID控制,本文方法下系統(tǒng)響應速度較快,其中第1階撓性模態(tài)振動衰減時間減少約250 s,控制力矩峰值減小約1.43N·m,姿態(tài)穩(wěn)態(tài)誤差減小約54%。
3)將智能優(yōu)化與現(xiàn)代控制理論結合,涉及的數(shù)學公式和參數(shù)較少,靈活性較強。本文方法彌補了經(jīng)典PID控制在不確定系統(tǒng)模型下魯棒性能以及抗干擾能力的不足;減小了帶極點配置的LMI方法在多指標約束凸優(yōu)化問題中的求解保守性;也能避免傳統(tǒng)方法將多目標轉化為加權指標函數(shù)時,由于多目標關系難以確定而導致的加權系數(shù)選擇困難。
4)本文智能優(yōu)化方法還可應用到更多復雜以及非線性系統(tǒng)控制參數(shù)優(yōu)化問題上,為大型撓性航天器魯棒控制器設計提供一種思路。