周偉靜,洪延姬,葉繼飛,李南雷,常浩
航天工程大學 激光推進及其應用國家重點實驗室,北京 101416
微納衛星具有重量輕、性能好、開發成本低、周期短、功耗低等特點,目前已經成為中型、大型等傳統衛星的有效補充,并帶來了全新的應用模式和應用理念。微納衛星的全面發展推動了微推力器的快速發展,冷氣微推力器[1]、脈沖等離子體微推力器[1]、場發射電推進[1]、膠體微推力器[1]、MEMS(Micro-Electro-Mechanical Systems)陣列[2]、激光燒蝕微推力器[3]、微波微推力器[4-5]等新型微推力器應運而生。
隨著這些微推力器的設計、研制和應用,微推力測量技術亦備受關注[6-7]。微推力測量一般分為穩態推力測量[8-11]、脈沖平均推力[9-10]和脈沖沖量測量[9-10,12-14],穩態推力范圍一般在幾百微牛至幾十毫牛,脈沖沖量范圍一般在幾十微牛秒至幾百微牛秒。穩態推力測量和脈沖沖量測量均采用直接測量方法,即將推力器與推力測量裝置固連,利用推力的反沖作用,將推力轉化為推力測量裝置的力學行為,如振動幅值或轉動位移,測量振幅或位移等信息可獲得推力或沖量大小。目前,水平旋轉的扭擺是最適宜測量微小推力或沖量的測量裝置。
對于脈沖式微推力器,脈沖沖量測量的實現相對簡單,只要力的作用時間小于欠阻尼扭擺測量系統振動周期的1/4[15],可以看做脈沖力瞬間作用于扭擺,根據扭擺響應的最大幅值即可計算脈沖沖量。對于穩態微推力器或高頻多脈沖微推力器(輸出等同于穩態推力),只要力或多脈沖的作用時間大于欠阻尼扭擺系統允許誤差范圍的調整時間,并且維持穩態一定時間,就可以根據扭擺響應的穩態幅值計算穩態推力大小[16]。
在機理研究階段和設計研制過程中,微推力器的推進性能需要測試跟蹤,以形成研究閉環,不斷改進和提升推進性能。目前,微推力器的研制多處于原理樣機或初樣機階段。脈沖式微推力器的單脈沖沖量測量容易實現。但是,由于攜帶推進劑受限(如激光燒蝕微推力器)、頻率控制受限(陰極微弧微推力器)、功率受限(如微波推力器)等原因,單次實驗條件下,穩態或多脈沖微推力器開機時長較短,在扭擺響應達到調整時間之前就需關機[17],無法根據穩態響應還原穩態推力。
通過增大阻尼比和固有頻率,可以減小扭擺系統的調整時間,使扭擺響應能盡快進入穩態。但是,一方面,為了能夠高精度地標定測量系統的系統參數,阻尼比最好在[0.1,0.4]范圍內;另一方面,微推力器的推重比較小,在承載較大質量推力器的基礎上,要能夠辨識微牛量級的推力,測量系統振動頻率均小于1 Hz。因此,通過增大阻尼比和固有頻率的方法,很難從根本上解決減小調整時間的問題,而且在實驗操作時,需要多次調節阻尼比和轉動慣量,并多次測試后方可達到測量要求,實施繁瑣,周期長。
本文針對微推力器開機工作時間與測量系統調整時間的匹配難題,提出了一種基于動態補償技術的微小穩態推力測量還原方法,利用動態補償技術改變測量系統的調整時間,并提出了一種補償器設計原則,通過實驗驗證了該方法的有效性。利用該方法可以從根本上縮短調整時間,不需要從硬件上做出修改,大大縮短實驗時間,并且具有很強的通用性。
水平扭擺基本結構如圖1所示。扭擺橫梁一端安裝推力器,另一端安裝配重,使得橫梁本身、推力器和配重的重心位于轉軸上。水平扭擺是二階質量-彈簧-阻尼系統,其基本運動方程為
(1)
式中:θ為扭擺橫梁相對于橫梁零點的轉角;J為測量系統轉動部件(包括橫梁、推力器和配重)相對于轉軸的轉動慣量;ζ為阻尼比;ωn為無阻尼自振角頻率;f(t)為推力器推力;b為推力器推力作用力臂。
因推力器輸出的穩態力上升時間非???,可以當做階躍力處理,即f(t)=F,則扭擺系統在階躍力作用下,系統響應為
(2)
(3)
式中:k為測量系統的扭轉剛度系數。
扭擺系統的典型響應如圖2所示。ts為扭擺的調整時間,t0為扭擺的1/4振動周期時刻。當穩態力作用時間Tf=t2,且穩態轉角采樣區間[ts,t2]大于振動周期的0.2倍時,就可以利用采樣區間內的穩態轉角均值計算穩態推力,且這種平均法帶來的相對誤差可以控制在一定水平,后文將詳細分析。但是當穩態力作用時間Tf=t1>t0,系統響應還沒達到穩態,此時無法利用穩態響應值計算推力。文獻[18]給出了一種基于杜哈梅積分法實時計算微推力的方法,但該方法的精度取決于采樣頻率和噪聲,且真值的收斂速度取決于測量系統的固有頻率。因此,本文提出利用動態補償技術改變測量系統的動態特性,減小測量系統的調整時間。
動態補償的原理[19]是:在原測試系統傳遞函數上增加一個補償環節,使得總的傳遞函數達到理想狀態,從而達到改善系統動態特性的目的。補償環節以辨識建模得到的模型為依據,設計出一種動態補償濾波器,與原來的測試系統相串聯,使級聯補償器后系統總的動態性能滿足使用要求。根據扭擺測量系統特點,基于動態補償技術的推力測量方案如圖3所示。
首先,將標定力加載到測量系統上,根據測量系統在標定力下的響應,進行系統動態模型辨識;然后,根據系統動態模型設計動態補償濾波器;最后,將推力器推力加載到測量系統上,將測量系統在推力作用下的響應,經過動態補償濾波器,還原出推力大小。
動態模型辨識的方法較多[19],這里不予贅述。本文采用的動態模型辨識方法的基本思路為:首先,根據模型階次估計準則,利用同時辨識模型階次和參數的方法確定模型階次;其次,根據輸入輸出數據,利用最小二乘法估計模型參數以作為迭代初值;最后,利用特殊白化濾波器的廣義最小二乘法估計模型參數。
圖3所示的測量過程中,動態補償濾波器是關鍵,進行軟件補償即可。常采用的設計方法有直接選擇等效系統法和零極點配置法,還有眾多研究人員提出了粒子群算法、神經網絡等復雜的補償器設計方法。本文所涉及的微推力測量系統屬于低階系統,考慮采用基于直接選擇等效系統法的動態補償濾波器設計方法。
直接選擇等效系統法設計動態補償數字濾波器的基本思想是:根據系統動態標定的瞬態響應曲線,由經典控制理論知識,直接得到描述系統特性的動態數學模型,對其模型進行動態分析,當其動態特性不滿足測量要求時,給系統串聯一個動態補償數字濾波器,使得串聯后構成的等效系統能夠滿足測量系統的動態特性。
取fr=kθ/b,則式(1)可改寫成
(4)
測量系統傳遞函數為
(5)
式(5)可串連一個阻尼比與固有頻率均相同的二階微分環節,為了減少高頻干擾,同樣可以加上一個二階低通濾波器,即將動態補償濾波器設計成具有下述傳遞函數的形式:
(6)

