摘 要：針對(duì)現(xiàn)代航天器因攜帶多類型充液貯箱及可機(jī)動(dòng)柔性附件而導(dǎo)致的動(dòng)力學(xué)建模困難，文中基于參數(shù)辨識(shí)理論和有限元方法實(shí)現(xiàn)了該類復(fù)雜航天器的模塊化建模。首先，考慮航天器運(yùn)行時(shí)貯箱內(nèi)液體的小幅晃動(dòng)問題，利用球形貯箱小幅晃動(dòng)的參數(shù)化模型推導(dǎo)出球形充液貯箱小幅晃動(dòng)的動(dòng)力學(xué)方程；其次，根據(jù)薄板小幅振動(dòng)理論，運(yùn)用有限單元法建立柔性附件四邊形板單元?jiǎng)恿W(xué)模型，將附件相對(duì)航天器變化的姿態(tài)角代入坐標(biāo)轉(zhuǎn)換矩陣，構(gòu)造柔性附件作大范圍機(jī)動(dòng)時(shí)的時(shí)變坐標(biāo)轉(zhuǎn)換矩陣。再次，基于凱恩方程推導(dǎo)航天器整體系統(tǒng)的動(dòng)力學(xué)狀態(tài)方程，利用MATLAB軟件編制出相應(yīng)的模塊化建模程序。最后，通過數(shù)值算例分析，研究典型構(gòu)型航天器中柔性附件以不同方式機(jī)動(dòng)的系統(tǒng)整體耦合動(dòng)力學(xué)性能，驗(yàn)證了該建模方法的通用性、適用性和準(zhǔn)確性。

關(guān)鍵詞：充液柔性航天器；參數(shù)化模型；有限單元法；凱恩方程；模塊化建模

中圖分類號(hào)：V412.4；V448.2 DOI：10.16375/j.cnki.cn45-1395/t.2024.01.001

0 引言

現(xiàn)代航天器不斷向多任務(wù)、模塊化和長壽命方向發(fā)展，攜帶諸多設(shè)備和大尺寸柔性附件的航天器是復(fù)雜的剛-液-柔耦合系統(tǒng)，準(zhǔn)確建立系統(tǒng)動(dòng)力學(xué)模型不僅能描述航天器系統(tǒng)的耦合動(dòng)力學(xué)特性，也能為系統(tǒng)控制律的設(shè)計(jì)提供理論依據(jù)和前期研究基礎(chǔ)。

航天器的典型構(gòu)型通常為中心剛體攜帶多個(gè)充液貯箱，外側(cè)連接太陽能帆板。太陽能電池板因具有較大長寬尺度和極小厚度的幾何尺寸特征，根據(jù)具體長寬比可簡化為Euler-Bernoulli梁模型[1-3]或Kirchhoff-Love薄板模型[4]。附件相對(duì)主剛體運(yùn)動(dòng)的耦合效應(yīng)[3，6-7]通過定義相對(duì)運(yùn)動(dòng)廣義坐標(biāo)[3]或構(gòu)造時(shí)變坐標(biāo)轉(zhuǎn)換矩陣[6]等途徑對(duì)航天器動(dòng)力學(xué)模型進(jìn)行改良?，F(xiàn)代航天器所攜帶的液體燃料占據(jù)總質(zhì)量的較大比重[8]，等效力學(xué)模型被廣泛應(yīng)用于液體的晃動(dòng)效應(yīng)的簡化和替代，例如彈簧-質(zhì)量-阻尼模型[9-11]、單擺模型[12-13]和復(fù)合擺模型[14-15]。區(qū)別于傳統(tǒng)等效模型參數(shù)與機(jī)械結(jié)構(gòu)參數(shù)相對(duì)應(yīng)的限制，液體晃動(dòng)參數(shù)化模型[16-18]為實(shí)驗(yàn)或仿真數(shù)據(jù)逆向辨識(shí)所得到，適用性更強(qiáng)，而且對(duì)仿真數(shù)據(jù)的辨識(shí)精度可達(dá)92%以上[19]，具有較好的研究價(jià)值和應(yīng)用前景。對(duì)于非完整約束系統(tǒng)動(dòng)力學(xué)方程的建立，Kane等[20-21]在推導(dǎo)系統(tǒng)動(dòng)力學(xué)方程時(shí)，提出一種基于偏速度和偏角速度為系統(tǒng)自由度的建模方法，又稱為凱恩方程。相較于傳統(tǒng)拉格朗日法，凱恩方程步驟程序化、物理意義清楚，廣泛應(yīng)用于航天器系統(tǒng)動(dòng)力學(xué)中[22]。

