何澤鵬，畢世華，傅德彬，李云峰

(1.北京理工大學(xué) 宇航學(xué)院，北京 100081；2.吉林江機特種工業(yè)有限公司技術(shù)部，吉林 132021)

0 引言

長期以來，火箭、導(dǎo)彈發(fā)射擾動和彈道計算作為獨立的領(lǐng)域各自發(fā)展，取得了較好的成果。實際應(yīng)用中兩者互相耦合的現(xiàn)象也主要是通過弱耦合的方式進(jìn)行關(guān)聯(lián)，即先計算初始擾動，再進(jìn)行有擾動情況下的飛行彈道計算，而兩者直接耦合計算的研究較少。

魏靖彪和鄭華利[1]結(jié)合動量理論和渦流理論，推導(dǎo)了槳葉環(huán)量方程，采用廣義尾流模型計算了有旋翼尾流干擾和無旋翼尾流干擾兩種條件下的導(dǎo)彈初始彈道。車競等[2]通過建立機彈分離的湍流模型，并嵌入導(dǎo)彈的六自由度彈道方程，運用CFD計算出導(dǎo)彈在載機流場干擾下的非定常氣動力和初始彈道。王潔瑤等[3]利用準(zhǔn)平衡滑翔條件建立滑翔射程和飛行時間的解析估算公式，研究高超聲速遠(yuǎn)程導(dǎo)彈的彈道性能。陳陽陽[4]等針對起飛質(zhì)量對小型導(dǎo)彈總體參數(shù)設(shè)計的制約問題，以最小起飛質(zhì)量為目標(biāo)函數(shù)，應(yīng)用復(fù)合形法，對小型飛航導(dǎo)彈縱向平面運動的方案彈道進(jìn)行了總體參數(shù)優(yōu)化設(shè)計研究。駱連珍等[5-6]通過建立炮彈的多體動力學(xué)運動模型，研究發(fā)射裝置的振動和炮彈間隙對火箭彈產(chǎn)生的初始擾動影響。徐偉國等[7]采用多體動力學(xué)仿真分析和實驗驗證的方法，研究了火箭彈發(fā)射裝置的支撐剛度對初始擾動的影響。趙憲斌等[8]采用物理模型和軟件分析方法研究發(fā)射裝置液壓缸振動、發(fā)射仰角和箭架配合間隙等因素對火箭發(fā)射擾動的影響。這些文獻(xiàn)在各自領(lǐng)域研究中較為深入，對于發(fā)射擾動與彈道計算的耦合關(guān)聯(lián)較少。

本文針對火箭導(dǎo)彈發(fā)射擾動與初始彈道互相耦合引起的彈道散布的問題，提出發(fā)射擾動與彈道解算相耦合的計算分析模型，并以滾轉(zhuǎn)彈為實例，對模型的計算結(jié)果進(jìn)行分析，為相關(guān)研究提供參考。

1 耦合計算模型

1.1 耦合計算思路

火箭導(dǎo)彈發(fā)射和飛行過程中，通常包括多個階段，在發(fā)射初期，通?？煞譃榘l(fā)射約束段、無控飛行段和有控飛行段。在本研究中，重點關(guān)注發(fā)射約束段和無控飛行段。

在發(fā)射約束段，火箭導(dǎo)彈與發(fā)射裝置相互作用，產(chǎn)生復(fù)雜的動力學(xué)響應(yīng)過程，是產(chǎn)生發(fā)射擾動和彈道散布的主要因素之一；在無控飛行段，彈體會保持發(fā)射段受到的擾動，并受到低速氣動載荷的影響，對彈道散布產(chǎn)生一定影響。為有效地將發(fā)射擾動與無控飛行段彈道解算相耦合，主要采用如圖1所示的狀態(tài)建立計算分析模型。

(1)發(fā)射約束段，采用包含發(fā)射裝置及彈體的多體動力學(xué)模型，對發(fā)射動力學(xué)響應(yīng)進(jìn)行求解計算。

(2)無控飛行段，在前述多體動力學(xué)模型基礎(chǔ)上，將彈體承受的氣動載荷轉(zhuǎn)換至彈體坐標(biāo)系，利用動力學(xué)方法進(jìn)行彈道計算。

圖1 彈體飛行三階段

圖1中，F(xiàn)推為彈體受到發(fā)動機的推力；F擾為彈體受到發(fā)射架的擾動力；F氣為彈體受到氣動力；G為彈體的重力；M氣為彈體受到氣動力矩；V為彈體的速度。

1.2 多體動力學(xué)模型