動態補償濾波器與測量系統串聯后的等效系統傳遞函數為
(7)

根據式(2),設在時間區間[T0,T1]內,扭擺系統在階躍推力作用下的響應均值為
(8)

利用穩態轉角均值計算穩態轉角的相對誤差為
|e-ζωnT1sin(ωdT1+2α)-e-ζωnT0·
(9)
設T0=k02π/ωd(k0為T0的周期倍數因子),T1=k12π/ωd(k1為T1的周期倍數因子),且ωd(T1-T0)≥x1≥1,ζωnT0≥x2?1時,式(9)變為
(10)

根據上述分析,當x1≥10和x2≥5,扭擺系統阻尼比不同時,采樣區間和相對誤差如表1所示。由表1可以看出,具有一定阻尼比的扭擺系統在推力作用下,只要進入穩態的時間區間至少能覆蓋某一區間[ts1,ts2],就可以將穩態轉角均值代替穩態扭轉角時的相對誤差控制在可忽略的水平。
表1不同阻尼比下的采樣區間和相對誤差
Table1Samplingrangesandrelativeerrorsfordifferentdampingratios

ζ采樣區間[ts1,ts2]δ/%0.2[3.898 5T,5.49T]≤0.076 10.3[2.530 4T,4.122T]≤0.070 30.4[1.823 3T,3.414 9T]≤0.068 20.5[1.378 3T,2.969 9T]≤0.067 60.6[1.061T,2.652 6T]≤0.067 40.7[0.811 9T,2.403 4T]≤0.067 4