現(xiàn)階段，柔性附件充液航天器的動(dòng)力學(xué)研究存在以下局限性：1）等效力學(xué)模型適用范圍存在一定局限，而參數(shù)化模型可以實(shí)時(shí)辨識(shí)并修改參數(shù)模型，描述更準(zhǔn)確。2）拉格朗日方程推導(dǎo)多自由度耦合系統(tǒng)動(dòng)力學(xué)方程較為困難且繁瑣。因此，本文基于文獻(xiàn)[19]的結(jié)論，構(gòu)造球形充液貯箱小幅振動(dòng)參數(shù)化模型；運(yùn)用有限單元法推導(dǎo)柔性附件動(dòng)力學(xué)方程，根據(jù)附件相對(duì)航天器的時(shí)變姿態(tài)角計(jì)算時(shí)變坐標(biāo)轉(zhuǎn)換矩陣；最終由凱恩方程模塊化地建立航天器系統(tǒng)的耦合動(dòng)力學(xué)方程，研究系統(tǒng)狀態(tài)變量在柔性部件做不同機(jī)動(dòng)和不同數(shù)量下的動(dòng)力學(xué)響應(yīng)。

1 航天器系統(tǒng)剛-液-柔耦合動(dòng)力學(xué)

方程推導(dǎo)

航天器系統(tǒng)的典型構(gòu)型通常簡化為中心剛體帶多個(gè)充液貯箱，并外接柔性附件，如圖1所示。定義以下參考系：慣性參考系N，[n1、n2、n3]為慣性參考系基矢量；與航天器本體固連的本體參考系B，其原點(diǎn)位于主剛體質(zhì)心O，[b1、b2、b3]為本體參考系的基矢量，初始時(shí)刻本體系基矢量與慣性系對(duì)應(yīng)基矢量方向一致；附件局部坐標(biāo)系S，其原點(diǎn)位于附件在航天器的安裝點(diǎn)[OS]，[s11、s12、s13]為附件局部參考系的基矢量；附件有限元模型中第i個(gè)四邊形板單元的局部坐標(biāo)系[Ki]，其原點(diǎn)位于四邊形板單元幾何中心，[ki1][、][ki2][、][ki3]為單元參考系基矢量，單元參考系基矢量與附件局部參考系S對(duì)應(yīng)基矢量方向一致。

1.1 航天器主體姿態(tài)動(dòng)力學(xué)建模

航天器的軌道機(jī)動(dòng)由質(zhì)心O相對(duì)于慣性系N的平動(dòng)進(jìn)行描述，姿態(tài)調(diào)整由本體系B相對(duì)于慣性系N的轉(zhuǎn)動(dòng)進(jìn)行描述。定義[rO=NrTON=rO1，rO2，rO3T]

[n1，n2，n3]為主剛體質(zhì)心O相對(duì)慣性系的矢徑，其中列矩陣[NrO=rO1，rO2，rO3]為矢徑[rO]投影到慣性坐標(biāo)系N在其基矢量方向上的分量陣，左上角標(biāo)N表示矢量向慣性參考系投影；[N=n1，n2，n3]為慣性系N的基矢量組合。同理，矢徑[rO]可由本體系B基矢量表示為：[rO=BrTOB]，其中[BrTO、B]的定義與上述N同理。航天器在軌運(yùn)行時(shí)其相對(duì)慣性系的矢徑快速變化，在研究小擾動(dòng)對(duì)航天器的影響時(shí)，將本體參考系下計(jì)算矢徑響應(yīng)結(jié)果即為相對(duì)穩(wěn)定運(yùn)行的影響。對(duì)主剛體質(zhì)心O相對(duì)慣性系N的絕對(duì)速度[vO]由矢徑[rO]對(duì)時(shí)間求導(dǎo)得到，并投影到本體坐標(biāo)系B中，則有

[vO=drOdt=BvTOB+Bω×0BrOTB，] （1）

式中：[Bω×0]為角速度在本體系投影的分量矩陣[Bω0]的反對(duì)稱矩陣，其中[ω0]為本體參考系B相對(duì)慣性參考系N旋轉(zhuǎn)的角速度。主剛體的絕對(duì)加速度[aO]可通過對(duì)速度求導(dǎo)后得到，即

[aO=BaTOB+Bω×0BvOTB+Bα×0BrOTB+]

[Bω×0Bω×0BrOTB，] （2）

式中：[Bα×0]為角速度在本體系投影的分量陣[Bα0]的反對(duì)稱矩陣，其中[α0]為本體參考系B相對(duì)慣性參考系N旋轉(zhuǎn)的角加速度。

1.2 柔性附件動(dòng)力學(xué)建模

