吳 靖,劉湘一,宋山松
(海軍航空大學(xué)航空基礎(chǔ)學(xué)院,山東 煙臺 264001)
直升機(jī)旋翼/機(jī)體耦合動態(tài)響應(yīng)仿真對于其動穩(wěn)定性分析[1],特別是非線性動穩(wěn)定性分析[2-4]等動力學(xué)問題具有重要意義,在通過飛行及試驗(yàn)難以獲得真實(shí)數(shù)據(jù)樣本的情況下,也可為直升機(jī)旋翼故障診斷研究提供大量仿真數(shù)據(jù)樣本[5-7],運(yùn)用軟件實(shí)現(xiàn)快速準(zhǔn)確的動態(tài)響應(yīng)仿真則是其中的關(guān)鍵。
在進(jìn)行直升機(jī)旋翼機(jī)體耦合動力學(xué)仿真時(shí),需要將位移量和速度量代入動力學(xué)模型中,求解加速度量,然后進(jìn)行迭代[8]。在不進(jìn)行簡化的情況下,槳葉一點(diǎn)的速度及加速度計(jì)算公式非常復(fù)雜[9],計(jì)算氣動力、慣性力以及氣動力矩、慣性力矩時(shí),需要進(jìn)行積分,給出其計(jì)算公式比較困難。運(yùn)用MATLAB等軟件進(jìn)行仿真時(shí),可以通過給加速度定義字符變量,在動力學(xué)模型中輸入位移量和速度量后,對加速度量進(jìn)行數(shù)值求解,但定義字符變量之后的求解速度很慢,每步求解的時(shí)間大約為數(shù)秒,步長為10-3s,計(jì)算20s的響應(yīng),需要花費(fèi)十幾個(gè)小時(shí)。對于直升機(jī)旋翼/機(jī)體耦合動態(tài)響應(yīng)仿真分析,特別是基于動力學(xué)模型的直升機(jī)旋翼故障診斷研究,需要大量數(shù)據(jù)樣本的情況,采用此種方法實(shí)現(xiàn)起來非常困難。
因此,在不進(jìn)行小角度線性簡化的前提下,保留運(yùn)動非線性,并通過分離加速度變量的方式給出一種快速的直升機(jī)旋翼/機(jī)體耦合動態(tài)響應(yīng)仿真方法,用于對其動態(tài)響應(yīng)進(jìn)行仿真分析。
基于矩陣運(yùn)算推導(dǎo)了槳葉任一點(diǎn)p在槳葉坐標(biāo)系xbybzb中的速度和加速度為
(1)
(2)
式中,r1為p點(diǎn)距揮舞/擺振鉸的距離,e為揮舞/擺振鉸外伸量,xc為槳轂中心距機(jī)體重心縱向距離,h為槳轂中心距機(jī)體運(yùn)動軸距離,L1~L6為變量矩陣,由基礎(chǔ)的坐標(biāo)變換矩陣運(yùn)算而來,直接展開形式復(fù)雜,不利于后續(xù)的積分計(jì)算氣動力和慣性力,一般都會進(jìn)行相應(yīng)的小角度線性簡化,為計(jì)入運(yùn)動的非線性,保留其如式(1)、(2)的矩陣運(yùn)算形式。
槳葉作用于揮舞/擺振鉸的力矩包括慣性力矩、彈簧力矩、結(jié)構(gòu)阻尼力矩及氣動力矩等,因此,第k片槳葉的揮舞運(yùn)動及擺振運(yùn)動方程為
(3)
(4)
(5)
(6)
式中,Ix、cx、kx分別是機(jī)體在滾轉(zhuǎn)方向上的慣性矩、阻尼和剛度;Iy、cy、ky分別是機(jī)體在俯仰方向上的慣性矩、阻尼和剛度。
作用在旋翼上的氣動力是非定常的,對于低頻振動的直升機(jī)來說,用動力入流模型能較好地描述非定常氣動力的作用。
用擴(kuò)展的Pitt-Peters動力入流模型[10]來描述非定常氣動力,其動力入流方程為
(7)
(8)
(9)
式中,L41、L51和L61為只含位移量和速度量的矩陣,不含加速度量,L42、L52和L62可以表示為
(10)
式中,L421、L422、L423、L424、L523、L524、L623和L624為只含位移量和速度量的矩陣。
在仿真時(shí),根據(jù)式(1)和式(9)中的非加速度變量項(xiàng)進(jìn)行積分計(jì)算氣動力、慣性力以及氣動力矩、慣性力矩,可用于計(jì)算式(8)中的b,根據(jù)式(9)中的加速度變量項(xiàng)進(jìn)行積分可計(jì)算式(8)中的A,再根據(jù)式(8)可計(jì)算加速度量,從而實(shí)現(xiàn)迭代求解,完成動態(tài)響應(yīng)仿真。
所用模型為美國NASA采用的無鉸旋翼模型,旋翼、機(jī)體模型的參數(shù)主要取自文獻(xiàn)[12],如表1所示。旋翼設(shè)定轉(zhuǎn)速為500r/min,槳葉初始安裝角為6°,來流角為0°。
表1 旋翼及機(jī)體模型參數(shù)
采用的仿真軟件為MATLAB,迭代算法為四階龍哥庫塔算法,采用分離加速度變量的方法和定義加速度字符變量的方法分別進(jìn)行仿真,步長為10-3s,分別計(jì)算20次。兩種方法每步計(jì)算時(shí)間分布如圖1所示。
圖1 計(jì)算用時(shí)對比
由圖1可知,仿真時(shí),采用分離加速度變量的方法每步計(jì)算平均用時(shí)為0.0033s,而定義加速度字符變量的方法每步計(jì)算平均用時(shí)為3.24s,前者僅為后者的1‰。采用分離加速度變量的方法計(jì)算,步長為10-3s,歷程為20s的動態(tài)響應(yīng)只需要1min左右,可以快速為基于動力學(xué)模型的旋翼故障診斷研究提供大量數(shù)據(jù)樣本。
以旋翼槳距不平衡為例,驗(yàn)證所提仿真方法的有效性。設(shè)定直升機(jī)第1片槳葉兩種槳距不平衡的情況,其安裝角分別對應(yīng)5.9°和6.1°,各槳葉揮舞運(yùn)動的動態(tài)響應(yīng)如圖2所示,為與槳距平衡的情況進(jìn)行對比,圖中給出了第1片槳葉槳距平衡時(shí)的揮舞響應(yīng)(旋翼槳距平衡時(shí)各槳葉揮舞角相同,旋翼共錐),機(jī)體響應(yīng)如圖3所示,槳距平衡時(shí),機(jī)體滾轉(zhuǎn)和俯仰角均為0。
圖2 槳葉揮舞響應(yīng)
圖3 機(jī)體響應(yīng)
由圖2可知,旋翼槳距不平衡時(shí),各槳葉揮舞不一致,引起旋翼不共錐,不平衡槳葉自身偏離平衡值最嚴(yán)重,約為2.76%。不平衡槳葉對其它槳葉也有影響,其它槳葉中后續(xù)第1片槳葉,即第4片槳葉所受影響最大,約為0.4%,后續(xù)第3片槳葉,即第2片槳葉所受影響最小,約為0.18%。分析可知,對于安裝角小于初始安裝角的情況,不平衡槳葉槳距減小,產(chǎn)生向下的誘導(dǎo)速度減小,后續(xù)槳葉在誘導(dǎo)速度減小的情況下,升力增加,揮舞值增加,如圖2(a)所示;而安裝角大于初始安裝角的情況,結(jié)果相反,使得后續(xù)槳葉揮舞值減小,如圖2(b)所示。
由圖3可知,槳距不平衡,旋翼作用在機(jī)體上的力矩不平衡,引起機(jī)體振動。由圖3(a)和(b)對比可知,第1片槳葉安裝角小于和大于初始安裝角兩種情況使得作用在機(jī)體上的不平衡力矩方向相反,因此引起機(jī)體的動態(tài)響應(yīng)相位相差180°。算例中引起機(jī)體的振動幅值較小,滾轉(zhuǎn)幅值僅為0.0019°,俯仰幅值為0.00037°。經(jīng)計(jì)算第1片槳葉槳距不平衡為10%,即安裝角為5.4°或6.6°時(shí),引起機(jī)體的滾轉(zhuǎn)幅值也僅為0.011°,俯仰幅值為0.0022°。
直升機(jī)旋翼機(jī)體耦合動力學(xué)模型中,為計(jì)入運(yùn)動的非線性,保留了槳葉任一點(diǎn)速度和加速度計(jì)算公式中的變量矩陣形式,使用MATLAB等軟件進(jìn)行仿真時(shí),采用定義加速度符號變量的方式,計(jì)算時(shí)間過長,而采用所提出的分離加速度變量的方法,計(jì)算時(shí)間可大大縮短,能快速地為基于動力學(xué)模型的旋翼故障診斷研究提供數(shù)據(jù)樣本。采用分離加速度變量的仿真方法,對槳距不平衡故障的情況進(jìn)行了計(jì)算,結(jié)果表明該方法能有效地將槳距不平衡對旋翼和機(jī)體動態(tài)響應(yīng)的影響規(guī)律仿真出來,可為直升機(jī)旋翼不平衡故障診斷研究提供參考。