(11)
綜上所述,當推力實際工作時間大于扭擺系統周期的0.25倍,且小于扭擺系統調整時間時,可通過動態補償方法建立新的等效測量系統,使得推力作用在等效測量系統上時,能夠利用穩態響應均值計算推力大小。由2.2節和3.1節可知,根據補償后輸出的均值計算穩態力的誤差主要來自3個方面:
1) 補償器輸入fr存在的誤差
由fr=kθ/b可知,fr的誤差來源于θ、k和b,根據文獻[20]可知,由于環境等干擾會產生測量噪聲,造成實際系統響應數據在平均位置曲線附近上下波動,可采用正交多項式局部滑動擬合方法,對實際系統響應測量值進行平滑降噪處理,確定系統響應的平均位置曲線。剛度系數k為具有置信區間的標定值,b存在誤差區間的測量力臂值。因此,fr的誤差可認為來自于k和b,且k越大、b越小,fr越大。
2) 補償器Hd(s)中系統參數ζ和ωn帶來的誤差
系統參數ζ和ωn均為原始測量系統具有置信區間的標定值。其中,ζ主要影響Hd(s)系統穩定的快速性,不影響穩態響應。因此,僅考慮系統參數ωn對補償器穩態輸出的影響,且ωn越大,穩態響應越大。
3) 補償后輸出在[0.811 9T′,2.403 4T′]區間內取均值帶來的誤差
根據3.1節的結論,在阻尼比為0.7的條件下,在[0.811 9T′,2.403 4T′]區間內取均值,利用穩態轉角均值代替穩態扭轉角時的相對誤差≤0.067 4%,該誤差可以忽略。
綜合上述分析,提出動態補償下的穩態推力還原步驟如下:
步驟1標定扭擺測量系統的系統參數[21],獲取扭擺測量系統參數及其置信區間ζ±Δζ、ωn±Δωn、T±ΔT、k±Δk,判斷測量系統的調整時間ts。
步驟2如果穩態力作用時間0.25T≤Tf≤ts,進入步驟3,否則直接根據扭擺響應的穩態幅值計算穩態推力大小。