考慮三軸姿態(tài)穩(wěn)定航天器在軌運(yùn)行時(shí)擾動(dòng)較小，帆板變形基于小變形理論描述，則無需考慮材料大幅變形時(shí)的動(dòng)力剛化效應(yīng)。構(gòu)造太陽能帆板四邊形薄板有限單元模型，將太陽能帆板劃分為多個(gè)較小四邊形板單元并對(duì)節(jié)點(diǎn)編號(hào)。根據(jù)有限元法，板單元各節(jié)點(diǎn)具有3個(gè)自由度，即單元變形時(shí)節(jié)點(diǎn)產(chǎn)生3個(gè)方向的位移：相對(duì)[ki1、ki2]的角位移[θ1、θ2]以及垂直于板平面的撓度w。撓度僅與該點(diǎn)在單元局部參考系的平面坐標(biāo)相關(guān)，即[w（ki1，ki2）]；根據(jù)Kirchhoff-Love板理論，薄板厚度方向各點(diǎn)具有相同撓度變形w。

考慮第[sl]個(gè)柔性附件上一點(diǎn)P，該點(diǎn)在慣性參考系的矢徑可由矢量疊加法得到，即

[rP=rO+rsl+rslP+uP=BrTOB+BrTslB+]

[CBSSrTslP+SuTPB，] （3）

式中：[rsl=BrTslB]為第[sl]個(gè)柔性附件安裝點(diǎn)[OSl]相對(duì)主剛體質(zhì)心的矢徑，并在本體參考系B下表示，[Brsl]為[rsl]在本體參考系B投影的分量矩陣；[rslP=CBSSrTslPB]為P點(diǎn)相對(duì)附件安裝點(diǎn)的矢徑，[SrslP]為[rslP]在附件局部參考系S投影的分量矩陣，[CBS]為附件局部參考系與本體系的坐標(biāo)轉(zhuǎn)換矩陣，轉(zhuǎn)換關(guān)系為[S=CBSB]；[uP=CBSSuTPB]為P點(diǎn)發(fā)生變形時(shí)相對(duì)初始位置的位移矢量，列矩陣[SuP]為[uP]在附件局部參考系S投影的分量矩陣。根據(jù)有限元法，第i個(gè)四邊形板單元上一點(diǎn)P垂直于板平面的撓度變形w由第i個(gè)單元的節(jié)點(diǎn)位移[δei]和該點(diǎn)在單元局部參考系坐標(biāo)[ki1，ki2]共同決定，即[SuP=0，0，wki1， ki2T=0， 0， Nki1， ki2?δeiT]，其中[Nx， y]為4節(jié)點(diǎn)板單元形函數(shù)。同理，質(zhì)點(diǎn)P在慣性參考系的矢徑[rP]對(duì)時(shí)間求導(dǎo)得到絕對(duì)速度[vP]，并將其投影到本體坐標(biāo)系后，有

[vP=drPdt=BvToB+Bω×0BroTB+Bω×0BrslTB+]

[CBSSω×0+Sω×1SrP+SuPTB+CBSSuTPB，] （4）

式中：[Sω×0和Sω×1]分別為主剛體角速度和附件相對(duì)角速度在局部參考系S投影的分量矩陣[Sω0和Sω1]的反對(duì)稱矩陣，[ω1]為柔性附件相對(duì)主剛體的旋轉(zhuǎn)角速度。

柔性附件變形后質(zhì)點(diǎn)相對(duì)慣性系的絕對(duì)加速度[aP]可描述為

[aP=dvPdt=BaToB+Bω×0BvoTB+Bα×0BrOTB+][Bω×0Bω×0BrOTB+Bα×0BrslTB+Bω×0Bω×0BrslTB+][Sα×0+Sα×1SrslP+SuPTS+2Sω×0+Sω×1SuPTS+][Sω×0+Sω×1Sω×0+Sω×1SrslP+SuPTS+SuTPS，] （5）

式中：[α1]為柔性附件相對(duì)主剛體的旋轉(zhuǎn)角速度；[Sα×0和Sα×1]分別為主剛體角加速度和附件相對(duì)角加速度在局部參考系S投影的分量矩陣[Sα0和Sα1]的反對(duì)稱矩陣。

有限單元法中將作用于板單元的分布載荷由虛功原理等效為單元節(jié)點(diǎn)載荷，由式（5）質(zhì)點(diǎn)加速度[aP]和形函數(shù)[Nx，y]在單元內(nèi)積分得到單元節(jié)點(diǎn)載荷矩陣：[Fe=NTρhaPdS]。按節(jié)點(diǎn)編號(hào)將單元節(jié)點(diǎn)載荷矩陣組裝得到附件總體載荷矩陣[Fsl]，總體質(zhì)量矩陣M和總體剛度矩陣K[23]。定義[δsl]為附件節(jié)點(diǎn)位移列矩陣，聯(lián)立得到柔性附件四邊形板單元振動(dòng)變形的動(dòng)力學(xué)方程為

