垂直著陸中直升機旋翼動力學行為研究

胡國才，吳靖，劉湘一，劉書巖

海軍航空大學 航空基礎學院，煙臺 264001

目前大、中型鉸接式或彈性軸承旋翼直升機廣泛采用剛度低、阻尼大的液壓阻尼器，提高旋翼擺振面阻尼以避免直升機發生地面共振[1]。液壓阻尼器的阻尼力與活塞速度或槳葉擺振速度的平方近似呈比例關系，為了限制使用中阻尼器的動態載荷，以滿足阻尼器和旋翼結構的強度及壽命要求，阻尼器內部一般設置定壓安全活門，當阻尼力超過活門開啟壓力時安全活門打開，從而限制載荷的進一步增加[2-3]。設置安全活門帶來的主要問題是直升機起降及地面運轉時，若旋翼發生大的擾動安全活門開啟后，阻尼器的等效阻尼會迅速下降，對地面共振產生不利影響[4]。因此，在正常使用中不會發生地面共振的直升機可能在粗暴著陸時因過大的初始擾動而發生地面共振，這在實際使用中仍時有出現。

直升機地面共振詳細設計要求GJB720.5—89規定，如果直升機采用顯著影響地面共振的非線性元件，則應給出不發生地面共振的初始擾動幅值，并在設計使用包線內可能遇到的擾動下應不發生地面共振。帶有安全活門的液壓阻尼器具有嚴重的非線性，采用等效線性系統，可方便地確定滿足地面共振穩定性的初始擾動幅值。但在設計使用包線內可能遇到的擾動幅值，在直升機設計時難以預估，這為直升機地面共振設計帶來了不確定性。雖然也可以通過半主動或主動控制的方法使直升機工作在全穩定狀態，如直接使用半主動式阻尼器磁流變阻尼器[5]，或者增加機體狀態反饋[6]等來達到對直升機地面共振的自適應控制[7]。但是，目前使用最多的還是結構相對簡單的被動式阻尼器，這就要求盡量提高阻尼器的等效阻尼避免直升機產生地面共振，但要滿足結構強度與壽命的要求，必然付出結構重量代價。要合理解決以上兩者之間的矛盾，前提是需要準確預估直升機使用包線內可能遇到的擾動大小。

經驗表明，直升機單輪粗暴著陸或著艦產生的擾動可能最為嚴重。現有關于直升機著陸和著艦方面的研究，重點關注起落架及機身結構的載荷、強度及疲勞壽命。許多學者采用多體動力學模型，用經典方法預估直升機著陸載荷[8]；用剛-柔混合有限元模型計算著陸載荷及機體和起落架組件的應力[9-10]。文獻[11]在全局有限元模型中改進了關鍵部位模型，研究粗暴著陸時關鍵部位的疲勞應力。但關于直升機著陸或著艦時激起的旋翼擾動及動力穩定性問題的研究報道很少。

本文通過直升機著陸動力學建模與仿真，從理論上分析垂直著陸時起落架撞擊地面對旋翼動力學行為的影響，希望為預估粗暴著陸時阻尼器的擾動幅值提供一種有效的方法。

1 分析模型

假設直升機限定在垂直著陸狀態，考慮質心垂向運動、機體繞質心的滾轉和俯仰等3個自由度，不考慮質心在水平面內的移動和機體的偏航自由度。因機體剛度比起落架剛度要大得多，可將機體看做剛體。對于鉸接式或彈性軸承旋翼來說，對阻尼器軸向運動及地面共振起主要作用的是槳葉繞揮舞鉸和擺振鉸的轉動，故本文只考慮剛硬槳葉的揮舞、擺振和變距自由度。起落架輪胎底面接觸地面后，輪胎和緩沖支柱壓縮產生彈性力和阻尼力并作用到機體上，對機體產生力和力矩，引起機體和旋翼的擾動運動。直升機垂直著陸動力學分析模型如圖1所示。

圖1 直升機垂直著陸模型Fig.1 Vertical landing model of helicopter

圖1中，Ogxgygzg為地面坐標系，也是慣性坐標系；Ofxfyfzf為機體坐標系，跟隨機體運動，原點Of位于直升機質心處，Ofxf指向機身尾部；Ohxhyhzh為槳轂不旋轉坐標系，跟隨機體運動，原點Oh位于槳轂中心，Ohxh指向機身尾部；Obxbybzb為槳葉坐標系，跟隨槳葉運動，原點Ob位于彈性軸承處，Obxb固定在槳葉上。設旋翼槳葉片數為Nb，旋翼轉速為Ω(俯視逆時針旋轉)，第k片槳葉的方位角為ψk、揮舞角為βk(向上為正)、擺振角為ζk(逆旋轉方向為正)、變距角為θk，彈性軸承的外伸量為e，機體滾轉和俯仰歐拉角為γ和?，機體滾轉和俯仰角速度為ωx和ωy，旋翼軸到質心距離為xc，槳轂中心到質心高度為h。

