陳竑宇 陳 提
(南京航空航天大學(xué)航空學(xué)院,南京 210016)
姿態(tài)控制系統(tǒng)是航天器的重要組成部分.用于描述姿態(tài)運(yùn)動(dòng)的旋轉(zhuǎn)矩陣的構(gòu)型空間為Lie 群SO(3),由于不是歐幾里得空間,為了便于姿態(tài)控制的研究,學(xué)者們提出了多種姿態(tài)表示方法[1].目前常用的姿態(tài)表示方法主要包括歐拉角,修正羅德里格參數(shù)(modified Rodrigues parameters,MRPs)和單位四元數(shù)等.但是這些姿態(tài)表示方法均有一定局限性,比如歐拉角和MRPs 存在奇異點(diǎn),不能夠全局描述姿態(tài)運(yùn)動(dòng);而在單位四元數(shù)表示中,對(duì)于每一種旋轉(zhuǎn)都存在兩組單位四元數(shù)與之相對(duì)應(yīng),因此在姿態(tài)表示上并不唯一.所以,基于單位四元數(shù)的連續(xù)反饋控制率可能會(huì)產(chǎn)生不期望的退繞(unwinding)現(xiàn)象[1-2].因此,已有許多研究嘗試直接在SO(3)上解決姿態(tài)控制問題.在SO(3)上剛體姿態(tài)動(dòng)力學(xué)由12 個(gè)非線性方程描述,但通常有3 個(gè)控制輸入,是典型的一個(gè)多輸入多輸出非線性系統(tǒng).Lee 定義了一個(gè)新的姿態(tài)誤差函數(shù),并設(shè)計(jì)了指數(shù)穩(wěn)定控制器[3].Sanyal等[4]討論了由旋轉(zhuǎn)矩陣描述的剛體姿態(tài)估計(jì)和控制問題.由于時(shí)不變光滑反饋律最多僅能保證SO(3)上姿態(tài)控制的近全局收斂,Berkane 等[5]提出了一種基于中心協(xié)同混合方法的全局指數(shù)穩(wěn)定控制方案.Wang 等[6]研究了角速度為分段正弦信號(hào)的SO(3)的姿態(tài)控制.在實(shí)際工程中,航天器通常帶有撓性構(gòu)件,因此不能忽略撓性振動(dòng)與姿態(tài)運(yùn)動(dòng)的耦合關(guān)系,龔柯杰[7]將航天器撓性附件振動(dòng)視作外界擾動(dòng)項(xiàng),并設(shè)計(jì)了自適應(yīng)滑??刂破?He 等[8]使用假設(shè)模態(tài)法建立SO(3)上撓性航天器姿態(tài)動(dòng)力學(xué)模型,提出了一種迭代學(xué)習(xí)觀測器用于觀測噪聲,并基于觀測器設(shè)計(jì)了一個(gè)復(fù)合控制器來穩(wěn)定姿態(tài).在這些研究中,控制器設(shè)計(jì)完全或部分取決于建模,但撓性航天器姿態(tài)動(dòng)力學(xué)通常具有強(qiáng)非線性、強(qiáng)耦合、參數(shù)不確定性,模態(tài)測量困難等特點(diǎn),此外,大尺寸撓性航天器存在模態(tài)集中和“全局運(yùn)動(dòng)疊加全局模態(tài)”等復(fù)雜現(xiàn)象[9],因此,建立撓性航天器精確動(dòng)力學(xué)模型是一個(gè)極其困難的任務(wù).
對(duì)于復(fù)雜的非線性系統(tǒng),為了分析的方便,往往會(huì)考慮將其線性化.典型的方法包括泰勒線性化、分段線性化和正交函數(shù)逼近線性化等.但這些方法存在明顯局限,泰勒線性化只在工作點(diǎn)附近進(jìn)行,分段線性化需要更多模型信息,正交函數(shù)逼近線性化計(jì)算量大且控制器設(shè)計(jì)復(fù)雜.1931 年,Koopman[10]提出非線性動(dòng)力系統(tǒng)可以通過作用于測量函數(shù)上的無限維線性算子表示,即提出了一種非線性系統(tǒng)的全局線性化方法.隨著近些年計(jì)算機(jī)技術(shù)和大數(shù)據(jù)的發(fā)展,Koopman 算子理論被越來越多地運(yùn)用于復(fù)雜非線性系統(tǒng)行為分析[11-14].一些研究已經(jīng)提出了一些低階非線性系統(tǒng)的解析觀測函數(shù)集[12,15],但對(duì)于更復(fù)雜的非線性動(dòng)力學(xué)或者未知?jiǎng)恿W(xué)系統(tǒng),則難以解析地得到Koopman 算子的表示.因此可以采取數(shù)據(jù)驅(qū)動(dòng)方法被辨識(shí)Koopman 算子,這種方法在建模過程中只需要仿真或者實(shí)驗(yàn)數(shù)據(jù),而不像傳統(tǒng)方法那樣依賴于物理參數(shù).由于Koopman 算子通常是無限維算子,在實(shí)際應(yīng)用中往往需要得到其有限維截?cái)?Schmid[16]在流體動(dòng)力學(xué)領(lǐng)域提出了動(dòng)態(tài)模態(tài)分解算法(dynamics mode decomposition,DMD),隨后Rowley 等[17]在DMD 和Koopman 理論之間建立了重要聯(lián)系.目前,辨識(shí)Koopman 算子常用的算法包括DMD 算法[16]和拓展動(dòng)態(tài)模態(tài)分解(EDMD)算法[18]等.Dahdah 等[19]針對(duì)現(xiàn)有算法辨識(shí)開環(huán)系統(tǒng)Koopman 算子的局限性,結(jié)合輸入輸出子空間辨識(shí)法和EDMD 算法得到了閉環(huán)系統(tǒng)Koopman 算子.由于DMD 算法和EDMD 算法直接計(jì)算Koopman算子有限維截?cái)?截?cái)嗾`差會(huì)不可避免影響建模精度.Brunton 等[20]提出非線性動(dòng)力學(xué)的稀疏性辨識(shí)(SINDY),通過稀疏回歸算法精確辨識(shí)非線性動(dòng)力學(xué)的廣義線性模型.Kaiser 等[21]在少量數(shù)據(jù)基礎(chǔ)上基于SINDY 算法辨識(shí)非線性動(dòng)力學(xué),驗(yàn)證了SINDY算法的準(zhǔn)確性,并基于模型使用非線性模型預(yù)測控制器(NMPC)實(shí)現(xiàn)系統(tǒng)控制.
近年來,許多研究人員開始關(guān)注基于Koopman算子理論的非線性系統(tǒng)控制問題,將成熟的線性控制理論用于解決非線性系統(tǒng)控制問題.Kaiser 等[22]基于一些低階非線性系統(tǒng),討論了基于Koopman 算子的控制問題.Ren 等[23]基于Koopman 理論設(shè)計(jì)了輪式機(jī)器人的滑模控制器.Korda 等[24]在Koopman算子理論的基礎(chǔ)上采用線性模型預(yù)測控制對(duì)非線性系統(tǒng)進(jìn)行控制.為提高閉環(huán)Koopman MPC 的魯棒性,Zhang 等[25]提出一種基于Tube 模型預(yù)測控制(tube-based MPC)的魯棒Koopman 控制器.Goyal 等[26]提出了一種基于Koopman 算子和數(shù)據(jù)驅(qū)動(dòng)方法的三自由度腕部康復(fù)機(jī)器人的軌跡跟蹤控制器.
目前也已有一些研究將Koopman 算子理論用于SO(3)上的姿態(tài)控制.比如,Chen 等[27]基于SO(3)構(gòu)型提出了一組非線性觀測函數(shù),得到了SO(3)上姿態(tài)動(dòng)力學(xué)Koopman 算子的解析表示,并在此基礎(chǔ)上設(shè)計(jì)了最優(yōu)控制器.Zinage 等[28]考慮了剛體位姿控制,將Koopman 建模的研究推廣到SE(3)群上的無人機(jī)動(dòng)力學(xué).這些研究主要是針對(duì)剛體姿態(tài)動(dòng)力學(xué)的建模與控制,但由于實(shí)際中航天器帶有撓性構(gòu)件,不能忽略撓性振動(dòng)與姿態(tài)運(yùn)動(dòng)的耦合關(guān)系,因此SO(3)上撓性航天器姿態(tài)動(dòng)力學(xué)的Koopman 辨識(shí)會(huì)更加復(fù)雜.
本文擬提出一組觀測函數(shù)用于Koopman 算子的數(shù)據(jù)驅(qū)動(dòng)辨識(shí),這組觀測函數(shù)僅依賴于狀態(tài)測量數(shù)據(jù),而不需要任何先驗(yàn)信息.在此基礎(chǔ)上,本文基于SINDY 算法得到SO(3)上撓性姿態(tài)動(dòng)力學(xué)的稀疏化廣義線性模型和Koopman 算子的有限維截?cái)?實(shí)現(xiàn)對(duì)撓性航天器姿態(tài)動(dòng)力學(xué)的全局線性化建模.基于截?cái)嗟腒oopman 算子設(shè)計(jì)一種高效率的線性LQR 控制器,在實(shí)現(xiàn)撓性航天器姿態(tài)控制的同時(shí)抑制撓性構(gòu)件振動(dòng),最后通過仿真算例驗(yàn)證所提出的數(shù)據(jù)驅(qū)動(dòng)建模方法和最優(yōu)控制器,并和傳統(tǒng)最優(yōu)控制方法進(jìn)行對(duì)比以論證所提方法的優(yōu)勢.
Koopman 算子理論可以將有限維非線性系統(tǒng)轉(zhuǎn)換為無窮維線性系統(tǒng).考慮如下非線性動(dòng)力系統(tǒng)
式中x∈S為光滑流形上的N維線性系統(tǒng)的狀態(tài),S?RN為狀態(tài)空間.
定義實(shí)值觀測函數(shù)g:M→R,它是Hilbert 空間上的元素.Koopman 算子 Kt是一個(gè)無限維線性算子,它對(duì)測量函數(shù)g的作用如下
其中 ? 為合成算子.
若式(1)中系統(tǒng)變?yōu)殡x散時(shí)間形式,那么對(duì)于時(shí)間步長為 Δt的離散時(shí)間系統(tǒng)x(k+1)=FΔt(x(k)),式(2)變?yōu)?/p>
即,Koopman 算子定義了一個(gè)無限維線性動(dòng)力系統(tǒng),并使?fàn)顟B(tài)觀測值gk=g(xk) 推進(jìn)到下一時(shí)間步長,即測量函數(shù)的演化[29].
Koopman 算子的線性屬性來源于觀測函數(shù)空間中加法運(yùn)算的線性,即
其中 α1,α2為常數(shù).
Koopman 算子定義了一個(gè)無窮維線性系統(tǒng)來表示非線性系統(tǒng),但無窮維的性質(zhì)會(huì)給Koopman 算子的應(yīng)用帶來問題.在實(shí)際應(yīng)用中,無窮維系統(tǒng)通常可以被截?cái)?如果不產(chǎn)生過大的精度損失,那么就可以通過有限維線性系統(tǒng)來近似Koopman 算子.
SO(3)上的剛性航天器姿態(tài)動(dòng)力學(xué)方程可表示為如下形式
其中R∈S O(3) 表示從航天器體坐標(biāo)系到慣性參考系的旋轉(zhuǎn)矩陣,ω ∈R3是體坐標(biāo)系的角速度,J∈R3×3為航天器慣性矩陣.對(duì)于任意向量a∈R3,a×定義為
在本文研究中,假設(shè)撓性航天器的結(jié)構(gòu)為中心剛體上固接梁,其構(gòu)型如圖1 所示.可采取假設(shè)模態(tài)法對(duì)撓性航天器進(jìn)行建模.SO(3)上的撓性航天器姿態(tài)動(dòng)力學(xué)方程可以寫為式(5)以及如下方程
圖1 撓性航天器示意圖Fig.1 Flexible spacecraft
考慮到撓性構(gòu)件模態(tài)坐標(biāo)一般不能直接測量,但可以測量撓性構(gòu)件上各點(diǎn)位移以及速度.撓性構(gòu)件位移w以及速度與模態(tài)坐標(biāo)q和之間的關(guān)系為
其中,Φ∈RN×m為各模態(tài)對(duì)應(yīng)固有振型矩陣,N為傳感器數(shù)量,且N>m.因此,可以使用替代模態(tài)坐標(biāo)參數(shù),兩者之間關(guān)系為
原始SINDY 算法[20]主要思想是利用非線性動(dòng)力學(xué)中通常只包含少數(shù)起作用項(xiàng)的性質(zhì),從人工選取的字典函數(shù)庫中尋找動(dòng)力學(xué)中的有效項(xiàng).考慮如下受控非線性動(dòng)力學(xué)系統(tǒng)
收集時(shí)序數(shù)據(jù)得到數(shù)據(jù)快照矩陣
建立候選非線性函數(shù)庫
其中 θi為人工選取的候選非線性函數(shù),以保證Θ(X,U) 中包含非線性動(dòng)力學(xué)f(x,u) 中的有效項(xiàng).通過以下廣義線性模型來逼近非線性動(dòng)力學(xué)f(x,u)
式中,Ξ 的列向量 ξk是一個(gè)系數(shù)向量,用于確定動(dòng)力學(xué)中的有效項(xiàng),且 ξk中的非零項(xiàng)要盡可能少.通過稀疏回歸算法辨識(shí)動(dòng)力學(xué)對(duì)應(yīng)有效項(xiàng)
將系統(tǒng)原始狀態(tài)映射到高維觀測函數(shù)空間,再辨識(shí)觀測空間中動(dòng)力學(xué)的有效項(xiàng),便可以得到在觀測函數(shù)空間上非線性動(dòng)力學(xué)的廣義線性模型,其中的線性部分便能作為Koopman 算子的有限維截?cái)?
選取構(gòu)成Koopman 提升空間的非線性觀測函數(shù) Ψ ∈RN和控制耦合函數(shù) Ψu∈RM以及高階補(bǔ)償觀測函數(shù)如下
基于數(shù)據(jù)可得到提升的數(shù)據(jù)快照矩陣為
通過稀疏回歸算法可以得到提升空間上的廣義線性模型如下所示
為了后續(xù)討論方便,將非線性系統(tǒng)的近似廣義線性模型,即式(25),寫為狀態(tài)空間形式,即
即得到非線性動(dòng)力學(xué)Koopman 算子的有限維截?cái)?
目前常用的稀疏回歸算法包括套索回歸(LASSO)或者序貫閾值最小二乘(sequential thresholded leastsquares,STLS)等[20].本文采取序貫閾值最小二乘算法(STLS)來計(jì)算SINDY 模型.STLS 算法基本思想是以稀疏提升系數(shù)為閾值將最小二乘解稀疏化,其具體計(jì)算過程如下:
(3)將 Ξ 中所有小于稀疏提升系數(shù) λ 的系數(shù)置零.矩陣Bi表示 Ξ 中非零元素位置,其中非零項(xiàng)對(duì)應(yīng)位置元素為1,其他位置元素為0;
(4)通過最小二乘法重新計(jì)算稀疏最優(yōu)擬合算子 Ξ,但是將僅使用非稀疏項(xiàng)所對(duì)應(yīng)的數(shù)據(jù),即ξk=,其中 ξk為 Ξ 的第k列,Bik為Bi的第k列,·為點(diǎn)乘運(yùn)算;
(5)重復(fù)步驟(2)~ (4),直到達(dá)到最大迭代步數(shù)或者解收斂.
對(duì)于線性系統(tǒng),其性能泛函是狀態(tài)變量和控制變量的二次型函數(shù)的積分,則這樣的最優(yōu)控制問題稱為線性二次型最優(yōu)控制問題,線性二次型調(diào)節(jié)器(LQR)可以得到狀態(tài)線性反饋的閉環(huán)最優(yōu)控制率[30].由于求解方便,LQR 的使用非常廣泛,是現(xiàn)代控制理論中最重要成果之一.Koopman 算子理論可以將本文所關(guān)心的撓性航天器的非線性控制問題映射到線性控制,因此可以采用LQR 方法設(shè)計(jì)相應(yīng)的控制策略.本節(jié)基于SINDY 辨識(shí)得到的Koopman 算子有限維近似模型,提出了一種LQR 控制器,在后續(xù)內(nèi)容中稱其為Koopman-LQR 控制器.
其中QLQR和RLQR為正定矩陣.Koopman-LQR 控制增益F可以通過求解代數(shù)Riccati 方程得到
Koopman-LQR 控制輸入便可以寫為
注意到Koopman-LQR 控制輸入 ΨuLQR需要進(jìn)行逆變換才能得到實(shí)際輸入到撓性姿態(tài)動(dòng)力學(xué)系統(tǒng)中的控制力矩,且根據(jù)非線性觀測項(xiàng)的選擇不同,Ψu可能為狀態(tài)依賴量.因此需要通過狀態(tài)依賴轉(zhuǎn)換矩陣T(x) 計(jì)算實(shí)際輸入的控制力矩u,如下所示
其中T(x) 為將u映射到 Ψu的變換矩陣,其形式基于Ψu的選取,上標(biāo) ? 表示偽逆.變換矩陣T(x) 通??梢曰谏弦徊蓸訒r(shí)刻的狀態(tài)測量計(jì)算得到.
Chen 等[27]提出了一組觀測函數(shù),并解析得到了SO(3)上剛性航天器姿態(tài)動(dòng)力學(xué)的Koopman 算子,但是所提出的觀測函數(shù)并不適用于撓性情況,且所提觀測函數(shù)中控制輸入的計(jì)算要求已知航天器轉(zhuǎn)動(dòng)慣量矩陣.本文根據(jù)撓性航天器姿態(tài)動(dòng)力學(xué)所包含的非線性成分,提出了一組用于SINDY 算法辨識(shí)SO(3)上撓性航天器姿態(tài)動(dòng)力學(xué)Koopman 算子有限維截?cái)嗟姆蔷€性觀測函數(shù).
對(duì)于帶撓性構(gòu)件的航天器姿態(tài)動(dòng)力學(xué),非線性觀測選取為如下
因此,基于Koopman 算子理論,在小角速度情況下,撓性航天器姿態(tài)動(dòng)力學(xué)可近似地表示為
由于同時(shí)需要足夠大的階數(shù)保證線性化系統(tǒng)有足夠的精度,以及避免觀測函數(shù)項(xiàng)過多而產(chǎn)生較大的回歸誤差,在之后的仿真中都將式(33)~ 式(37)中的觀測函數(shù)階數(shù)n取為3.
在仿真算例中,假設(shè)有3 個(gè)模態(tài)與航天器轉(zhuǎn)動(dòng)耦合,其他模態(tài)和轉(zhuǎn)動(dòng)的耦合系數(shù)均為零,并且假設(shè)撓性構(gòu)件上有3 組位移和加速度傳感器,即模態(tài)振型矩陣 Φ 可逆,則1.2 節(jié)中式(11)成立.
SINDY 是一種數(shù)據(jù)驅(qū)動(dòng)算法,因此需要采集訓(xùn)練數(shù)據(jù).基于1.2 節(jié)中的撓性航天器姿態(tài)動(dòng)力學(xué)數(shù)學(xué)模型,在MATLAB 中使用ODE45 進(jìn)行數(shù)值仿真,將所得仿真數(shù)據(jù)作為訓(xùn)練數(shù)據(jù).假設(shè)有3 個(gè)模態(tài)與航天器轉(zhuǎn)動(dòng)耦合,只給出前3 階模態(tài)參數(shù).撓性航天器的參數(shù)如表1 所示.
表1 航天器參數(shù)Table 1 Parameters of spacecraft
數(shù)據(jù)的采樣間隔取為0.01 s,初始條件為: 隨機(jī)的初始旋轉(zhuǎn)矩陣R(0) 和隨機(jī)初始角速度ω(0)∈[-0.01,0.01] rad/s,撓性構(gòu)件初始模態(tài)坐標(biāo)參數(shù)為q(0)=0,(0)=0 .將控制力矩設(shè)為如下形式的正弦激勵(lì)
其中aj∈[0,1],bj∈[0,2],cj∈[0,2π],j=1,2,3,Amp=0.1.
將在以上隨機(jī)初始條件和控制力矩下得到的160 組,時(shí)間為4 s,數(shù)據(jù)長度為401 的狀態(tài)軌跡用作訓(xùn)練數(shù)據(jù).SINDY 算法所需要的導(dǎo)數(shù)數(shù)據(jù)由五點(diǎn)法數(shù)值微分得到,因此實(shí)際上訓(xùn)練數(shù)據(jù)為160 組數(shù)據(jù),每組包含396 組采樣.
取SINDY 算法中稀疏提升系數(shù)為λ=1×10-5.基于數(shù)據(jù)得到SINDY 模型,并使用均方根誤差來評(píng)估模型對(duì)訓(xùn)練數(shù)據(jù)的擬合程度
其中‖ ‖2表示矩陣的2 范數(shù).
SINDY 模型對(duì)訓(xùn)練數(shù)據(jù)的均方根誤差為RMSEtrain≈1.2×10-9.另外將40 組長度4 s 的隨機(jī)初始條件出發(fā)的采樣數(shù)據(jù)作為驗(yàn)證數(shù)據(jù)集,模型對(duì)驗(yàn)證數(shù)據(jù)集的均方根誤差為RMSEverify≈1.5×10-9.結(jié)果表明SINDY 模型對(duì)訓(xùn)練數(shù)據(jù)的擬合精度高,且具備良好的泛化能力.
為驗(yàn)證數(shù)據(jù)驅(qū)動(dòng)模型的預(yù)測精度,進(jìn)行20 組隨機(jī)初始條件下的仿真,其結(jié)果與數(shù)值解的誤差表示了模型預(yù)測精度.由于需要驗(yàn)證撓性構(gòu)件存在初始振動(dòng)的情況下的數(shù)據(jù)驅(qū)動(dòng)模型的預(yù)測精度,將初始模態(tài)坐標(biāo)設(shè)為非0,其范圍如下
角速度預(yù)測誤差由角速度中各個(gè)元素的平均誤差Errω來表示
旋轉(zhuǎn)矩陣預(yù)測誤差由旋轉(zhuǎn)矩陣中各個(gè)元素的平均誤差ErrR來表示
撓性航天器姿態(tài)動(dòng)力學(xué)SINDY 模型的50 s 內(nèi)的預(yù)測誤差如圖2 所示.圖2(a)顯示角速度預(yù)測誤差在大約 ±1.6×10-5rad/s 范圍內(nèi),圖2(b)顯示50 s內(nèi)旋轉(zhuǎn)矩陣預(yù)測誤差在 ±1×10-3范圍內(nèi).
圖2 SINDY 模型預(yù)測誤差Fig.2 Predict error of SINDY model
結(jié)果表明SINDY 算法準(zhǔn)確辨識(shí)到了姿態(tài)動(dòng)力學(xué)系統(tǒng)在提升空間中對(duì)應(yīng)的非線性成分,證明了數(shù)據(jù)驅(qū)動(dòng)建模的有效性.而對(duì)于撓性航天器姿態(tài)動(dòng)力學(xué)Koopman 算子有限維截?cái)?即式(27),也需驗(yàn)證其預(yù)測能力.Koopman 算子有限維截?cái)嗟念A(yù)測誤差如圖3 所示,結(jié)果表明由于忽略了更高階非線性項(xiàng)成分,誤差相對(duì)SINDY 模型有明顯增大.
圖3 Koopman 模型預(yù)測誤差Fig.3 Predict error of Koopman model
由于回歸誤差和軟件計(jì)算誤差等因素,SINDY模型不能完全精確地得到非線性動(dòng)力學(xué)中有效項(xiàng)的系數(shù),特別是其對(duì)應(yīng)系數(shù)小的情況,這是影響模型精度的主要原因.事實(shí)上,非線性函數(shù)項(xiàng)越多,在回歸過程中產(chǎn)生的誤差越大.注意到無論是SINDY 模型還是截?cái)嗟腒oopman 算子,其角速度的預(yù)測精度都要好于旋轉(zhuǎn)矩陣的預(yù)測精度,這是由于表示旋轉(zhuǎn)矩陣的非線性觀測項(xiàng)要多于表示角速度的項(xiàng),并產(chǎn)生了相對(duì)更大的回歸誤差.因此,簡單地增加更多觀測函數(shù)項(xiàng)并不一定能夠提高模型精度.
基于4.2 節(jié)得到的線性化模型,根據(jù)第2 章內(nèi)容設(shè)計(jì)LQR 控制器.本節(jié)中仿真的采樣間隔設(shè)為0.2 s,在Koopman-LQR 作用下,進(jìn)行20 組時(shí)長為200 s的隨機(jī)初始條件出發(fā)的仿真.
Koopman-LQR 控制器的代價(jià)函數(shù)參數(shù)可以寫為如下形式
其中權(quán)重矩陣QR,Qω,Qw對(duì)應(yīng)旋轉(zhuǎn)矩陣R,角速度ω和撓性件位移 δ,QRω對(duì)應(yīng)R和 ω 耦合項(xiàng),RΨu對(duì)應(yīng)控制耦合觀測 Ψu.根據(jù)第3 章中的觀測函數(shù)的選取,令QR=100I9,Qω=20 000I3,QRω=0,Qw=1000I6,RΨu=I138.
定義姿態(tài)誤差指標(biāo)函數(shù)為
定義撓性構(gòu)件振動(dòng)的指標(biāo)函數(shù)為
其中,q為模態(tài)坐標(biāo),m為模態(tài)坐標(biāo)的數(shù)量.
20組仿真結(jié)果如圖4 所示,圖4(a)為姿態(tài)跟蹤誤差,結(jié)果表明在不同初始條件下Koopman-LQR 控制器均能在120 s 內(nèi)將撓性航天器姿態(tài)誤差V1控制到 1×10-3范圍內(nèi).圖4(b)顯示Eq在40 s 內(nèi)快速減小到0.1 內(nèi),其中子圖表明撓性構(gòu)件振動(dòng)在大約100 s前被抑制到很小范圍內(nèi).仿真結(jié)果表明在SINDY 辨識(shí)的Koopman-LQR 控制器能夠?qū)⒑教炱髯藨B(tài)有效控制到目標(biāo)姿態(tài),并且抑制住撓性振動(dòng).
圖4 隨機(jī)初始條件下的20 組仿真Fig.4 20 simulation results with random initial condition
仿真使用的計(jì)算機(jī)處理器為Intel i7-11700.在仿真過程中,Koopman-LQR 控制器平均單步所需計(jì)算時(shí)間大約為1.7 ms,實(shí)際上可以采用更高的采樣頻率.
此外,為了解代價(jià)函數(shù)參數(shù)對(duì)Koopman-LQR 控制器性能的影響,從同一初始條件出發(fā),對(duì)比不同代價(jià)函數(shù)參數(shù)的Koopman-LQR 作用下的結(jié)果.共7 種不同的QLQR和RLQR的選取方案,如表2 所示,其中情況3 與先前的仿真中的選取一致.從相同的初始條件出發(fā)進(jìn)行仿真,姿態(tài)跟蹤誤差V1如圖5(a)所示,撓性構(gòu)件振動(dòng)指標(biāo)函數(shù)Eq如圖5(b)所示.
表2 不同的 QLQR 和 RLQR 選取Table 2 Different chosen of QLQR andRLQR
圖5 不同 QLQR 和 RLQR 對(duì)Koopman-LQR 姿態(tài)控制的影響Fig.5 Influence of different QLQR andRLQR
觀察圖5(a) 和圖5(b),表2 中情況3 與情況6 比較表明,增大QLQR中旋轉(zhuǎn)矩陣元素對(duì)應(yīng)權(quán)重系數(shù)可以減小控制時(shí)間;情況3 和情況5 對(duì)比表明增大QLQR中角速度 ω 對(duì)應(yīng)權(quán)重系數(shù),能夠限制控制器的響應(yīng)速度,獲得更平滑的控制效果,但意味著需要更多時(shí)間才能控制到期望姿態(tài).通過選取不同的RLQR可以調(diào)節(jié)控制時(shí)間,但如圖5(a)中情況6 和情況7 所示,RLQR過大可能會(huì)導(dǎo)致控制器出現(xiàn)超調(diào)現(xiàn)象.從圖5(c)可以看出,情況1,2,3 在100 s 時(shí)撓性構(gòu)件振動(dòng)顯著大于情況4,這意味著增大QLQR中撓性構(gòu)件模態(tài)坐標(biāo)對(duì)應(yīng)參數(shù)能增強(qiáng)Koopman-LQR的振動(dòng)抑制能力.
除了進(jìn)行姿態(tài)控制,在協(xié)同或同步任務(wù)中,撓性航天器往往需要與目標(biāo)航天器或相鄰航天器的運(yùn)動(dòng)同步,這就需要控制器既能夠高精度地跟蹤預(yù)定軌跡,也能夠抑制撓性構(gòu)件的振動(dòng).因此對(duì)本文所提控制器跟蹤時(shí)變期望姿態(tài)的情況進(jìn)行仿真.
假設(shè)期望姿態(tài)的軌跡為滿足初始時(shí)刻姿態(tài)為R(0)=I3的剛性航天器運(yùn)行軌跡,其角速度滿足如下方程
撓性航天器初始姿態(tài)矩陣為I3,初始角速度為rad/s,初始模態(tài)坐標(biāo)為q=0,=0 .實(shí)際上,姿態(tài)跟蹤性能也受到代價(jià)函數(shù)參數(shù)的影響.將LQR 代價(jià)函數(shù)參數(shù)選取為表2 中的其中情況,分別進(jìn)行長度為500 s 的仿真.
跟蹤時(shí)變期望姿態(tài)的仿真結(jié)果如圖6 所示,圖6(a)為角速度 ω1的跟蹤誤差,而圖6(b)為旋轉(zhuǎn)矩陣的跟蹤誤差.結(jié)果表明,在不同的代價(jià)函數(shù)參數(shù)選取情況下,Koopman-LQR 能夠很好地抑制振動(dòng)并跟蹤期望姿態(tài),姿態(tài)跟蹤誤差大約在 [0,2.5×10-3] 范圍內(nèi),而角速度 ω 中第一個(gè)元素 ω1的跟蹤誤差大約在±2×10-3rad/s范圍內(nèi).
圖6 不同代價(jià)函數(shù)參數(shù) QLQR 和 RLQR 對(duì)Koopman-LQR 姿態(tài)跟蹤效果的影響Fig.6 Influence of QLQR and RLQR on attitude tracking
將圖6 中情況1,2,3 和情況4 的結(jié)果相對(duì)比,發(fā)現(xiàn)撓性構(gòu)件參數(shù) δ 在代價(jià)函數(shù)中權(quán)重越高,振動(dòng)抑制效果越好,但姿態(tài)跟蹤誤差增大.將情況3 與情況5 以及情況6 的對(duì)比發(fā)現(xiàn),QLQR中旋轉(zhuǎn)矩陣和角速度 ω 對(duì)應(yīng)權(quán)重系數(shù)越大,姿態(tài)跟蹤誤差越小.
以上姿態(tài)控制和姿態(tài)跟蹤的兩個(gè)算例說明了Koopman-LQR 控制器的計(jì)算效率較高,且具有不錯(cuò)的控制精度.但由于代價(jià)函數(shù)中參數(shù)的選取對(duì)所提控制器性能的影響比較大,需要根據(jù)需求選取合適的代價(jià)函數(shù).
為了評(píng)估本文所提控制器的控制器性能,將Koopman-LQR 控制器與狀態(tài)依賴Riccati 方程(SDRE)控制器進(jìn)行對(duì)比.SDRE 控制器最初用于解決仿射非線性系統(tǒng)控制問題,作為線性LQR 控制理論在非線性系統(tǒng)上的推廣,其簡化了非線性優(yōu)化問題中HJB (Hamilton-Jacobi-Bellman)方程的求解,避免了求解兩點(diǎn)邊值問題[31].
對(duì)于SO(3)上撓性姿態(tài)動(dòng)力學(xué),通過在每一時(shí)間步迭代求解離散時(shí)間形式的狀態(tài)依賴仿射系統(tǒng)的Riccati 方程得到全狀態(tài)反饋控制率.為了對(duì)比,令SDRE 代價(jià)函數(shù)與Koopman-LQR 代價(jià)函數(shù)形式相同.SDRE 控制器的具體設(shè)計(jì)見附錄.
從同一初始條件出發(fā),采樣間隔為0.2 s,進(jìn)行時(shí)長200 s 的仿真.SDRE 控制器與Koopman-LQR 的前120 s 仿真結(jié)果對(duì)比如圖7 和表3 所示.其中圖7(a)為姿態(tài)跟蹤誤差V1,圖7(b)為式(44)中定義的撓性振動(dòng)指標(biāo)函數(shù)Eq.另外,將跟蹤時(shí)變期望姿態(tài)情況下Koopman-LQR 和SDRE 的性能進(jìn)行對(duì)比,初始條件和期望姿態(tài)取為與4.2 節(jié)中一致,過程中的最大姿態(tài)跟蹤誤差如表3 所示.
表3 Koopman-LQR 與SDRE 結(jié)果對(duì)比Table 3 Comparison between Koopman-LQR and SDRE
圖7 Koopman-LQR 與SDRE 的對(duì)比Fig.7 Comparison between Koopman-LQR and SDRE
圖7(a)和圖7(b)表明在兩種控制器代價(jià)函數(shù)J形式相同的情況下,SDRE 控制器與Koopman-LQR的效果比較接近,SDRE 控制器相對(duì)更快地將撓性航天器控制到期望姿態(tài).表3 結(jié)果表明SDRE 跟蹤時(shí)變姿態(tài)的姿態(tài)誤差要小于Koopman-LQR.而圖7(c)表明Koopman-LQR 對(duì)撓性振動(dòng)的抑制要好于SDRE.在計(jì)算效率方面,SDRE 控制器共用時(shí)約17.9 s,Koopman-LQR 共用時(shí)約1.7 s,這是由于SDRE 控制器在每一采樣間隔內(nèi)都需要求解代數(shù)Riccati 方程,使得其計(jì)算效率相對(duì)更低.表3 中第3 行是跟蹤時(shí)變期望姿態(tài)情況下的最大姿態(tài)跟蹤誤差,結(jié)果表明Koopman-LQR 和SDRE 均能夠使以較高精度使航天器跟蹤期望姿態(tài).由于在線求解最優(yōu)控制的緣故,SDRE 的姿態(tài)跟蹤的精度相對(duì)更高.
綜合來看,不像SDRE 那樣需要對(duì)航天器姿態(tài)動(dòng)力學(xué)進(jìn)行精確建模,Koopman-LQR 只基于數(shù)據(jù)驅(qū)動(dòng)辨識(shí)模型進(jìn)行設(shè)計(jì),雖然相較于SDRE 控制器控制精度有所不如,但是Koopman-LQR 具有計(jì)算效率更高和對(duì)撓性振動(dòng)的抑制效果更好的優(yōu)勢.
本文提出了一組觀測函數(shù)用于數(shù)據(jù)驅(qū)動(dòng)辨識(shí)得到SO(3)上撓性航天器姿態(tài)動(dòng)力學(xué)的Koopman 算子的有限維截?cái)?觀測函數(shù)構(gòu)建為旋轉(zhuǎn)矩陣與姿態(tài)角速度及其部分耦合項(xiàng),還包括撓性構(gòu)件測量點(diǎn)位移與加速度和它們與姿態(tài)角速度和旋轉(zhuǎn)矩陣的耦合項(xiàng).基于SINDY 算法,通過拓展觀測項(xiàng)和引入高階補(bǔ)償項(xiàng)來得到提升空間上的撓性航天器姿態(tài)動(dòng)力學(xué)的廣義線性模型和截?cái)郖oopman 算子.使用SINDY模型,可以以較高精度預(yù)測航天器狀態(tài)軌跡,而其線性部分即截?cái)嗟腒oopman 算子也能在一定時(shí)間內(nèi)以一定精度預(yù)測系統(tǒng)未來狀態(tài).此外,基于Koopman算子理論設(shè)計(jì)了一個(gè)LQR 控制器,能夠?qū)崿F(xiàn)撓性航天器姿態(tài)控制以及振動(dòng)抑制,且該控制器能夠以較高的計(jì)算效率實(shí)現(xiàn)接近最優(yōu)的控制效果.
未來將考慮數(shù)據(jù)存在噪聲和航天器存在外界擾動(dòng)(如重力梯度力矩)的情況,并嘗試將數(shù)據(jù)驅(qū)動(dòng)方法與Koopman 理論用于航天器集群控制等領(lǐng)域.
附錄A SDRE 控制器
考慮如下非線性系統(tǒng)
將上式寫為控制仿射形式,即
其中,f(x)=A(x)x.上述控制仿射形式方程具有如下形式的無線時(shí)間性能函數(shù)
為了確保局部穩(wěn)定,對(duì)于任意x,Q(x) 須為半正定矩陣,R(x)為正定矩陣.上式可以看作線性LQR 方法在非線性控制里的推廣.若P(x) 是滿足如下狀態(tài)依賴代數(shù)Riccati 方程(SDRE)的正定解
則控制輸入可以表示為
在每一采樣間隔內(nèi),將狀態(tài)依賴系統(tǒng)視為線性時(shí)不變系統(tǒng),便可通過求解時(shí)不變代數(shù)Riccati 方程得到控制量.
其中u∈R3為控制力矩,,J為航天器轉(zhuǎn)動(dòng)慣量矩陣.在A(x) 和B(x) 中
其中rij為旋轉(zhuǎn)矩陣R的元素.
其中
其中Jij為航天器轉(zhuǎn)動(dòng)慣量矩陣J的元素.
其中m為模態(tài)坐標(biāo)數(shù)量,K為模態(tài)剛度矩陣,C為模態(tài)阻尼矩陣.
SDRE 控制代價(jià)函數(shù)形式為
為保證SDRE 的代價(jià)函數(shù)與Koopman-LQR 一致,且注意到Koopman-LQR 的控制量是在觀測函數(shù)空間的,因此對(duì)應(yīng)的控制代價(jià)參數(shù)需要進(jìn)行相應(yīng)變換,即令其中為從u到 Ψu的狀態(tài)依賴變換矩陣.在小角速度情況下,有另外,由于在Koopman 建模中采用的是撓性構(gòu)件位移w,而本節(jié)中狀態(tài)依賴控制仿射系統(tǒng)使用模態(tài)坐標(biāo)q,因此也需進(jìn)行變換以保證代價(jià)函數(shù)x形式一致.由于有,則
其中P(x) 是滿足如下方程的正定解