[Mδsl+Cδsl+Kδsl=Fsl.] （6）

1.3 液體晃動(dòng)等效力學(xué)模型

航天器液體燃料貯箱安裝位置如圖1所示。由于三軸姿態(tài)穩(wěn)定航天器運(yùn)行時(shí)產(chǎn)生的擾動(dòng)較小，航天器控制系統(tǒng)消耗的燃料忽略不計(jì)，即假設(shè)液體燃料質(zhì)量不變。由文獻(xiàn)[19]的結(jié)論，基于勢(shì)流理論和傅里葉-貝塞爾函數(shù)展開系數(shù)條件，在忽略高階項(xiàng)時(shí)球形貯箱液體晃動(dòng)狀態(tài)方程如式（7）所示。顯然該狀態(tài)方程的結(jié)構(gòu)與傳統(tǒng)物理等效力學(xué)模型動(dòng)力學(xué)方程相似，例如單擺等效力學(xué)模型等。

[ddtb1（t）a1（t）=01-ω2-2ωζb1（t）a1（t）+0Aax（t），] （7）

式中：[b1（t）， a1（t）T]為系統(tǒng)狀態(tài)向量；[ζ]為一階模態(tài)的阻尼比；[ω]為系統(tǒng)固有頻率；[A]為輸入矩陣系數(shù)；[ax（t）]為橫向晃動(dòng)加速度。

在研究液體晃動(dòng)效應(yīng)時(shí)大多不需要推導(dǎo)流體域各處情況，通常將晃動(dòng)液體視為整體參與晃動(dòng)，因此根據(jù)狀態(tài)方程，設(shè)參與晃動(dòng)的液體等效質(zhì)心m在外激勵(lì)下做小幅振動(dòng)時(shí)，質(zhì)心相對(duì)貯箱球心以固定長度晃動(dòng)；未參與晃動(dòng)的靜止液體的等效質(zhì)心設(shè)為[m1]。由矢量疊加法，質(zhì)心m在慣性系中的矢徑如下

[rm=rO+rsP+rl+um=BrTOB+BrTsPB+BrTlB+BuTmB，]

（8）

式中：[rsP=BrTsPB]為球形貯箱球心相對(duì)主剛體質(zhì)心的矢徑，列矩陣[BrsP]為矢徑[rsP]在本體參考系B投影的分量矩陣；[rl=BrTlB]為充液貯箱靜置時(shí)質(zhì)心m相對(duì)球形貯箱球心的矢徑，列矩陣[Brl]為矢徑[rl]在本體參考系B投影的分量矩陣；[um=BuTmB]為質(zhì)心m晃動(dòng)時(shí)偏離初始位置的位移矢量，列矩陣[Bum]是[um]在本體參考系B投影的分量矩陣。

在文獻(xiàn)[19]中，參數(shù)化模型對(duì)球形充液貯箱受橫向外激勵(lì)晃動(dòng)的輸出響應(yīng)辨識(shí)效果較好，而且小幅晃動(dòng)下質(zhì)心m擺動(dòng)的角位移可近似等效為沿[b1、b2]方向的線位移，即[Bum=[um1，um2，0]T]。同理，晃動(dòng)液體等效質(zhì)點(diǎn)m在慣性系下的速度和加速度通過對(duì)時(shí)間求導(dǎo)后，得到

[vm=drmdt=BvTOB+Bω×0BrsP+Brl+BumTB+BuTmB，] （9）

[am=dvmdt=BaTOB+Bω×0BvOTB+][Bα×0BrsP+Brl+BumTB+Bω×0Bω×0BrsP+Brl+BumTB+][Bω×0BumTB+BuTmB.] （10）

最終，球形充液貯箱小幅晃動(dòng)（包括待辨識(shí)參數(shù)）的等效動(dòng)力學(xué)狀態(tài)方程可整理為

[m00mum1um2+c00cum1um2=-mBam（1）-mBam（2）]， （11）

式中：c為液體晃動(dòng)等效阻尼系數(shù)；[Bam（1）、Bam（2）]為質(zhì)心m的加速度本體參考系的投影在[b1、b2]方向上的分量。

2 基于凱恩方程的航天器耦合系統(tǒng)

動(dòng)力學(xué)方程推導(dǎo)

凱恩將偏速度和偏角速度作為系統(tǒng)自由度，并聯(lián)立廣義主動(dòng)力和廣義慣性力在偏速度和偏角速度方向的投影，最終得到系統(tǒng)動(dòng)力學(xué)方程。

2.1 偏速度和偏角速度

偏速度和偏角速度統(tǒng)一用[v（r）i]表示，其中右上角標(biāo)（r）表示第r個(gè)偏速度。對(duì)于航天器主剛體相對(duì)慣性系的平動(dòng)和轉(zhuǎn)動(dòng)，定義偏速度和偏角速度：