任一片槳葉上距離彈性軸承r處的葉素在機體坐標系中的矢徑為

ρ=[xf,yf,zf]if

(1)

式中：if為機體坐標系單位向量，且

(2)

葉素在地面坐標系中的矢徑為

ρ=[xg,yg,zg]ig

(3)

式中：ig為地面坐標系單位向量，且

(4)

葉素在地面坐標系中的速度和加速度可由式(3)分別對時間求一次和二次導數得到。把葉素在慣性坐標系中的速度及加速度轉換到槳葉坐標系中，可以計算葉素的氣動力及慣性力。

著陸時起落架撞擊地面會引起機體和旋翼的擾動運動，為了提高分析精度，應計入非定常氣動力影響。對于旋翼低頻擾動運動，Pitt/Peters動力入流模型可反映非定常氣動力的影響[12]，這在旋翼/機體耦合動穩定性研究中得到了驗證[13]。考慮到下降速度對旋翼誘導速度的影響，采用擴展的動力入流模型[14]。設槳盤平面上任一點的無因次誘導速度為

v=v0+vs(r+e)sinψ+vc(r+e)cosψ

(5)

式中：v0為平均誘導速度；vs和vc分別為橫向和縱向入流系數。Pitt/Peters動力入流方程為

(6)

式中：M為顯在空氣質量矩陣；L為入流增益矩陣；CT、CL及CM分別為旋翼總的氣動升力系數、對槳轂中心的氣動滾轉及俯仰力矩系數。

前飛狀態時矩陣L為

(7)

式中：

直升機以速度vz垂直下降時，前進比μ=0，入流角as=-90°，上述系數改寫為：vT=v0-vz，vm=2v0-vz，α=90°。那么垂直下降時的矩陣L為

(8)

直升機接近地面時產生的地面效應，會降低旋翼的誘導速度，增加旋翼升力[15]。旋翼總距一定的情況下，旋翼升力增加會增大旋翼平均揮舞角和擺振角，但對旋翼周期揮舞和周期擺振的影響相對較小，因此可不考慮地面效應對旋翼擾動運動的影響。

采用彈性軸承旋翼系統，液壓阻尼器安裝時一般將一端用球鉸與槳轂連接，另一端用球鉸與槳葉或過渡件連接，這樣槳葉繞彈性軸承的揮舞、擺振和變距運動都會影響到阻尼器速度，故計算阻尼器速度時需計入幾何耦合的影響[16]。液壓阻尼器的力-速度曲線由試驗數據得到，但為了仿真計算的方便，一般采用雙折線來擬合試驗數據獲得理論計算公式。

起落架觸地時機體姿態角不大的情況下，沿緩沖支柱的軸向載荷對機體起主要作用，本文暫不考慮輪胎觸地后產生的側向力和縱向力。

根據達朗貝爾原理導出直升機垂直著陸的動力學方程，以一階向量表示為

(9)

機體的運動學方程為

(10)

圖2 直升機起落架模型Fig.2 Model of helicopter landing gear

為便于計算起落架著陸載荷，建立起落架參考坐標系。起落架觸地后壓縮或伸展時，緩沖支柱和輪胎均會產生非線性彈性力和阻尼力，分別以非線性彈簧S1、S2及非線性黏壺D1、D2表示，起落架模型如圖2所示。前起落架用下標“N”表示，左側主起落架用下標“ML”表示，右側主起落架用下標“MR”表示。

起落架參考點Or與直升機質心處于同一鉛垂線上，參考點垂直位移與直升機質心垂直位移之間的關系為zr=-z。起落架的垂直位移取決于參考點的垂直位移zr和機體姿態角。假定參考點到前起落架的距離為xN，到主起落架連線的距離為xM，主起落架輪距為2yM。

則起落架的垂向位移為

(11)

起落架的垂向位移對時間求導，可得起落架的垂向速度。

(12)

式中：Nw為單個起落架的輪胎數量。