在本研究中，發(fā)射擾動的計算和彈道解算均以多體動力學(xué)模型[9-11]為基礎(chǔ)。對于含N個物體作空間運動的多體系統(tǒng)，取一個慣性參考基e和每個物體的連體基ei，質(zhì)心相對于慣性基的坐標(biāo)(xi，yi，zi)與連體基相對于慣性基的歐拉四元數(shù)(λ0i，λ1i，λ2i，λ3i)構(gòu)成描述物體βi的笛卡爾坐標(biāo)陣為

qi=(xi，yi，zi，λ0i，λ1i，λ2i，λ3i)T

(1)

對于該多體系統(tǒng)，每個物體均可用式(1)描述位形，則描述整個多體系統(tǒng)位形的坐標(biāo)總數(shù)為7N。引入拉格朗日乘子，得到整個系統(tǒng)的動力學(xué)方程為

(2)

式中Z為廣義質(zhì)量陣；z為力陣；zI為慣性力陣；σ=(σ1…σs)T為約束方程對應(yīng)的拉格朗日乘子陣。

2 發(fā)射擾動與氣動載荷處理

2.1 發(fā)射擾動處理方法

發(fā)射過程中，發(fā)射裝置與彈體之間、發(fā)射裝置與平臺或地面之間，除傳統(tǒng)約束外，往往還包含大量的接觸關(guān)系，而這些接觸關(guān)系和狀態(tài)，又是引起發(fā)射擾動的主要因素之一。在研究模型中，利用運動副模擬物體間的相互約束關(guān)系，利用基于碰撞函數(shù)的接觸算法[12-14]模擬接觸關(guān)系，如圖2所示。

上述碰撞模型的接觸力計算公式可表示為

(3)

圖2 碰撞接觸模型

為滿足裝配要求，彈架配合通常留有一定的間隙，如圖3所示。在發(fā)射過程中，彈架配合間隙同樣是引起初始擾動的重要因素。本文為模擬發(fā)射約束段定向?qū)к壟c彈上滑塊間的配合間隙作用，引入隨機過載模型[15-16]，在多體動力學(xué)模型中將受間隙尺寸影響的隨機擾動力施加在彈體上，并通過接觸狀態(tài)模擬隨機擾動作用下的彈架動態(tài)響應(yīng)。

圖3 導(dǎo)軌與定向件的配合間隙

作用于彈體的隨機過載表示為

ny=A(H)×R(t)

(4)

式中A(H)為與間隙尺寸有關(guān)的過載幅值；R(t)為取值范圍在0和1之間的隨機數(shù)。

過載幅值采用式(5)計算：

(5)

式中h為彈體支腳與發(fā)射架之間的配合間隙；λ為導(dǎo)軌平度的波長；V1為彈體速度。

2.2 氣動載荷與彈道解算

本研究以多體動力學(xué)的模型為基礎(chǔ)，將彈體受到的氣動載荷均轉(zhuǎn)換到彈體坐標(biāo)系下，并利用動力學(xué)方法[17]進(jìn)行彈道計算。

為通過飛行狀態(tài)和氣動參數(shù)進(jìn)行氣動載荷的計算，并將其轉(zhuǎn)換到彈體坐標(biāo)系下，引入彈道解算的如下常用坐標(biāo)系：地面坐標(biāo)系Oxyz、彈體坐標(biāo)系Ox1y1z1、彈道坐標(biāo)系Ox2y2z2、速度坐標(biāo)系Ox3y3z3、適用于低速滾轉(zhuǎn)彈的準(zhǔn)彈體坐標(biāo)系Ox4y4z4和準(zhǔn)速度坐標(biāo)系Ox*y*z*。在多體動力學(xué)模型中，常用坐標(biāo)系間的關(guān)系如表1所示。

表1 氣動載荷參數(shù)的定義與求解

表1中，Vx、Vy、Vz分別為速度在地面坐標(biāo)系Ox、Oy、Oz軸上的分量；Vxz為導(dǎo)彈速度在Vx與Vz的矢量和。

利用上述坐標(biāo)和坐標(biāo)轉(zhuǎn)換角，可將不同坐標(biāo)系下的氣動載荷轉(zhuǎn)換到彈體坐標(biāo)系下。對于滾轉(zhuǎn)彈體，其坐標(biāo)轉(zhuǎn)換過程較為復(fù)雜，這里主要針對其進(jìn)行介紹，對于無滾轉(zhuǎn)彈體，可參考這類方法進(jìn)行轉(zhuǎn)換。

2.2.1 準(zhǔn)速度坐標(biāo)系下的氣動力

在準(zhǔn)速度坐標(biāo)系下，彈體受到的阻力X*，升力Y*以及側(cè)向力Z*可以表示為

(6)

2.2.2 準(zhǔn)彈體坐標(biāo)系下的氣動力矩

為求解準(zhǔn)彈體坐標(biāo)系下的力矩，需要用到準(zhǔn)彈體坐標(biāo)系下的姿態(tài)角速度的如下轉(zhuǎn)換關(guān)系式：