[vji=BvO（j）?bj， " " " " "j=1，2，3，] （12）

[v3+ji=Bω0（j）?bj]， [j=1，2，3]. （13）

設(shè)航天器共安裝[sl]個(gè)附件，每個(gè)柔性附件劃分n個(gè)節(jié)點(diǎn)，四邊形板單元每個(gè)節(jié)點(diǎn)含3個(gè)自由度。對(duì)于第s個(gè)柔性附件的第j個(gè)節(jié)點(diǎn)的振動(dòng)變形，進(jìn)一步定義偏速度和偏角速度：

[v6+3n（s-1）+3（j-1）+1i=θ1j?b1，]

[s=1，2，…，sl;j=1，2，…，n]， （14）

[v6+3n（s-1）+3（j-1）+2i=θ2j?b2，]

[s=1，2，…，sl;j=1，2，…，n]， （15）

[v6+3n（s-1）+3（j-1）+3i=wj?b3，]

[s=1，2，…，sl;j=1，2，…，n.] （16）

另設(shè)航天器總共攜帶q個(gè)部分充液貯箱。對(duì)于第k個(gè)充液貯箱的參與晃動(dòng)液體質(zhì)心m相對(duì)主剛體的運(yùn)動(dòng)，可定義如下偏速度和偏角速度：

[v6+3nsl+2（k-1）+1i=um1?b1， k=1，2，…，q]， （17）

[v6+3nsl+2（k-1）+2i=um2?b2][， k=1，2，…，q]. （18）

2.2 廣義慣性力和廣義主動(dòng)力

對(duì)于第j個(gè)偏速度方向，航天器耦合系統(tǒng)所受的廣義慣性力[Fj]的計(jì)算公式為

[Fj=-mBBaO?v（j）i-IB?NαB0+NωB0×IB?NωB0?v（j）i-]

[s=1slSNaP?v（j）idmP-k=1qmkBamk?v（j）i，] （19）

式中：[mB]為主剛體的質(zhì)量；[IB]為主剛體關(guān)于本體系的二階慣性張量；[mP]為附件上任意質(zhì)點(diǎn)P的質(zhì)量；[mk]為第k個(gè)貯箱內(nèi)液體晃動(dòng)等效質(zhì)心的質(zhì)量；[amk]為第k個(gè)充液貯箱晃動(dòng)液體等效質(zhì)點(diǎn)的晃動(dòng)加速度。

航天器系統(tǒng)所受的廣義作用力包括液體晃動(dòng)的反作用力，以及柔性附件振動(dòng)變形的彈性作用力。對(duì)于第j個(gè)偏速度方向，航天器耦合系統(tǒng)所受的廣義作用力[Fj*]的計(jì)算公式為

[Fj*=FsP?v（j）i+Fsl?v（j）i]， （20）

式中：[FsP=k=1qmk（um1b1+um2b2）+c（um1b1+um2b2）]為液體晃動(dòng)等效質(zhì)心的作用力；[Fsl]為柔性附件振動(dòng)變形對(duì)支座的反作用力，計(jì)算公式為[Fsl=s=1slCSBKj，：sδsls]，其中[CSB]為附件局部參考系S到本體系B的坐標(biāo)轉(zhuǎn)換矩陣，（j，：）表示第j行所有列的全部元素。

2.3 基于凱恩方程的系統(tǒng)動(dòng)力學(xué)方程建模

將相同偏速度方向下的廣義主動(dòng)力和廣義慣性力分量聯(lián)立得到系統(tǒng)動(dòng)力學(xué)方程。在第j個(gè)偏速度方向下的凱恩方程可寫為

[Fj=F*j]. （21）

將式（12）—式（18）中求得的偏速度代入式（19）和式（20），可建立描述系統(tǒng)運(yùn)動(dòng)的動(dòng)力學(xué)方程組。將航天器系統(tǒng)所有自由度的響應(yīng)按偏速度的順序組合成位移列矩陣u，最后將得到的動(dòng)力學(xué)方程組整理成如下形式，

[Mu+Cu+Ku=F]， （22）

式中：M、C、K分別為系統(tǒng)質(zhì)量矩陣、阻尼矩陣和剛度矩陣。質(zhì)量矩陣M由主剛體、柔性附件和液體晃動(dòng)參數(shù)模型的質(zhì)量矩陣對(duì)角排列構(gòu)成，且包含了廣義主動(dòng)力中系數(shù)為廣義速率導(dǎo)數(shù)的項(xiàng)；同理，阻尼矩陣C和剛度矩陣K的排列構(gòu)成類似，且包含了廣義主動(dòng)力中系數(shù)為廣義速率的項(xiàng)和系數(shù)為自由度位移的項(xiàng)?？紤]MATLAB微分方程ODE求解器僅適用于一階方程，需對(duì)系統(tǒng)方程改寫為等效一階方程組。系統(tǒng)廣義位移和廣義速度記為系統(tǒng)的狀態(tài)向量[u，uT]，將上述動(dòng)力學(xué)方程降階為