緩沖支柱和輪胎的彈性力取決于壓縮量，為提高預估精度，采用靜壓縮試驗數據。緩沖支柱的阻力包含油液流經阻尼孔產生的阻力fd1h和活塞摩擦力fd1c。設ρ為油液密度、Cd為流量系數、A0為節流閥座的承壓面積、As為壓縮行程的阻尼孔面積、Ar為伸展行程的阻尼孔面積，則流體阻力為

(13)

活塞摩擦力fd1c大小與活塞的摩擦系數及緩沖支柱內部空氣壓力有關，其方向與活塞運動方向相反。影響摩擦力的因素較多，需要由試驗確定。由于缺少摩擦力試驗數據，本文暫不考慮摩擦力的作用。輪胎的阻尼力采用經驗公式[17]：

(14)

式中：Cw為輪胎當量阻尼系數。

起落架作用在機體上的垂向力Fz、滾轉力矩Mx和俯仰力矩My，即各起落架緩沖支柱的彈性力和阻尼力之和。僅考慮起落架的垂向載荷時，作用在機體上的力和力矩分別為

(15)

起落架過載系數定義為

ngi=(fs1i+fd1i)/p0ii=N,ML,MR

(16)

式中：p0i為直升機在零升力停機狀態時各個起落架的地面支反力。

2 著陸仿真與結果分析

2.1 計算方法

采用4階Runge-Kutta法對直升機著陸運動方程進行數值積分，模擬直升機垂直著陸過程中旋翼槳葉、機體及起落架的動態響應。給定不同的高度和旋翼總距，可獲得不同的著陸速度；設置不同的周期變距可使直升機產生不同的姿態角。首先將直升機高度固定，設置旋翼總距和周期變距，對運動方程進行初值積分，經過一段時間后旋翼運動和機體姿態趨于平衡狀態。直升機趨于平衡后，放松質心垂向自由度，使直升機垂直下降，當起落架觸地時，機體在起落架支反力作用下產生擾動運動，同時激起旋翼槳葉的揮舞和擺振運動。對全部槳葉揮舞和擺振響應進行多槳葉坐標轉換，得到旋翼整體模態響應，對其進行快速傅里葉變換(FFT)及濾波運算可以辨識各個模態。

旋翼和機體運動方程的積分變量為方位角ψ，起落架動力學方程的積分變量為時間t，兩者之間的關系為ψ=Ωt。最終的載荷及動態響應統一以時間t或方位角ψ為自變量的形式輸出。

2.2 結果分析

以某彈性軸承旋翼直升機為例，旋翼和機體的主要參數如表1所示。

不穩定流場、或駕駛員操縱不當等因素都會導致直升機單輪粗暴著陸，表現為直升機以一定傾斜角和較大速度接地。作為算例，假定讓右起落架最先接地，接地速度大約為2 m/s。為模擬這種著陸狀態，直升機需從較高的高度開始下降，駕駛桿需向右、向后拉桿。設直升機從高度為1 m(輪胎底面距離地面的平均高度)開始垂直下降著陸，旋翼總距θ0=6°、周期變距θc=-2°及θs=0.5°，重心位于中位xc=0。圖3為機體滾轉和俯仰角響應，圖4為起落架過載系數隨時間的變化曲線。

圖3顯示，開始階段機體呈現向右傾斜約3.5°、抬頭1°的姿態，該姿態角與右壓桿(θc=-2°)及后拉桿(θs=0.5°)的配平特性相吻合。運動方程積分時，旋翼周期揮舞和周期擺振的初值設置為零，故要經過一定時間才能達到配平狀態，在此過程中機體姿態角隨時間會存在一定波動，機體滾轉角的波動更為明顯。

表1 旋翼及機體主要參數Table 1 Main parameters of rotor and body

圖3 機體滾轉和俯仰響應Fig.3 Body roll and pitch responses

從t=1.57 s開始下降，大約在t=2.34 s時右起落架輪胎觸地，隨后起落架載荷迅速增加(見圖4)，機體在該載荷作用下產生左滾和低頭力矩，使機體迅速向左滾轉并低頭，如圖3所示。大約經過0.11 s，左起落架和前起落架先后(幾乎同時)觸地。右起落架最大過載系數為1.44，左、前起落架最大過載系數分別為1.58和1.53，比右起落架分別高9.7%和6.2%。著陸后期，機體略有向右傾斜和抬頭的姿態，左、右起落架載荷也不相同，主要是因為直升機著陸后旋翼槳距沒有回零的緣故。

圖4 起落架過載系數隨時間變化曲線Fig.4 Curves of overload factor of landing gears vs time

圖5 槳葉擺振角隨方位角的變化曲線Fig.5 Curves of blade lag angle vs azimuth

