小行星轨道问题解

小行星轨道问题学号：专业：班级：姓名：摘要：要确定一颗小行星绕太阳运行的轨道，需要在轨道平面内建立以太阳为原点的空间直角坐标系，然后在不同时刻对小行星进行观测，以确定其轨道。由于观测所得到的数据较多，但毫无疑问是存在误差的。因此，应该尽量使用到这些数据来求最优解。关键字：椭圆，最优解，最小二乘法^MatlabLingo问题描述：要确定一颗小行星绕太阳运行的轨道，需要在轨道平面内建立以太阳为原点的空间直角坐标系，然后在不同时刻对小行星进行观测，以确定其轨道。已知在5个不同时刻对某颗小行星进行了5次观测，表B给出了相应的观测数据。表1：某小行星的5次观测数据（单位：天文单位）12345X坐标5.7646.2866.7597.1687.480Y坐标0.6481.2021.8322.5263.360注：一个天文单位等于地球到太阳的平均距离，即L4959787x1011米。你所要作的工作是：确定这颗小行星的轨道，如椭圆的半长轴、半短轴、半焦距、近日点、远日点，以及椭圆轨道的周长等。题解：由开普勒第一定律：太阳系中的所有行星围绕太阳运动的轨道都是椭圆，太阳处在所有椭圆的一个焦点上。可知：此小行星的轨道为椭圆，太阳是这个椭圆的一个焦点。下面以此为前提来求解小行星的轨道方程。由于题目提供了5组观测数据，这里有两种可行的方法：1、仅考虑椭圆上的这5个点，不考虑太阳的位置（不一定是焦点），设椭圆轨道方程为一般方程ax2+2ax^y+ay2+2ax+2ay+1=012345然后代入5个点的坐标，使用matlab算出结果。2、在考虑太阳的位置（必须是焦点）基础上，对另外五个点用最小二乘法的方法进行曲线拟合，求出最优解，从而得出方程，即认为这五个点都有一定的误差，拟合后得出最优解。

详细步骤：1、第一种方法需要确定系数a，顷，2,3,4,5；利用已知的数据，不妨设qP顷，2,3,4,5;欲确定系数七•等价于求解一个线性方程组:TOC\o"1-5"\h\zax2+2axy+ay2+2ax+2ay+1=011211314151ax2+2axy+ay2+2ax+2ay+1=012222324252*ax2+2axy+ay2+2ax+2ay+1=013233334353ax2+2axy+ay2+2ax+2ay+1=014244344454ax2+2axy+ay2+2ax+2ay+1=0"1525535455533.22377.47010.419911.52801.296039.513815.11151.444833.22377.47010.419911.52801.296039.513815.11151.444812.57202.404045.684124.64333.323313.51803.646051.380236.21276.380714.33605.052054.878549.781811.289614.81606.7200即:「a]「—「1a-1a=-13a.—1a4-1L5」1——1x212xy11y212x12y1「a]1「-「x222x2y2y222x22y2a2TA=x232x3y3y232x32y3X=a3b=Tx22x4y4y242x42y4a4T其中，x2L52x5y5y252x52y5」，aL5」，_-1_可见，解答上述问题就是对线性方程进行求解。模型求解：1.1利用matlab求解矩阵代码如下：A=[33.2237,7.4701,0.4199,11.5280,1.2960;39.5138,15.1115,1.4448,12.5720,2.4040;45.6841,24.6433,3.3233,13.5180,3.6460;51.3802,36.2127,6.3807,14.3360,5.0520;54.8785,49.7818,11.2896,14.8160,6.7200]B=[-1;-1;-1;-1;-1]inv(A)*B结果：A=33.22377.47010.419911.52801.296039.513815.11151.444812.57202.404045.684124.64333.323313.51803.646051.380236.21276.380714.33605.052054.878549.781811.289614.81606.7200B=-1-1-1-1-1ans=0.0507-0.03510.0381-0.22650.1321所以，轨道方程是：0.050力-2*0.035必+0.038^2-2*0.2265+2*0.132^+1=01.2画出轨迹图：代码如下：>>x=-3.5:0.01:15.7;>>y=(-24.92+5.66*x+sqrt((5.66*x-24.92).A2-4*2.58*(4.23*x<2-41.02*x+1)))/5.16;>>y1=(-24.92+5.66*x-sqrt((5.66*x-24.92).A2-4*2.58*(4.23*x.A2-41.02*x+1)))/5.16;>>plot(x,y,x,y1),gridon结果：显然，太阳（坐标（0,0））不在椭圆的焦点，而在小行星的轨道附近，这是不符合实际的。因此，使用这种方法得到的结果，误差较大，舍弃。

2.第二种方法根据椭圆上的点到两焦点的距离和为定值，设另一个焦点的坐标为F(m,n),半长轴为a。使用最小二乘法求最优解。列出式子如下：i22:~222+\Xi+Yi-2a]Minz[«(Xi-m)+(Yi+\Xi+Yi-2a]model:sets:VOR/1..5/:x,y;Endsetsdata:x,y=5.7640.6486.2861.2026.7591.8327.1682.5267.483.36;enddatamin=@sum(VOR:@sqr(@sqrt(xA2+yA2end结果为：Localoptimalsolutionfound.Objectivevalue:Objectivebound:Infeasibilities:Extendedsolversteps:Totalsolveriterations:VariableCostM0.000000N-0.1280014E-08A0.2761780E-08X(1))+@sqrt((x-m)A2+(y-n)A2)-2*a));0.3372721E-040.0000000.0000001493708552)+@sqrt((x-m)A2+(y-n)A2)-2*a));0.3372721E-040.0000000.0000001493708552ValueReduced6.3880464.8504135.0241995.764000X(2)X(3)6.7590000.0000000.000000X(4)7.1680000.0000000.000000X(5)7.4800000.0000000.000000Y(1)0.64800000.0000000.000000Y(2)1.2020000.0000000.000000Y(3)1.8320000.0000000.000000Y(4)2.5260000.0000000.000000Y(5)3.3600000.0000000.000000RowSlackorSurplusDualPrice10.3372721E-04-1.000000-1.000000即M=6.388046N=4.850413A=5.024199另一焦点坐标为F(6.388046,4.850413)，*左rh#-rCnoai半长轴为5.024199.所以，小行星的轨迹方程为:22(X-6.388)+(Y-4.850)+:X+Y=10.048轨迹图:

经过简单的计算可得到以下结果：半长轴a为：5.024119个天文单位；半短轴b为：2.972619个天文单位；半焦距c为：