[IOOMddtuu=OI-K-Cuu+OF]. （23）

3 數(shù)值仿真

對(duì)上述建立的含柔性附件充液柔性航天器的系統(tǒng)動(dòng)力學(xué)方程，編寫相應(yīng)MATLAB程序。在算例1中，航天器各部件參數(shù)和構(gòu)型參照文獻(xiàn)[24]，航天器所受外激勵(lì)如下所示

[M02（t）10， 0lt;tlt;10， " " " " " " " " " " " " " -10， 20lt;tlt;30， " " " " " " " " " " 0， t≤0， 10≤t≤20， t≥30. ] （24）

編程計(jì)算航天器系統(tǒng)在外激勵(lì)作用下的耦合動(dòng)力學(xué)響應(yīng)如圖2—圖5所示。其中，航天器速度與角速度響應(yīng)如圖2與圖3所示，柔性附件自由邊中點(diǎn)的位移響應(yīng)如圖4所示，液體晃動(dòng)等效質(zhì)心位移響應(yīng)如圖5所示。

編程求解結(jié)果與文獻(xiàn)[24]相比較，由于忽略了帆板小幅振動(dòng)時(shí)的阻尼效應(yīng)，在外激勵(lì)作用結(jié)束后帆板仍進(jìn)行小幅的振動(dòng)變形，而航天器響應(yīng)與附件位移響應(yīng)曲線的趨勢(shì)和數(shù)值大小與參考文獻(xiàn)結(jié)果相接近，初步驗(yàn)證了推導(dǎo)和程序的準(zhǔn)確性。

太陽能帆板屬于大撓度附件，當(dāng)帆板相對(duì)航天器本體運(yùn)動(dòng)時(shí)對(duì)其姿態(tài)角產(chǎn)生較大影響。在算例2、3中，考慮柔性附件進(jìn)行機(jī)動(dòng)時(shí)，計(jì)算航天器動(dòng)力學(xué)響應(yīng)。算例2的航天器構(gòu)型如圖1所示，初始時(shí)刻各部件處于靜止?fàn)顟B(tài)。航天器系統(tǒng)包含4個(gè)球形充液貯箱，貯箱參數(shù)參考文獻(xiàn)[25]，其球心相對(duì)本體坐標(biāo)系原點(diǎn)的坐標(biāo)如表1所示。

航天器本體安裝一個(gè)機(jī)動(dòng)太陽能帆板，該柔性附件安裝點(diǎn)為[rsl=[0，0.15，0]T]，安裝角為[θ01=[0，0，π 2]T]。柔性附件的尺寸和材料參數(shù)如表2所示。

機(jī)動(dòng)附件的驅(qū)動(dòng)類型采用Bang-Bang驅(qū)動(dòng)[26]，柔性附件繞局部參考系的[s12]軸相對(duì)航天器旋轉(zhuǎn)。其相應(yīng)的角加速度的變化如下所示：

[α（t）=0.021， 0lt;tlt;5， " " " " " "-0.021， 5lt;tlt;10， " " " 0， t≤0， t=5， t≥10. " ] （25）

編程計(jì)算柔性附件機(jī)動(dòng)時(shí)航天器動(dòng)力學(xué)響應(yīng)如圖6—圖7所示，充液貯箱晃動(dòng)質(zhì)心位移響應(yīng)如圖8—圖9所示，柔性附件自由邊中點(diǎn)位移響應(yīng)如圖10所示。

附件繞坐標(biāo)軸[s12]進(jìn)行旋轉(zhuǎn)機(jī)動(dòng)時(shí)，振動(dòng)變形的反作用力矩為本體參考系的[b1]方向，因此結(jié)合反作用力對(duì)質(zhì)心的力矩，在附件機(jī)動(dòng)時(shí)角速度響應(yīng)[ω01]較大。附件的機(jī)動(dòng)形式和驅(qū)動(dòng)形式相關(guān)，振動(dòng)趨勢(shì)與文獻(xiàn)[26]中算例相接近。由于部分充液貯箱等效晃動(dòng)質(zhì)心一般并不位于球心，因此主剛體運(yùn)動(dòng)產(chǎn)生的牽連載荷各不相同，產(chǎn)生復(fù)雜的晃動(dòng)現(xiàn)象。