利用文獻[22]中構建的水平扭擺測量系統進行實驗驗證,該測量系統由水平面調節機構、水平扭轉機構、位移測量裝置、電磁標定裝置和阻尼裝置組成。水平面調節機構用于提供整個測量機構的安裝底座支撐以及參考水平面。水平扭轉機構的主要機構是水平橫梁,用于安裝推力器。水平面調節機構和水平扭轉機構之間利用無摩擦樞軸連接,從而推力器工作時水平橫梁能夠在水平面內發生轉動。位移測量裝置用于測量水平橫梁的位移,在扭轉角小于5°時,可以利用水平橫梁的位移量代替扭轉角度。電磁標定裝置為非接觸電磁力產生裝置,由永磁鐵部分和線圈部分組成,永磁鐵部分安裝在水平橫梁上,線圈安裝在固定底座上,當線圈通特定電流時,即可產生電磁力使水平橫梁發生轉動,標定裝置的電磁力采用電子天平計量,作為測量系統的標定源。阻尼裝置用于提供測量系統的阻尼。對該測量系統的彈性軸、標定力裝置和位移測量裝置進行了調整,量程擴展到10 mN,放置在真空模擬環境中的測量系統如圖5所示。
首先對測量系統進行標定。真空度為3.1×10-4Pa,真空艙內溫度為20 ℃,位移傳感器采樣頻率為520.8 Hz,電磁阻尼電流設置為0.4 A。對扭擺測量系統分別加載和卸載0.139、0.283、0.601、0.890、1.179 mN的電磁力,電磁力作用力臂為0.5 m,傳感器測量力臂為0.59 m,測量系統響應如圖6所示。由于真空系統振動的影響,從局部放大圖圖7可以看出,測量系統響應存在高頻噪聲,由頻譜分析可知,噪聲集中在47 Hz和 35 Hz。對數據進行濾波處理后,根據文獻[21]的處理方法,可得扭擺系統的系統參數及95%置信區間如表2所示。
在測量系統標定后,利用標定力模擬實際推力器輸出的穩態力,設定標定力幅值為500.469 μN,控制穩態力工作時間分別為3、10、20 s。按3.3節的步驟,設計了相應的補償濾波器,取均值獲取了相應的推力大小,測量系統的實際響應(取的是位移的相對值)以及補償后的輸出如圖8所示(圖8中相關參數的取值為k=0.381 38 N·m/rad、b=0.5 m、ωn=0.783 63 rad/s)。根據3.3節的步驟,獲得的推力大小范圍及與輸入的偏差如表3所示。

表2 系統參數及95%置信區間Table 2 System parameters and 95% confidence interval

輸入推力/μN工作時間/s系統周期倍數補償后輸出/μN相對誤差/%500.46930.365 4[492.313 5,508.101 2][-1.63,1.53]500.469101.218 1[492.441 7,506.817 1][-1.6,1.27]500.469202.436 2[492.225 8,507.952 4][-1.65,1.5]
由圖8(a)可以看出,在穩態力作用下,由于力作用時間以及測量系統的特點,測量系統均未達到穩定狀態,無法利用扭擺的穩態位置相對變化量換算推力。由圖8(b)可以看出,利用動態補償技術后,等效系統最終的輸出能夠很快達到穩定狀態,并利用穩態平均值還原出推力大小。從表3可以看出,根據本文步驟能夠給出還原推力大小范圍,且與輸入的偏差能保證在2%以內。因此,當穩態推力實際工作時間大于扭擺系統周期的0.25倍,且小于扭擺系統調整時間時,可以利用本文提出的補償方法還原出穩態推力大小。
新型微推力器的力學性能測量對推力器的工作機理、設計和應用都具有很強的工程價值。由于新型推力器研究和設計階段的特殊性,有時推力器穩態推力輸出時間較短,與機械式測量裝置的穩態調整時間匹配困難,無法利用穩態推力與測量系統穩態響應成線性關系這一原理還原推力大小。本文提出一種基于動態補償技術的方法,詳細分析了扭擺穩態響應與穩態推力測量的關系,提出了補償器的設計原則,總結了推力還原步驟,并通過實驗驗證了該方法的有效性。具體研究結論如下:
1) 穩態推力長時間作用在水平扭擺式微小推力測量系統時,系統響應會進入穩態,如果在穩態響應某一區間內采樣穩態轉角值,可以利用穩態轉角均值代替穩態扭轉角,從而通過穩態扭轉角與穩態推力的線性關系計算穩態推力大小,而且這種計算方法帶來的相對誤差不超過0.1%。
2) 實驗結果表明,當推力實際工作時間大于扭擺系統周期的0.25倍,且小于扭擺系統調整時間時,可利用本文提出的基于動態補償技術的微小推力還原方法,建立新的等效測量系統,使得推力作用在等效測量系統上時,能夠利用穩態響應均值計算推力范圍。
3) 基于動態補償技術的微小推力還原方法適用于推力實際工作時間大于扭擺系統周期的0.25倍,且小于扭擺系統調整時間的情況。所以,一方面在設計測量系統時,在兼顧量程等指標的前提條件下,盡量使得扭擺系統的周期小;另一方面采用動態補償方法還原推力不需要修改硬件或進行復雜的硬件調整,適用于類似水平扭擺的二階測量系統,具有很強的通用性。