(7)

式中ωx1、ωy1、ωz1分別為彈體坐標(biāo)下的滾轉(zhuǎn)角速度、偏航角速度和俯仰角速度。

在準(zhǔn)彈體坐標(biāo)系下，彈體受到的滾動力矩Mx4、偏航力矩My4以及俯仰力矩Mz4可以表示為

(8)

2.2.3 彈體坐標(biāo)系下力與力矩轉(zhuǎn)換

在彈體坐標(biāo)系下，可利用如下關(guān)系式計算其受到的氣動力：

(9)

其中

O=(sinβ*)X*+(cosβ*)Z*

M=(cosα*cosβ*)X*+(sinα*)Y*-

(cosα*sinβ*)Z*

N=(-sinα*cosβ*)X*+(cosα*)Y*+

(sinα*sinβ*)Z*

在彈體坐標(biāo)系下，彈體承受的力矩可由準(zhǔn)彈體坐標(biāo)系下的力矩的關(guān)系式獲得：

(10)

3 基于滾轉(zhuǎn)彈的實例分析

3.1 實例模型及參數(shù)

為應(yīng)用多體動力學(xué)模型進(jìn)行發(fā)射擾動和初始彈道的耦合計算，這里以滾轉(zhuǎn)彈為例進(jìn)行計算分析。該滾轉(zhuǎn)彈簡化模型包括地面、發(fā)射架導(dǎo)軌、彈體支腳和彈體四部分，各部件間的關(guān)系見圖4。

圖4 多體動力學(xué)簡化模型

在計算模型中，以地面坐標(biāo)系為全局坐標(biāo)系，并在彈體質(zhì)心位置建立多個局部坐標(biāo)系，彈體上的局部坐標(biāo)系方向與彈體坐標(biāo)系相同。所有坐標(biāo)系均為右手坐標(biāo)系。滾轉(zhuǎn)彈的質(zhì)量為11.7 kg，彈體的赤道轉(zhuǎn)動慣量為0.032 kg·m2，極轉(zhuǎn)動慣量為0.318 kg·m2，彈體質(zhì)量為11.7 kg，彈長1 m，質(zhì)心位于距彈頭640 mm的軸線上，發(fā)射軌的有效長度為200 mm，發(fā)射角度為6°，發(fā)射前已將彈體放置在導(dǎo)軌上，并采用同時離軌方案。計算用到的接觸參數(shù)設(shè)置參考文獻(xiàn)[18]。

3.2 耦合模型典型響應(yīng)

利用商用多體動力學(xué)軟件ADAMS對耦合實例模型進(jìn)行求解計算，獲得了耦合模型作用下彈體典型的響應(yīng)如下所述。

圖5給出發(fā)射2 s內(nèi)彈體俯仰角、彈道傾角和準(zhǔn)攻角變化曲線。點火初期，由于彈體軸向速度很低，支腳與導(dǎo)軌接觸力引起的Y向速度相對較大，彈體姿態(tài)角出現(xiàn)瞬時震蕩，但對后續(xù)計算結(jié)果影響很小；隨著彈體軸向速度增加，初始俯仰角和彈道傾角穩(wěn)定在5°～7.5°之間，初始攻角為零；彈體離軌后，俯仰角先逐漸增加，在受重力產(chǎn)生的彈體下沉速度影響下，彈體彈道傾角逐漸減小，準(zhǔn)攻角快速增加至1.9°左右；此后，由于彈體受到氣動載荷的作用，彈體俯仰角和準(zhǔn)攻角先逐漸減小后增加，在一定范圍內(nèi)成波動變化，彈道傾角也先增加后減小，在一定范圍內(nèi)成波動變化。

整個考察時間范圍內(nèi)，準(zhǔn)攻角與彈道傾角之和等于俯仰角，與理論狀態(tài)一致，表明模型響應(yīng)是合理的。

3.3 不同模型結(jié)果對比

為考察耦合模型對發(fā)射過程和初始彈道的作用，這里分別對表2所示的模型進(jìn)行求解計算，并對其典型響應(yīng)進(jìn)行對比分析。表2中，在下文中“有/無彈架間隙作用”記為“有/無間隙”；“有/無氣動載荷作用”記為“有/無氣動”；在有彈架間隙作用的狀態(tài)中，采用1 mm的間隙。

圖5 俯仰角、彈道傾角及準(zhǔn)攻角的變化曲線

表2四種組合模型

Table 2The four types of combined model

影響因素彈架間隙作用氣動載荷作用模型一有有模型二有無模型三無有模型四無無

3.3.1 俯仰角對比