算例3：柔性附件繞局部參考系的[s11]軸相對(duì)航天器旋轉(zhuǎn)，驅(qū)動(dòng)和其余參數(shù)與算例2對(duì)應(yīng)參數(shù)相同。編程計(jì)算柔性附件機(jī)動(dòng)時(shí)航天器動(dòng)力學(xué)響應(yīng)如圖11—圖12所示，充液貯箱晃動(dòng)質(zhì)心位移響應(yīng)如圖13—圖14所示，柔性附件自由邊中點(diǎn)位移響應(yīng)如圖15所示。

附件繞坐標(biāo)軸[s11]進(jìn)行旋轉(zhuǎn)機(jī)動(dòng)時(shí)，振動(dòng)變形的反作用力矩為本體參考系的[b2]方向，相較于前文算例，附件機(jī)動(dòng)時(shí)產(chǎn)生一定的角速度響應(yīng)[ω02]。由于附件相對(duì)[s11]的慣性矩更小，因此在相同驅(qū)動(dòng)下各部件的動(dòng)力學(xué)響應(yīng)也更小。

4 結(jié)論

本文將三軸姿態(tài)穩(wěn)定充液柔性航天器為研究對(duì)象，考慮航天器運(yùn)行過程中的小幅擾動(dòng)，基于球形貯箱液體晃動(dòng)參數(shù)化模型推導(dǎo)充液貯箱晃動(dòng)的動(dòng)力學(xué)方程，運(yùn)用有限元法建立板類柔性附件的四邊形單元模型，基于Kirchhoff-Love板理論描述單元位移。由凱恩方程建立航天器剛-液-柔耦合系統(tǒng)動(dòng)力學(xué)模型。仿真對(duì)比分析充液柔性航天器在單個(gè)板類附件沿不同轉(zhuǎn)軸機(jī)動(dòng)下的動(dòng)力學(xué)響應(yīng)，結(jié)果顯示航天器姿態(tài)角和軌道狀態(tài)參數(shù)、液體晃動(dòng)效應(yīng)、附件機(jī)動(dòng)與柔性變形的相互耦合作用，驗(yàn)證了該方法建立的耦合系統(tǒng)動(dòng)力學(xué)模型的準(zhǔn)確性。體現(xiàn)了該方法對(duì)不同附件驅(qū)動(dòng)類型、不同附件數(shù)量航天器的模塊化建模特性。

參考文獻(xiàn)

[1] GHORBANI H，VATANKHAH R，F(xiàn)ARID M. General planar motion modeling and control of a smart rigid-flexible satellite considering large deflections[J]. Nonli-near Dynamics，2022，108：911-939.

[2] GHORBANI H，VATANKHAH R，F(xiàn)ARID M. Adaptive nonsingular fast terminal sliding mode controller design for a smart flexible satellite in general planar motion[J]. Aerospace Science and Technology，2021，119. Doi：10.1016/j.ast.2021.107100.

[3] MANCINI M，CAPELLO E. Adaptive sliding mode-based control system for flexible spacecraft[J]. Proceedings of the American Control Conference，2021（2021-May）：2968-2973.

[4] DA FONSECA I M，RADE D A，GOES L C S，et al. Attitude and vibration control of a satellite containing flexible solar arrays by using reaction wheels，and piezoelectric transducers as sensors and actuators[J]. Acta Astronautica，2017（139）：357-366.

[5] HE G Q，CAO D Q，CAO Y T，et al. Dynamic modeling and orbit maneuvering response analysis for a three-axis attitude stabilized large scale flexible spacecraft installed with hinged solar arrays[J]. Mechanical Systems and Signal Processing，2022，162：108083.

[6] WANG J，LI D X，JIANG J P. First order coupled dynamic model of flexible space structures with time-varying configurations[J]. Acta Astronautica，2017，132：117-123.

[7] ZHANG G W，WANG F，CHEN J，et al. Fixed-time sliding mode attitude control of a flexible spacecraft with rotating appendages connected by magnetic bearing[J]. Mathematical Biosciences and Engineering，2022，19（3）：2286-2309.

[8] BURK T A. Cassini orbit trim maneuvers at saturn-overview of attitude control flight operations[M]. AIAA Guidance，Navigation，and Control Conference，2011.

[9] ALLARD C，RAMOS M D，SCHAUB H，et al. Spacecraft dynamics integrating hinged solar panels and lumped-mass fuel slosh model[C]//AIAA/AAS Astrodynamics Specialist Conference 2011，2011.

[10] CHIBA M，MAGATA H. Coupled pitching dynamics of flexible space structures with on-board liquid sloshing[J]. Acta Astronautica，2021，181：151-166.

[11] COLAGROSSI A，LAVAGNA M. Integrated vibration suppression attitude control for flexible spacecrafts with internal liquid sloshing[J]. Multibody System Dynamics，2020，51（2）：123-157.