圖5顯示了各片槳葉的擺振響應隨第1片槳葉方位角的變化曲線。可以看到槳葉的平均擺振角為6°左右，定常周期擺振的幅度約為0.18°。著陸過程中各片槳葉產生了不同程度的擺振擾動，擾動幅度大約在0.5°～1.2°之間。

從仿真結果看到，右起落架觸地前，各片槳葉按照固定的相位差(60°)以轉速Ω作定常擺振運動。一個有意思的現象是，全部起落架觸地引起的擾動擺振響應與定常擺振響應疊加后，至第1片槳葉處于330°方位角時，各片槳葉均近似處于平均擺振角位置。當然，這個現象與各起落架觸地時機、機體姿態等諸多因素有關，不具有普遍意義，但為下面的力學分析提供了便利。

圖6表示第1片槳葉處于330°方位(坐標軸xf為0°、方位角按順時針確定)，前、左起落架接地后對機體產生向后及向右的力矩，機體擾動使槳轂中心產生向后及向右的加速度，對槳葉施加向前和向左的慣性力作用。

從圖6可以看出，此時1號和2號槳葉在慣性力矩作用下將逆旋轉方向擺振，擺振角將增加(因規定了逆旋轉擺振為正)，綜合兩個慣性力的作用，1號槳葉應比2號槳葉向后擺振更劇烈；4號和5號槳葉將順旋轉方向擺振，擺振角將減小，并且4號槳葉應比5號槳葉向前擺振更劇烈。對于3號和6號槳葉來說，由滾轉和俯仰運動產生的慣性力矩起相互抵消的作用，擺振幅值小，而擺振方向與兩者的大小有關，圖5顯示的響應曲線與上述分析結果吻合。隨后的槳葉擺振運動因與機體運動、槳葉旋轉及定常擺振運動等因素存在復雜關系，難以一一分析清楚。

圖6 機體擾動引起的槳葉慣性力Fig.6 Blade inertial force caused by fuselage perturbation

全部起落架從2.34 s至3.0 s完成第1次壓縮-回彈運動，持續時間大約為0.66 s (見圖4)，之后起落架和機體處于小幅震蕩。因此，在著陸擾動期間難以將擺振后退型的瞬態響應分離出來，因而也就不能根據包絡線來識別其模態阻尼[18]，但可知著陸撞擊激起的擺振擾動大小。

對旋翼擺振余弦分量進行2 Hz低通、4 Hz高通和3～4 Hz帶通濾波可知，除了擺振后退型模態響應外，還存在其他頻率的響應成分，具體如圖8所示。

圖7 旋翼擺振的余弦分量及2～3 Hz帶通濾波后的響應Fig.7 Cosine component and 2-3 Hz band pass filtered rotor lead-lag response

很明顯，幅度最大、頻率最低的響應成分是由起落架觸地后壓縮-回彈引起的(實線)；高通濾波后僅存在頻率大約為4.5 Hz的響應成分(虛線)，因擺振前進型模態頻率ωζ P=4.46 Hz，可以確定這是擺振前進型響應，其幅值約為0.17°；頻率為3.5 Hz持續振動的響應成分是定常周期擺振，著陸時其幅值也稍有增加(點線)。

從以上仿真結果及分析可知，單輪著陸會激起旋翼擺振后退型響應，雖然其幅值并不大，但由于存在各種頻率成分的背景振動，尤其是低頻、大幅值的沖擊響應，液壓阻尼器將處于大速度狀態而打開定壓活門，進而使其等效阻尼下降。圖9為各片槳葉的液壓阻尼器軸向速度響應，圖10為該阻尼器的等效阻尼隨速度幅值的變化曲線[19]，定壓活門的開啟速度約為0.01 m/s。

由圖9看到，定常狀態下阻尼器的速度大約為0.012 m/s，著陸時各阻尼器速度峰值高達0.035～0.05 m/s，已遠超定壓活門的開啟速度，由圖10看到，阻尼器的等效阻尼將隨著速度的增加而急劇下降。

圖8 旋翼擺振響應中的各種成分Fig.8 Components in rotor lead-lag response

圖9 各阻尼器的軸向速度響應Fig.9 Axial velocity response of each damper

圖10 等效阻尼隨軸向速度的變化曲線Fig.10 Curves of equivalent damping vs axial velocity

3 垂直著陸對載荷及旋翼動力學影響

