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

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

0 引言

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

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

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

1 耦合計算模型

1.1 耦合計算思路

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

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

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

(2)無控飛行段，在前述多體動力學模型基礎上，將彈體承受的氣動載荷轉換至彈體坐標系，利用動力學方法進行彈道計算。

圖1 彈體飛行三階段

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

1.2 多體動力學模型

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

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

(1)

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

(2)

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

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

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

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

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

(3)

圖2 碰撞接觸模型

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

圖3 導軌與定向件的配合間隙

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

ny=A(H)×R(t)

(4)

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

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

(5)

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

2.2 氣動載荷與彈道解算

本研究以多體動力學的模型為基礎，將彈體受到的氣動載荷均轉換到彈體坐標系下，并利用動力學方法[17]進行彈道計算。

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

表1 氣動載荷參數的定義與求解

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

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

2.2.1 準速度坐標系下的氣動力

在準速度坐標系下，彈體受到的阻力X*，升力Y*以及側向力Z*可以表示為

(6)

2.2.2 準彈體坐標系下的氣動力矩

為求解準彈體坐標系下的力矩，需要用到準彈體坐標系下的姿態(tài)角速度的如下轉換關系式：

(7)

式中ωx1、ωy1、ωz1分別為彈體坐標下的滾轉角速度、偏航角速度和俯仰角速度。

在準彈體坐標系下，彈體受到的滾動力矩Mx4、偏航力矩My4以及俯仰力矩Mz4可以表示為

(8)

2.2.3 彈體坐標系下力與力矩轉換

在彈體坐標系下，可利用如下關系式計算其受到的氣動力：

(9)

其中

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

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

(cosα*sinβ*)Z*

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

(sinα*sinβ*)Z*

在彈體坐標系下，彈體承受的力矩可由準彈體坐標系下的力矩的關系式獲得：

(10)

3 基于滾轉彈的實例分析

3.1 實例模型及參數

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

圖4 多體動力學簡化模型

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

3.2 耦合模型典型響應

利用商用多體動力學軟件ADAMS對耦合實例模型進行求解計算，獲得了耦合模型作用下彈體典型的響應如下所述。

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

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

3.3 不同模型結果對比

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

圖5 俯仰角、彈道傾角及準攻角的變化曲線

表2四種組合模型

Table 2The four types of combined model

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

3.3.1 俯仰角對比

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

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

3.3.2 彈道傾角結果對比

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

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

圖6 俯仰角對比曲線

圖7 彈道傾角對比曲線

3.3.3Y向位移量結果對比

圖8為四種模型下彈體在發(fā)射2 s內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 結論

(1)采用此方法能夠將彈架發(fā)射初始擾動與初始彈道進行相互耦合，且能有效模擬滾轉彈的發(fā)射狀態(tài)。

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

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

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