[12] SUN Y，ZHOU D，WANG J D，et al. Liquid sloshing in a cylindrical tank with multiple baffles under horizontal and pitching motions[J]. International Journal of Applied Mechanics，2020，12（7）：190-197.

[13] GASBARRI P，SABATINI M，PISCULLI A. Dynamic modelling and stability parametric analysis of a flexible spacecraft with fuel slosh[J]. Acta Astronautica，2016，127：141-159.

[14] LIU F，YUE B Z，TANG Y，et al. 3DOF-rigid-pendulum analogy for nonlinear liquid slosh in spherical propellant tanks[J]. Journal of Sound and Vibration，2019，460：114907.

[15] LU Y，YUE B Z，MA B L. Improved moving pulsating ball equivalent model for large-amplitude liquid slosh[J]. AIAA Journal，2022，60（8）：5004-5011.

[16] 周鳳霞.帶貯液罐類重型車輛操縱動(dòng)力學(xué)與穩(wěn)定性研究[D].柳州：廣西科技大學(xué)，2020.

[17] 黃中烈，吳文軍，周鳳霞，等.罐車動(dòng)力學(xué)建模與剛-液耦合特性研究[J].現(xiàn)代制造工程，2020（8）：20-26.

[18] 李超，吳文軍，周鳳霞，等.重型罐車內(nèi)液體晃動(dòng)參數(shù)化數(shù)學(xué)模型辨識(shí)[J].廣西科技大學(xué)學(xué)報(bào)，2021，32（2）：65-71.

[19] 王佐，吳文軍，姜丙超，等.球形貯箱內(nèi)液體橫向晃動(dòng)問題的系統(tǒng)辨識(shí)與參數(shù)化等效建模[J].廣西科技大學(xué)學(xué)報(bào)，2022，33（3）：1-7.

[20] KANE T R，LEVINSON D A. Formulation of equations of motion for complex spacecraft[J]. Journal of Guidance，Control，and Dynamics，1980，3（2）：99-112.

[21] KANE T R，RYAN R R，BANERJEE A K. Dynamics of a cantilever beam attached to a moving base[J]. Journal of Guidance，Control，and Dynamics，1987，10（2）：139-151.

[22] LIU F，YUE B Z，BANERJEE A K，et al. Large motion dynamics of in-orbit flexible spacecraft with large-amplitude propellant slosh[J]. Journal of Guidance，Control，and Dynamics，2020，43（3）：438-450.

[23] 韓清凱，孫偉，王伯平，等.機(jī)械結(jié)構(gòu)有限單元法基礎(chǔ)[M].北京：科學(xué)出版社，2013.

[24] 閆玉龍，申云峰，岳寶增.含板類柔性附件的充液航天器動(dòng)力學(xué)研究[J].宇航學(xué)報(bào)，2020，41（5）：531-540.

[25] 苗楠，李俊峰，王天舒.橫向激勵(lì)下液體大幅晃動(dòng)建模分析[J].宇航學(xué)報(bào)，2016，37（3）：268-274.

[26] 吳文軍.帶多充液貯箱及機(jī)動(dòng)柔性附件航天器耦合動(dòng)力學(xué)研究[D].北京：北京理工大學(xué)，2015.

Modular modeling and coupled dynamical study for liquid-filling spacecraft with multiple maneuverable flexible appendages

LUO Houlin， WU Wenjun*， WANG Zuo

（School of Mechanical and Automotive Engineering， Guangxi University of Science and Technology，

Liuzhou 545616， China）

Abstract： This paper focuses on the problem of dynamic modeling for modern spacecraft with multi-type liquid-filling tanks and maneuverable flexible appendages. Based on parameter identification theory and finite element method， the dynamic modeling for complex spacecraft is deduced. Firstly， in view of small disturbance during spacecraft operation， the dynamic equations of small-amplitude liquid sloshing in spherical tank are derived by using the parametric model for liquid sloshing. Secondly， based on the theory for small deformation of plate， the dynamic quadrilateral plate element model of flexible appendage is established using finite element method. The time-varying coordinate transformation matrix of flexible appendage relative to spacecraft is constructed by substituting the attitude angle of the appendage into the transformation matrix. Then， the dynamic state equations of the spacecraft system are derived based on Kane's equation， and the corresponding modular modeling program is compiled by MATLAB software. Finally， through numerical example analysis， the coupling dynamic performance of the typical configuration spacecraft with flexible appendages maneuvering in different ways is studied， and the universality， availability and accuracy of the modeling method are confirmed.

Keywords： liquid-filling flexible spacecraft; parametric model; finite element method; Kane's equation; modular modeling

（責(zé)任編輯：黎 婭、于艷霞）