圖6為四種模型下彈體在發(fā)射2 s內(nèi)俯仰角對比圖。從圖6可看出，四種模型影響下的俯仰角變化趨勢一致，在發(fā)動機點火初期即在0.041 4 s之前，導(dǎo)彈在導(dǎo)軌上滑行，俯仰角會有波動，但在離軌后，俯仰角逐漸增加后減小呈波動式變化；彈架間隙產(chǎn)生的擾動對彈體的飛行俯仰角影響較大，在有間隙作用的影響下，俯仰角離軌后的增加趨勢變小，且波動變化的范圍在6°～7.5°之間，而無間隙作用的影響下，俯仰角離軌后的增加趨勢較大，波動變化的范圍在6.2°～11.5°之間；氣動載荷的作用對彈體的俯仰角影響較小。

由此可見，彈架間隙擾動作用對彈體俯仰角影響較大，存在彈架間隙的擾動作用會使得俯仰角變化幅值范圍減小，而氣動載荷的影響對彈體俯仰角的影響較小。

3.3.2 彈道傾角結(jié)果對比

圖7為四種模型下彈體在發(fā)射2 s內(nèi)彈道傾角的對比圖。由圖7可知，在發(fā)射0.041 4 s彈體離軌前，彈體在軌道上滑行，四種模型下的彈道傾角都會產(chǎn)生波動，出現(xiàn)一定的下降變化；在彈體離軌后，推力發(fā)動機關(guān)機之前(0.5 s)，在無彈架間隙擾動作用的影響下，彈體的彈道傾角增加較快，并在氣動載荷作用的影響下，彈道傾角增加的幅度更快；在0.5 s時刻，推力撤去，在無氣動載荷作用的影響下，彈道傾角出現(xiàn)明顯的拐點；在0.5 s之后，彈體在氣動載荷作用影響下，彈道傾角成波動式變化，而在無氣動載荷的影響作用下，彈道傾角出現(xiàn)較為明顯的下降趨勢。

由此可見，彈架間隙擾動與氣動載荷作用對彈體彈道傾角影響較大，存在彈架間隙擾動作用，會使得彈道傾角變化幅值范圍減小，而存在氣動載荷作用，會使得彈道傾角成波動式變化。

圖6 俯仰角對比曲線

圖7 彈道傾角對比曲線

3.3.3Y向位移量結(jié)果對比

圖8為四種模型下彈體在發(fā)射2 s內(nèi)Y向位移量對比圖。由圖8可知，四種模型下的Y向位移量變化趨勢基本一致，在彈體離軌后，都成逐漸增加的變化，但變化增加速率不同。在有間隙擾動而無氣動載荷作用的影響下，彈體的Y向位移量最低，且在0.5 s后無推力作用下，Y向位移量增加速率緩慢，在1.5 s左右后會有下降的趨勢；在無間隙擾動而有氣動載荷作用的影響下，彈體的Y向位移量最高，且在0.5 s后無推力作用下，Y向位移量增加速率成波動式變化，Y向位移量持續(xù)增加；在有間隙擾動和有氣動載荷作用的影響下，彈體的Y向位移量增加速率基本保持一致變化，且在0.5 s后無推力作用下，Y向位移量持續(xù)平穩(wěn)增加；在無間隙擾動和無氣動載荷作用的影響下，彈體的Y向位移量增加速率成逐漸減小的變化趨勢，且在1.5 s后無推力作用下，Y向位移量增加緩慢。

圖8 Y向位移量對比曲線

由此可見，彈架間隙擾動和氣動載荷作用對彈體Y向位移量都有影響，存在彈架間隙擾動作用，會使得彈體Y向位移量增加速率變大，而存在氣動載荷作用，會使得彈體Y向位移量在撤去推力后持續(xù)的平穩(wěn)增加。

4 結(jié)論

(1)采用此方法能夠?qū)椉馨l(fā)射初始擾動與初始彈道進(jìn)行相互耦合，且能有效模擬滾轉(zhuǎn)彈的發(fā)射狀態(tài)。

(2)發(fā)射約束段的初始擾動作用，主要體現(xiàn)在彈架間隙引起的隨機擾動作用方面，當(dāng)彈架間隙為1 mm時，此擾動會使得彈體在飛行過程中的俯仰角和彈道傾角的幅值范圍減小4°左右，也使得彈體在飛行過程中的Y向位移量在1.5 s時刻減小6 m左右。

(3)氣動載荷擾動作用，主要與彈體的速度相關(guān)，在彈體離軌后影響較大，當(dāng)離軌速度為9.7 m/s時，存在此擾動會使得彈體在0.5 s撤去推力后的姿態(tài)角成波動式變化，使得滾轉(zhuǎn)彈穩(wěn)定飛行，也會使得彈體Y向位移量在撤去推力后持續(xù)的平穩(wěn)增加。

文中采用的方法可為研究發(fā)射擾動與飛行彈道的影響提供有效的途徑，也可為同類導(dǎo)彈初始擾動的分析提供參考和借鑒。