飛行方案設(shè)計大作業(yè)(1)

1、實用文檔航天飛行動力學(xué)大作業(yè)韓謹陽20213004641、方案飛行2、彈道設(shè)計3、衛(wèi)星攝動與機動第一局部飛行方案第三局部衛(wèi)星的攝動與機動第二局部彈道設(shè)計飛行方案大作業(yè)問題描述在導(dǎo)彈質(zhì)量、轉(zhuǎn)動慣量、發(fā)動機推力等參數(shù)的情況下,導(dǎo)彈分為三個飛行方案,即三個階段飛行.階段一：飛行距離在x9100m,采用追蹤法,其中方案高度與距離的關(guān)系、方案彈道傾角與高度的關(guān)系如下：*(1)H2000cos(0.0003141.1x)5000*、O八.o*、k(H-H)+&(H-困)階段二：飛行距離在24000mx9100m,采用追蹤法,其中方案高度與距離的關(guān)系、方案彈道傾角與高度的關(guān)系、導(dǎo)彈因燃料消耗而質(zhì)量改變參數(shù)如

2、下：3050mk(H-H*)+&由ms0.46kg/s(3)階段三:飛行方案x24000m&ay0,而最終目標位置為xm30000m采用比例導(dǎo)引法rdqVmsinVTsindtmT.yTymtanqxTxm*(4)八kdqdtdt*.,、0k(qq)zk(-)+k&(&)要求:1計算縱向理想彈道,給出采用瞬時平衡假設(shè)mzmzzz0時所有縱向參數(shù)隨時間的變化曲線.2不考慮氣動力下洗影響,計算飛行器沿理想彈道飛行時,你認為可以作為特性點的5個以上點處的縱向短周期擾動運動的動力系數(shù),并分析其在特性點處的自由擾動的穩(wěn)定性,以及計算在各個特性點處彈體傳遞函數(shù)Ws,W1ys,Ws建立模型基于“瞬時平衡假設(shè)

3、,導(dǎo)彈在鉛垂平面內(nèi)運動的質(zhì)心運動方程組為:dVmdtPcosbmVdtPsinbYbmgsinmgcosdxdT(5)VcosdydTVsin由于階段一不考慮導(dǎo)彈質(zhì)量隨時間的變化,其中攻角可根據(jù)瞬時平衡假設(shè)因此階段一的模型需要聯(lián)立公式1、公式5;從而可得到導(dǎo)彈攻角與彈道傾角之間的關(guān)系mrzm：(6)其中Cxq&efXCyqGef(7)其中假設(shè)公式(1)的zk(-*)+k&(&)中的k=-9,k&=-0.5；又由于階段二需要考慮導(dǎo)彈質(zhì)量隨時間的變化,因此階段二的模型需要聯(lián)立公式2公式5、公式6、公式7最后一階段,由于利用了比例導(dǎo)引法公式4的k=2,可得導(dǎo)彈到達目標的相對微分方程為drdTdqrd