直升機著陸具有一定隨機性，其構型、著陸場地、大氣環境以及駕駛員操縱等因素都會影響著陸狀態[20-21]。為模擬著陸隨機性，針對本文的直升機構型，給定高度范圍為0.3～1.5 m、間隔為0.3 m，在一定范圍內隨機給出開始下降時間、旋翼總距和周期變距，開始下降時間范圍為1～2 s、總距范圍為6°～7°、周期變距范圍為-2°～2°。在每個高度，隨機著陸5次，輸出模擬參數為各起落架著陸速度和最大過載系數、槳葉揮舞和擺振動幅值以及阻尼器速度峰值，以散點圖的形式顯示。起落架觸地瞬間的機體姿態角及角速度等動態參數，取決于著陸高度、開始下降時間及旋翼變距等因素，仿真結果如圖11所示。

圖11(a)和圖11(b)顯示了直升機不同高度隨機著陸時起落架的著陸速度及過載系數，反映了著陸的嚴重程度。

從圖11(a)和圖11(b)可以看到，著陸速度及過載與著陸高度具有一定的相關性。說明著陸高度對起落架著陸速度和過載有重要影響。同時也看到，在同一高度隨機著陸時，著陸速度和過載系數具有相當大的分散性，表明除了著陸高度外，觸地時機體的姿態角及角速度也是影響著陸速度和過載的重要因素。主起落架著陸速度的分散性大約為0.8 m/s，基本上不隨著陸高度而變。而過載系數的分散性卻隨著陸高度的增加而有所擴大。

圖11(c)和圖11(d)分別表示著陸時槳葉揮舞角和擺振角的最大動幅值。

圖11 仿真結果圖Fig.11 Results of simulation

從圖11(c)和圖11(d)可見，隨著陸高度的增加，揮舞角和擺振角的最大動幅值的分散性明顯增大。說明著陸高度越大，隨機著陸時槳葉可能會遇到更大的擾動。與揮舞相比，槳葉擺振擾動的分散性更為顯著，原因是擺振面氣動力遠小于揮舞面氣動力，因此著陸時機體擾動引起的擺振氣動阻尼也相應地小得多，而擾動產生的慣性力對槳葉擺振運動的作用要大于對揮舞運動的作用。從圖11(d)看到，高度為0.3 m隨機著陸時，擺振角幅度在0.41°～0.61°之間，而在1.5 m高度隨機著陸時，擺振幅度范圍擴大至0.67°～1.7°。

圖12表示隨機著陸時液壓阻尼器的軸向速度峰值隨著陸高度的變化。

圖12 液壓阻尼器軸向速度峰值隨著陸高度的變化Fig.12 Maximum of hydraulic damper axial velocity vs landing altitude

由圖12看到，在不同高度隨機著陸時，阻尼器軸向速度峰值也存在很大的分散性。由于阻尼器與槳葉揮舞、擺振及變距自由度存在幾何耦合，其軸向速度主要取決于槳葉的擺振速度外，還與槳葉揮舞和變距運動有關。因此阻尼器的軸向速度峰值與圖14所示的槳葉最大擺振動幅值之間并不是嚴格的對應關系。不同高度隨機著陸時，阻尼器軸向速度峰值的最小值處于24～30 mm/s范圍內，最大值處于45～74 mm/s范圍內，說明粗暴著陸時，阻尼器的軸向速度峰值會顯著增加。

4 結 論

通過建立直升機垂直著陸動力學模型并進行數值仿真，研究了直升機垂直著陸中起落架、機體及旋翼的動力學行為，分析了著陸狀態對動態響應參數的影響規律，得到以下主要結論：

1) 直升機著陸時，起落架撞擊地面引發的機體滾轉和俯仰擾動是激起旋翼動態響應的主要原因。

2) 起落架觸地后的壓縮-回彈周期中，將激起各片槳葉不同的擺振幅值，并形成旋翼擺振后退型響應。

3) 隨機著陸時，起落架著陸速度及過載系數與著陸高度具有明顯的相關性，而機體姿態角及角速度的不同使得同一高度著陸時起落架的著陸速度及過載具有一定分散性。

4) 隨著陸高度增加，著陸姿態對槳葉揮舞和擺振動幅值的影響增加，尤其對槳葉擺振動幅值的影響更為顯著，隨著著陸高度增加，擺振動幅值的分散性也增加。

5) 不同高度隨機著陸時，阻尼器軸向速度峰值呈現很大的分散性，且由于阻尼器與槳葉揮舞、擺振及變距自由度存在幾何耦合，使得其軸向速度峰值的分散性隨著著陸高度的變化規律不明顯。