4、tVcosVsin而導(dǎo)引率ddtk四、其中k=2;由于第三階段的初始參數(shù)及終點坐標均為直角坐標系,由下列圖可知將30000rcosqyrsinq代入到公式4,得到直角坐標系下的微分方程組dxdTdydt另外補充方程法向平衡方程：drcosdtdr.sindt一dqrsinq-dtdqrsinq-dtmV-dtPsinYmgcos三、算法實現(xiàn)編程使用MATLA瞅件,并運用歐拉方程解微分方程,即ode45函數(shù)；四、程序源代碼*階段一*functiondy=jieduan1(t,y)dy=zeros(4,1);m=320;g=9.8;P=2000;q=0.5*1.2495*(288.15-0.006

5、5*y(4)/288.15).A4.2558*y(1).A2;k=-9;dk=-0.5;Hi=2000*cos(0.000314*1.1*y(3)+5000;dHi=-2000*0.000314*1.1*sin(y(3);delta=k*(y(4)-Hi)+dk*(dy(3)-dHi);alpha=0.34*delta;Xb=(0.2+0.005*alphaA2)*q*0.45;Yb=(0.25*alpha+0.05*delta)*q*0.45;dy=zeros(4,1);dy(1)=P*cos(alpha)/m-Xb/m-g*sin(y(2);dy(2)=P*sin(alpha)/m/y(1

6、)+Yb/m/y(1)-g*cos(y(2)/y(1);dy(3)=y(1)*cos(y(2);dy(4)=y(1)*sin(y(2);end*階段二*functiondy=jieduan2(t,y)dy=zeros(4,1);m=320-0.46*t;g=9.8;P=2000;q=0.5*1.2495*(288.15-0.0065*y(4)/288.15),A4.2558*y(1).A2;k=-0.25;Hi=3050;delta=k*(y(4)-Hi);alpha=0.34*delta;Xb=(0.2+0.005*alphaA2)*q*0.45;Yb=(0.25*alpha+0.05*de

7、lta)*q*0.45;dy(1)=P*cos(alpha/180*pi)/m-Xb/m-g*sin(y(2)/180*pi);dy(2)=P*sin(alpha/180*pi)/m/y(1)+Yb/m/y(1)-g*cos(y(2)/180*pi)/y(1);dy(3)=y(1)*cos(y(2)/180*pi);dy(4)=y(1)*sin(y(2)/180*pi);end*階段三*functiondy=jieduan3(t,y)v=y(4);k=10;m=285.04-0.46*t;q0=-atan(3050/6000);g=9.8;q1=0.5*1.2495*(288.15-0.006

8、5*y(2)/288.15).A4.2558*y(4).A2;k1=10;dk1=0.05;dy=zeros(4,1);r=sqrt(y(1)A2+y(2)A2);q=atan(y(2)/(y(1)-30000);elta=q-y(3);dr=-v*cos(elta);tht=q0+k*(q-q0);dq=v/r*sin(elta);dtht=k*dq;delta=k1*(y(3)-tht)+dk1*(dy(3)-dtht);alpha=0.34*delta;dy(1)=-dr*cos(q)+r*sin(q)*dq;dy(2)=-dr*sin(q)-r*cos(q)*dq;Yb=(0.25*a

9、lpha+0.05*delta)*q1*0.45;dy(3)=(2000*sin(alpha)/m+Yb/m-g*cos(y(3)/v;y(4)=v;endmain*導(dǎo)彈質(zhì)量m(1)=287.2204;%P=2000;%發(fā)動機推力g=9.8;k=5;det(1)=0.045;a(1)=0.6186;sit(1)=-0.000002024;V(1)=217.2867;x(1)=24000;H(1)=3071;%H1(1)=3050;S=0.45;%L=2.5;%k1=-0.14;k2=-0.06;sit1(1)=sit(1);p0=1.2495;%初始速度%初始位置初始高度參考面積參考長度T0=

10、288.15;T(1)=T0-0.0065*H(1);p(1)=p0*(T(1)/T0)A4.25588;q(1)=1/2*p(1)*V(1)A2;%Cx(1)=0.2+0.005*a(1)A2;Cy(1)=0.25*a(1)+0.05*det(1)*180/pi;大氣密度計算公式%升力系數(shù)Y(1)=Cy(1)*q(1)*S;X(1)=Cx(1)*q(1)*S;SIT(1)=(P*sind(a(1)+(Y(1)-m(1)*g*cos(sit(1)/m(1)/V(1);Q(1)=atan(-H(1)/(30000-x(1)+pi;r(1)=6708.2039;R(1)=-V(1)*cos(Q(1

11、);n(1)=Q(1)+pi;SIT1(1)=k/r(1)*(V(1)*sin(n(1);mza=-0.1;%mzdet=0.024;%t=0;i=0;dt=0.01;ms=0.46;%whileH0&H10%i=i+1;t=t+dt;俯仰力矩系數(shù)對攻角的偏導(dǎo)數(shù)俯仰力矩系數(shù)對舵偏角的偏導(dǎo)數(shù)質(zhì)量秒消耗量運用迭代法求解det(i+1)=k1*(sit(i)-sit1(i)+k2*(SIT(i)-SIT1(i);a(i+1)=-mzdet/mza*det(i)*180/pi;Cy(i+1)=0.25*a(i)+0.05*det(i)*180/pi;Cx(i+1)=0.2+0.005*a(i)A2;Y

12、(i+1)=Cy(i)*q(i)*S;X(i+1)=Cx(i)*q(i)*S;m(i+1)=m(i)-ms*dt;sit(i+1)=sit(i)+(P*sind(a(i)+(Y(i)-m(i)*g*cos(sit(i)/m(i)/V(i)*dt;V(i+1)=V(i)+(P*cosd(a(i)-(X(i)+m(i)*g*sin(sit(i)/m(i)*dt;x(i+1)=x(i)+V(i)*cos(sit(i)*dt;H(i+1)=H(i)+V(i)*sin(sit(i)*dt;Q(i+1)=atan(-H(i)/(30000-x(i)+pi;sit1(i+1)=k*(Q(i)-Q(1);H1

13、(i+1)=H(i)+V(i)*sin(sit1(i);SIT(i+1)=(sit(i+1)-sit(i)/dt;r(i+1)=(H(i)A2+(30000-x(i)A2)A(1/2);R(i+1)=(r(i+1)-r(i)/dt;n(i+1)=acos(-R(i)/V(i)+pi;SIT1(i+1)=k/r(i)*(V(i)*sin(n(i);T(i+1)=T0-0.0065*H(i+1);p(i+1)=p0*(T(i+1)/T0)A4.25588;q(i+1)=1/2*p(i+1)*V(i+1)A2;endplot(x,H);holdont,y=ode45(jieduan1,039.056

14、4,250007000);plot(y(:,3),y(:,4);holdont,y=ode45(jieduan2,39.0564115,192.768-0.00991002998.71);plot(y(:,3),y(:,4);其中每一段的初始值,均為上階段的結(jié)束值所以每一階段計算結(jié)束后,需要再給出所有數(shù)據(jù)的結(jié)果,找到每一段距離相對應(yīng)的數(shù)據(jù),即為初始值.五、結(jié)果分析制出導(dǎo)彈三個階段的飛行軌跡如圖1圖1圖2是第一階段縱向參數(shù)隨時間的變化曲線;/flmnhnyM,1,-imV隱喪it由繳.器6隨談化曲線-e湍t變化曲線xT-曲步證曲線一gj；%、.丁4230n二.,卜、=*,-39-口圖2圖3時第二

15、階段縱向飛行參數(shù)隨時間的變化曲線由圖1導(dǎo)彈在第一階段,從初始高度7000m,開始下降飛行,在距離9100m時,開始變?yōu)榈歉唢w行,距離到達24000m至目標30000m這一階段為導(dǎo)彈的下降尋找目標階段；由圖2得,第二階段的飛行速度先增加后減小,在第一階段末尾階段速度減小至192.768m/s;彈道傾角先減小后增加,海拔高度隨時間的增加而減小;JM)一施凌枇曲嵯*M/于可仆J!_T:雁t變化曲葡,0度和3050m之由圖3得,第三階段為登高飛行,所以彈道傾角和海拔高度分別在間振蕩,而速度也根本在140m/s至150m/s之間徘徊；六、特性點的動力系數(shù)、傳函分別取特性點1:x=0時；特性點2:x=91

16、00時;特性點3:x=24000時;特性點4:x=30000時由縱向自由擾動的穩(wěn)定性條件a,4a22a340即縱向自由擾動運動穩(wěn)定.根據(jù)以下公式:M0=mqSL0(MZ)0a22a24a25a34a35(mzqSL)0(W)0Jz0(Ma).JZ0(q)0JZ0(PYv)(mV)0(Yz)0(mV)0得到以下值:特性點170000.6012250-2113.912507.33特性點291000.9273192-3676.815882.4357特性點3240000.92311206-1974.445473.86694特性點4300001.2495195-2672.563641.4152特性點10

17、.006349216.7108-1.61010.051420.0052848-0.01760.3860.0115特性點20.00634911.672-2.80140.070960.009191-0.026220.29270.01131特性點30.0063496.26808-1.504340.049680.004936-0.016680.399410.011189特性點40.006349218.48431-2.0360.0735810.0084-0.026080.34330.01372特性點1的傳遞函數(shù):W(s)=-10.017625(0.386013a1)s(0.149s20.0086s1)W(s)1wy(s)0.017625(s0.0063492)s(0.149s20.0086s1)0.4496(10.0012667(s0.006349)(0.149s2_0.00861)特征點2的傳遞函數(shù):W(s)=-10.02622(0.2927s1)s(0.0857s20.0066s1)W(s)1w/s)0.02622(1-0.0026899(s0.0063492)Z2s(0.0857s0.0066s1)0.71035(10.0022667(s0.006349)(0.857s20.0061)特性點3的傳遞函數(shù):W(s)=-10.016