空间碎片轨道预报中动力学模型基准的选取.pdf

第6 4卷第4期2023年7 月doi:10.15940/ki.0001-5245.2023.04.005天文学报ACTA ASTRONOMICA SINICAVol.64 No.4Jul.,2023空间碎片轨道预报中动力学模型基准的选取胡松华1,2,3汤靖师1,2,3十(1南京大学天文与空间科学学院南京2 10 0 2 3)(2现代天文与天体物理教育部重点实验室南京2 10 0 2 3)(3南京大学空间环境与航天动力学研究所南京2 10 0 2 3)摘要在不同的轨道预报场景中，使用的动力学模型也不同。例如，在低轨空间碎片的预报中大气阻力是十分重要的摄动力，而到了中高轨，大气阻力就可以忽略不计.如何为不同轨道类型的空间碎片选择最优(满足精度要求下的最简)动力学模型还没有系统、详尽的研究.将对不同精度需求、不同轨道类型下的大批量轨道进行预报，通过比较不同动力学模型下的预报结果，给出各种预报场景的最优动力学模型建议.可以为不同轨道类型的空间碎片在轨道预报时选择基准动力学模型提供参考或标准.关键词天体力学，地球，方法：数值，行星和卫星：动力学演化和稳定性中图分类号：P135；文献标识码：A1引言截至2 0 2 1年，太空中有超过36 0 0 0 个尺寸大于10cm的空间碎片 1，它们包括废弃的卫星、火箭残骸等，而撞击事件、自爆事件更会使空间碎片的数目爆炸式地增长大量聚集的空间碎片会威胁在轨航天器的安全，并且严重影响天文学观测 1太空中的碰撞事件已发生多次，比如2 0 0 9 年俄罗斯的退役卫星宇宙2 2 31号和美国的商业通讯卫星铱星33发生了碰撞事故，二者在碰撞后解体成为30 0 0 多个空间碎片 2 因为空间碎片数量巨大,且分布十分广泛，所以在对它们进行轨道预报的过程中，如何选取适当的动力学模型至关重要，精细的摄动模型可以带来精度更高的结果，但随之而来的是更大的计算量需求.对于不同高度、不同偏心率的空间碎片，使用的动力学模型差别很大.图1展示了典型面质比的空间碎片在不同高度受到的摄2022-04-18收到原稿，2 0 2 2-0 6-0 6 收到修改稿*国家自然科学基金项目(118 7 30 31)资助动力大小，从图中我们可以看到，随着轨道高度的变化，各项摄动力的量级有着巨大的变化，起主导作用的摄动力也在不断变化.在保证预报精度的前提下，尽可能地简化动力学模型可以很大地提高大批量空间碎片轨道预报的效率，如何在保证精度和保证速度之间取得平衡正是这篇文章研究的重点.北美防空司令部(The North American Aero-space Defense Command,NORAD)通过两行根数(Two-LineElement Set，T LE)编目空间目标,并使用SGP4(Simplified General Perturbations4)/SDP4(Simplified Deep-Space Perturbations4)算法对目标库进行应用.这种算法是一种基于正则变换的轨道分析解，由于SGP4和SDP4适用的轨道不同，它们考虑的动力学模型也不同，前者考虑的是地球非球形引力位的低阶近似以及符合幂律分布的大气阻力，而后者考虑的是地球的非球形引42-164卷力位的低阶近似和日月引力 3然而这种算法使用的动力学模型仅能保证特定的精度(在一周内的轨道预报可以保证位置误差在10 km量级 4-5)，并且依赖于其发布的TLE6.针对空间碎片，需要提出一种针对不同预报场景(即不同轨道类型，不同精度需求)的动力学模型.102GM10J2102J33104J3106MoQn.sun1081010102图1不同高度的空间碎片受到的摄动力.其中GM表示地球的二体引力，J或J*表示地球的非球形引力摄动的主要带谐项和主要田谐项，计算大气阻力和光压摄动的相关参数见2.1节.Fig.1 The perturbation force of space debris at differentheights.Among them,GM represents the two-bodygravitation of the earth,J or J*represents the main bandand tesseral harmonic terms of the aspheric gravitationalperturbation of the earth.The relevant parameters forcalculating the atmospheric drag and solar radiationpressure perturbation are listed in Section 2.1.我们使用不同的动力学模型对大量的空间碎片轨道进行轨道预报，通过比较不同模型下轨道预报的结果，以精细动力学模型下的结果为参考，给出不同的动力学模型下轨道预报的精度.参照预报结果的精度，我们可以针对不同的预报场景，给出最优化的动力学模型.本文第2 节将会描述轨道预报的过程，第3节会展示不同预报场景和不同动力学模型下的预报精度，第4节通过分析第3节的结果给出动力学模型建议，第5节是总结.天文学报2轨道预报方法2.1空间碎片的动力学空间碎片在地球的引力下运动，在理想情况下为一个二体运动，它是一个可积系统，存在解析解.而在现实中，地球质量和形状的不均匀、太空中稀薄大气的存在等因素使得空间碎片绕地球运动成为一个受摄动的二体问题.我们研究的摄动力主要包括地球的非球形引力、大气阻力、日月的第三体引力、太阳辐射压力(Solar RadiationPressure,SRP)、地球固体潮汐力、后牛顿效应和太阳系其他大天体的第三体引力.本文中使用的地球引力场模型是GGM03C(Global Gravity Model 03C)7,它是由GRACE(Gravity Recovery and Climate Experiment)卫星SolarRadiationPressure数据反演的引力场之一，最高支持36 0 36 0 阶次.Atmosphere在计算第三体引力摄动时，计算日月位置使用的是JPL(Jet Propulsion Laboratory)数值历表和一种103104Height/km4期平均化的分析历表，我们还会在下文中对二者做一个简单的比较.在大气阻力摄动的计算中，相对麻烦的是大气密度的计算，因为大气受太阳辐射强度和地磁活动等影响较大，大气密度模型的精度有限 8 本文中使用的是一种参考了NRLMSISE-00(Naval Research Laboratory Mass Spectrometerand Incoherent Scatter Radar Exosphere 2000)9大气模型的扁球指数平均大气模型 8.大气阻力和光压摄动的强度受面质比等因素的影响，在没有特别说明的情况下，本文采用表1所示参数.表1中列举的4项参数分别为大气阻力系数、迎风向等效面质比、太阳方向等效面质比和反射系数，它们将被应用到大气阻力摄动和太阳辐射压摄动的计算中。轨道预报采用数值方法，即使用变步长的Runge-Kutta-Felhberg7（8)积分器积分高斯型摄动运动方程，截断误差小于10-13（长度单位为地球半径，角度单位为弧度)，轨道根数采用第1类无奇点根数,即a,i,=ecosw,n=esinw,入=w+M,其中a为轨道半长径、i为轨道倾角、为升交点经度、e为轨道偏心率、w为近地点辐角、M为平近点角.42-264卷表1大气阻力和太阳辐射压相关参数取值表Table 1 Table of parameters related toatmospheric resistance and solar radiationParametersDrag coefficientWindward effectivesurface-to-mass ratioEffective surface-to-massratio towards the SunReflection coefficient2.2试验轨道的初值选择为了模拟真实的空间碎片分布，轨道初值将会根据TLE数据进行选择，一方面使其符合真实的空间碎片分布(主要是半长径和偏心率的分布)，另一方面需保证均匀且全面地取值如图2 所示，我们画出了目前TLE所有在编空间目标的轨道高度与轨道偏心率分布，并且使用4种不同的颜色代表上述4种轨道类型.我们对空间目标的轨道进行划分，共分为以下4类，这里轨道高度的定义为轨道的半Fig.2 The distribution of TLE objects with orbit height and orbit eccentricity胡松华等：空间碎片轨道预报中动力学模型基准的选取pressureValues2.27.08 10-3 m.kg-17.08 10-3 m.kg-10.51口Type10.9Type2Type3Type40.80.70.60.50.40.30.20.101024期长径与地球赤道半径之差，即轨道近地点与远地点的平均高度.类型一：低轨(Low Earth Orbit,LEO)，轨道高度在2 0 0 0 km以下，其中比较具有代表性的是太阳同步轨道和星链卫星所在轨道.我们研究的低轨目标总数为159 2 1个.类型二：中轨(Medium Earth Orbit,MEO)小偏心率，轨道高度高于2 0 0 0 km，低于36 0 0 0 km，e0.1.我们研究的中轨小偏心率目标总数为58 4个.类型三：中轨大偏心率，轨道高度高于2000km，低于36 0 0 0 km，e 0.1，包括闪电轨道(Molniya).我们研究的中轨大偏心率目标总数为2264个.类型四：高轨(High Earth Orbit，H EO),一部分为地球同步轨道(GeosynchronousOr-bit，G SO)，其中包括倾斜地球同步轨道(InclinedGeosynchronousOrbit，I G SO)、地球静止轨道(Geostationary Orbit，G EO)，另一部分为高度高于GSO的轨道(包含大偏心率轨道和小偏心率轨道)和高度等于GSO的大偏心率轨道.我们研究的高轨目标总数为156 0 个。OO103104Height/km图2 TLE在编目标随轨道高度与轨道偏心率的分布10542-364卷3车轨道预报的误差分析在以下小节中，我们比较的都是一天内的位置偏差.假设待研究的动力学模型为A，基准动力学模型为B.基准动力学模型中的摄动力计算是足够精确的，譬如在研究10 10 非球形引力位对预报精度的影响时，动力学模型A中的引力场阶次为10 10,而基准动力学模型B中的引力场阶次为8 0 8 0.分别使用这两种动力学模型预报同一条轨道，每隔一分钟记录一次，得到两条星历用位置表示两条星历分别为TA1，T A 2，T A 1440 和TB1，T B2，TB1440.定义位置偏差d为-1440IFAi-TB:l2d11440这个偏差就表征了使用动力学模型A在一天内积累的位置误差.3.1引力场阶次对轨道预报精度的影响作为地球形状和质量分布不均匀导致的摄动力，非球形引力摄动的效应会随轨道高度的增加而减弱.另外，考虑的引力场阶次越高，摄动力的计算越精确，但是轨道预报的速度也越慢不同轨道类型对引力场的阶次的敏感程度不同，为了比较不同引力场阶次下的位置偏差，我们使用不同精细程度的引力场对空间碎片进行轨道预报，除了地球的非球形引力摄动外，这里没有考虑别的摄动力.用以计算的初始轨道是2.2 节中提到的全部轨道.不同引力场阶次下的轨道预报误差如图3-4所示，其中图3对应低轨目标，图中每一个点表示相应轨道高度和引力场阶次下的位置误差，位置误差的大小用颜色表示线条为特定精度的分界线，例如黑色虚线表示的是精度为10 m的分界线，分界线上方区域的精度高于10 m，分界线下方区域的精度低于10 m.图4对应中高轨目标，图中的表示方法与图3类似从图3-4中可以看出，地球非球形引力摄动的影响随着轨道高度的增加而降低，这就意味着预报的轨道越高，采用引力场模型就可以越简单，例如，我们希望轨道预报的结果在一天内的位置误差不超过10 0 m，如果预报的对象是一个高度为30 0 km的低轨空间碎片，我们需要至少采用30 30 阶的引力场相同精度要求下，如果这是天文学报一个处于中轨的小偏心率轨道的空间碎片，那么引力场模型就可以简化到10 10 阶，对于地球同步轨道的空间碎片，更是可以简化到55阶.对于低轨空间碎片的轨道预报，引力场的准确性要求较高，随着引力场阶次的提高，轨道预报位置偏差的减小十分缓慢，当引力场采用8 0 8 0 阶时，绝大部分空间碎片的轨道预报位置偏差减小到了1m以下，但这样高阶次引力场的计算量已经十分大了，再提高轨道预报的精度十分困难.当轨道较低(例如近地点低于40 0 km的轨道)时，随着引力场阶次的提高，精度提升十分缓慢.如果我们希望这些空间碎片的轨道预报达到厘米级精度，这将会是十分困难的.另外，这样的精度要求也并不现实，(1)74期因为在低轨空间碎片的轨道预报中，大气阻力摄动的影响显著，且精确计算十分困难，所以由大气模型造成的误差一般远大于1cm.因此，对于轨道高度低于40 0 km的轨道，将引力场的精度定在厘米或以下意义不大.事实上，大多数低轨空间碎片不会处在极低的轨道上，随着轨道高度的增加，对引力场精确度的依赖逐渐在减小.同样很难使精度提高到厘米级的是中轨大偏心率轨道，其中一些轨道的近地点高度甚至低于2 0 0 km，所以对他们的处理与低轨类似.使用不同的引力场阶次需要的算力有很大区别，因此选取合适的引力场是节省轨道预报时间的关键如图5所示我们展示了使用不同精细程度的引力场预报一天轨道所需的时间.使用不同的硬件或软件环境时，所需时间会有所差异，我们这里使用的CPU是intel core i5-9500，使用的编译器是in-tel visual fortran，使用release方式编译.图中展示的4条轨道的初始开普勒根数如表2 所示.图5中纵坐标表示轨道预报一天使用的CPU时间，取10 次测试的平均值一般来说轨道越高，轨道预报的速度越快，这是因为轨道高度越高，变步长积分器可以取到的步长越长，同样的预报时间内步数可以越少.另一方面，对于同一高度的轨道来说，动力学模型越简单，轨道预报速度越快.从图5中的LEO轨道来看，如果将引力场阶次从10 10 提升到50 50，轨道预报的用时增加到了原来的5倍左右.因此，如果结合图3与图4给出某种预报场景下的最佳引力场阶次，就可以大大地节省算力.42-464卷胡松华等：空间碎片轨道预报中动力学模型基准的选取4期80807070606050504040103Deviation=100mDeviation=10mDeviation=1mDeviation=0.1mDeviation=0.01m10210130302020101020图3低轨空间碎片在不同阶次的引力场下预报结果的位置误差.横坐标为轨道的平均高度，纵坐标为引力场的阶次，其中2 0 表示仅考虑J2项，线Fig.3 The position error of the LEO space debris projected in different orders of gravitational field.The abscissa is the averageheight of the orbit,and the ordinate is the order of the gravitational field,where 2x0 means that only J2 is considered,and the10-110210-3400600line is the dividing line of a certain position error.800Height/km条为某一个位置误差的分界线.1000120014002525Deviation=100m-Deviation=10mDeviation=1 m2020Deviation=0.1mDeviation=0.01 m1515101055202图4中高轨空间碎片在不同阶次的引力场下预报结果的位置误差.横坐标为轨道的平均高度，纵坐标为引力场的阶次，其中2 0 表示仅考虑J2项，Fig.4 The position error of the projection results of MEO/HEO space debris under different orders of gravitational field.Theabscissa is the average height of the orbit,and the ordinate is the order of the gravitational field,where 2x0 means that only J2345678910Height(10km)线条为某一个位置误差的分界线.is considered,and the line is the dividing line of a certain position error.203042-564卷天文学报6LEO-MEO.Molniya54期GEO4320图5轨道预报用时随引力场阶次的变化图.当横坐标为0 时，表示二体引力，当横坐标为n时，表示nn阶次引力场图中Molniya表示闪电轨道，Fig.5 The orbit prediction time as a function of the order and the degree of the gravitational field.When degree or order is 0,it represents two body gravity.When degree or order is n,it represents the gravitational field of order nxn.Molniya in thefigure represents Molniya orbit.Molniya orbit is a special orbit with large eccentricity,and the orbital period is about 12 h.表2 4条轨道在初始时刻的开普勒根数Table 2 Keplerian elements of the four orbits atthe initial momentOrbita/kmLEO71780.01MEO103780.01Molniya266000.7463GEO423780.013.2大气阻力对轨道预报精度的影响本节对低轨试验轨道进行一天的轨道预报，对比是否考虑大气阻力的预报结果，研究动力学模型中忽略大气阻力造成的误差.大气阻力摄动的计算中，大气密度的计算尤其重要现有研究表明，目前的大气密度模型在比较算例中各有优劣，模型之间内符合误差最大可达2 0%-30%10-12，这里用来计算的轨道是类型一的轨道，即低地球轨道，总数16 0 0 条.随着太阳活动强度、地球磁场等环境的变化，大气密度的改变也相当显著，在没有特别说10e2W/M/453045303027003020闪电轨道是一种特殊的大偏心率轨道，轨道周期约为12 h.4050405050405030DegreeorOrder明的情况下，我们取太阳活动指数F107=130，地球磁场指数Ap=10.在这种环境下，高度440 km的大气密度约为1.6 2 7 10-12 kgm-3.在本小节我们也会适当研究不同F107和Ap下的位置偏差的区别轨道预报的结果如图6 所示，这里横坐标取近地点高度，这是因为即使是半长径很大的轨道，如果偏心率较大，那么由于它短暂地穿过较厚重的大气，使得轨道受到大气的影响依旧十分明显，所以使用近地点高度衡量大气阻力的地位更为合理.黑色散点表示在一般情况下(F107=130，A p=10)试验轨道的误差情况，实线是对黑色散点的拟合：这里使用的拟合公式是y=e-，使用这种公式构型的原因是这里的位置偏差与大气密度存在密切的关系，而大气密度随高度的变化正是这种形式，即一方面大气密度随高度指数式下降，另一方面下降的速度会越来越慢(标高越来越高).在一天的轨道预报中，40 0 km左右高度的空间碎片可以积累几公里到几十公里的误差，而当轨道高度低405060708042-664卷于40 0 km时，大气阻力在一天内对预报结果的影响会呈现指数式的上升，甚至可以达到百公里和千公里的量级同时，这部分轨道的高度下降也是十分快的，因此空间碎片或航天器所处的轨道高度一般不会低于30 0 km.而随着轨道稍稍抬高，大气胡松华等：空间碎片轨道预报中动力学模型基准的选取4期阻力的影响开始骤减几百公里高的轨道积累的误差只有几米到几十米之间，对于轨道高度大于1200km的空间碎片，位置误差均小于1m且迅速衰减，到了16 0 0 km，则一般场景下几乎可以不用考虑.口Ap=10,F107=130106Ap=10,F107=130Ap=5,F107=70Ap=15,F107=190SSOor StarLink1044102101200图6 大气阻力对轨道预报精度的影响随轨道高度的变化图方形散点为相应情况下不同轨道的位置偏差的分布情况，实线是对方形散点的拟合虚线和点划线是其他两种情况下的结果的拟合.圆形散点表示的是低轨目标中的一些具有代表性的轨道，太阳同步轨道(Sun-SynchronousOrbit，Fig.6 The effect of atmospheric resistence on the accuracy of orbit prediction as a function of orbital altitude.The squarescatters are the distribution of the positional deviations of different orbits in the corresponding cases,and the solid line is thefitting of the square scatter points.The dashed and dot-dashed lines are the fittings of the results in the other two cases.TheSSO卫星和星链(StarLink)卫星等是在低轨目和地磁活动会影响大气密度，从而间接地影响大气标中的一些典型目标，我们在图6 中标记了这些目阻力的大小，并且影响程度比较显著就轨道高度标.从图中我们可以看出，这些目标相比其他目标为8 0 0 km的空间碎片而言，虚线与点划线的结果有并没有特殊之处，即决定大气阻力造成的误差情况一个量级左右的差别.因此，在考虑大气阻力对轨主要在于轨道的平均高度或近地点高度.道预报结果的影响时，有必要考虑当日的太阳活动我们还在图6 中展示了其他两种太阳活动情况情况和地磁活动情况.和地磁活动情况下大气阻力对轨道预报结果的影响情况，这两种情况分别是F107=70，A p=5和F107=190，A p=15.因为过多的散点会重叠在一起使图像过于混乱，这里没有画出这两种F107和Ap下相应的散点图而仅展示了拟合结果，拟合结果分别用点划线与虚线表示.拟合公式同样是取=e-薛的形式，拟合参数如表3所示.太阳活动400circle scatters represent SSO and starlink satellite orbits.600SSO)和星链(StarLink)卫星轨道.800Perigeeheight/km3.3日月第三体引力对轨道预报精度的影响第三体引力摄动(这里指太阳和月球)是空间碎片动力学中的主要摄动力来源之一，在绝大多数的预报场景中都不可忽略对于任意高度的空间碎片，太阳和月球的引力不仅在同一个量级，而且大小接近，那么我们就不需要对他们分别研究，而1000120014001600180042-764卷是考虑它们的共同作用.第三体引力摄动的影响会随着空间碎片轨道高度的增加而增加，对于低轨空间碎片而言，使用不考虑第三体引力摄动的动力学模型预报一天的弧段会产生百米级的位置误差，而对于地球同步轨道空间碎片而言，位置误差会达到几十公里.第三体引力摄动的计算并不复杂，关键在于获得准确的太阳、月球位置我们将会使用两种方法来获得太阳和月球位置，分别是分析历表 13和数值历表(即JPL历表 14).在这一小节我们需要研究两个问题，其一是分析历表和数值历表的精度对轨道预报的影响有多大？其二是如果使用不同版本的JPL历表，它们的差别有多大？分析历表中的太阳位置可以使用太阳轨道的平根数获得，而月球位置的周期项较大，不可以简单地使用平根数获得，在本文中使用的是带有部分周期项的解析轨道 13，相对精度(即位置误差与位置大小之比)约为|rl/r=O(10-3).分析历表不仅计算速度快,而且不依赖于大体量的拟合参数数据文件，可以减小程序的内存占用，提高运行速度相对于JPL历表，分析历表精度较低，但依然能满足大多数轨道CircularMEo天文学报预报场景的需求.表33种结果的拟合参数Table 3 Fitting parameters for the threeoutcomesF107Ap130107404-4.323201.47051901510930-5.817423.8图7 展示了使用日月分析历表带来的位置误差，即相对于使用精度十分高的数值历表结果的位置偏差.这里计算的轨道分为4部分，分别是低地球轨道、中轨小偏心率轨道、闪电轨道和地球同步轨道.显然，轨道高度越高的轨道受到第三体引力摄动的影响越大，分析历表带来的误差也越大.对于处于低轨的空间碎片，一天内积累的位置误差在1m左右，对于处于中轨的空间碎片，一天内积累的位置误差在10 m量级，对于闪电轨道和地球同步轨道的空间碎片，位置误差则在百米量级.LEOMolniya4期a64678-3.05048.12GEOC10210110010-1500100015002000Height/km图7 使用日月分析历表产生的位置误差随轨道类型的变化图。从左到右按轨道高度排序，依次是低轨、中轨近圆、闪电轨道和地球同步轨道.Fig.7 The position error caused by the analytical ephemeris of the sun and the moon is plotted against the type of orbit.Fromleft to right,they are sorted by orbital height,followed by low orbit,medium orbit near circular,Molniya orbit,and104Height/kmgeosynchronous orbit.2Height/10*km2.23Height/10*km3.5442-864卷JPL历表是由美国喷气推进实验室研究并发布的基于多项式拟合的行星历表.最早获得广泛应用的是DE9615，后续JPL历表发展了许多版本，比如19 8 5年的DE125根据旅行者二号飞掠天王星的数据对历表精度进行了修正 16.本文中轨道预报默认使用的是DE43014，这是目前比较稳定且精度比较高的版本，而最新的JPL历表是DE440.如图8 展示的几个典型版本JPL历表的日月位胡松华等：空间碎片轨道预报中动力学模型基准的选取4期置精度，以最新的DE440为基准.图8（左)为月球的位置精度，图8(右)表示太阳的位置精度纵坐标为日月位置的相对误差，即e=lDE*一TDE440l/rDE440l，其中DE*是相应历表下的日月位置，例如DE440，T D E2 0 0 月球的相对位置误差从DE200的10-6 降低到了到DE430的10-8 太阳的相对位置误差从DE200的10-7 降低到了DE430的10-9,这也得益于月球激光测距的发展 17.Moon105DE200-.*DE413DE403DE421-DE405DE430106sun10-6DE200-.*DE413-DE403DE421DE405DE43010-710-710810810-9109Feb1Feb15Mar1Mar15Dates in2021图8 主要版本的JPL历表相对DE440（D e v e l o p me n t Ep h e me r i s 440)的位置偏差示意图Fig.8 Schematic diagram of the relative position deviation of the major version of the JPL planetary ephemerides compared日月的位置误差间接导致的轨道预报位置误这些轨道从TLE编目中随机选取，包含了2.2 节所述差如图9 所示空间碎片的初始轨道是一条高度为的各个类型的轨道图10 画出了4条具有代表性轨36000km的近圆轨道，在半个月的轨道预报中，位道的预报误差情况，它们分别是低轨空间碎片、中置误差与预报时间有接近线性的关系如果采用轨空间碎片、闪电轨道空间碎片和地球同步轨道历表DE200，每天积累的误差为厘米级，而使用空间碎片，初始轨道开普勒根数如表2 所示.DE430,每天积累的位置误差不到1mm.结果表明，SRP导致的位置误差有明显的周期特性，周期即使不使用比较新版本的JPL历表，轨道预报的精正是轨道周期.根据预报的结果，我们发现虽然太度也不会受到明显影响阳辐射压力的大小随轨道高度的变化并不显著，但是SRP对高轨空间碎片的摄动效果依然远大于低3.4太阳辐射压(SRP)对轨道预报精度的影响轨空间碎片.对于中低轨空间碎片，不考虑SRP的本文采用的典型面质比如表1所示.这里使用预报结果在一天内的位置误差约为10 m，对于地球不同的动力学模型预报了8 0 0 条轨道来进行比较，10-10LApr1Jan5Mar16May 25 Aug3Dates in2021with DE440同步轨道和闪电轨道，误差约为10 0 m.Oct12Dec2142-964卷天文学报4期DE20010一DE403DE413DE421DE43010-110210310*40Fig.9 The influence of different versions of JPL ephemeris on the accuracy of orbit prediction1041031022图9不同版本的JPL历表对轨道预报精度的影响LEO-MEO-MolniyaGEO46810Timeld1214161820101101010Fig.10 The infuence of solar radiation pressure on the accuracy of the four orbits在计算SRP摄动力时，地影的处理十分关键，当空间碎片处于地影之内时，不会受到SRP的影响.一般情况下，为了简化计算过程，我们将太阳考虑为一个足够远且足够小的点，那么地影呈圆柱形.在以上的计算中，地影模型均是如此处理.而现实中，因为太阳与地球的距离为有限值，地影的形状应为圆锥形.另外，考虑到太阳并非一个点而是一个圆球，地影会被分成本影和半影两个部分.在本10图10 太阳辐射压对4条轨道预报精度的影响2030Time/h影中，SRP摄动不起作用，在半影中，太阳被地球部分遮挡，由未被遮挡的部分提供光压摄动力，即光压摄动力的大小由未被遮挡部分的太阳视圆面所占整体视圆面的比例决定 13这里将会研究不同地影模型的使用会对轨道预报结果产生多大的影响因为地影主要影响低轨空间碎片的轨道而对中高轨的空间碎片的轨道影响很小，我们这里以低轨空间碎片作为研究目4050607042-1064卷标.目标轨道的轨道高度分别为40 0 km(LEO1)、700 km(LEO2)、10 0 0 k m(LEO3)、130 0 k m(LEO4)，偏心率取0.0 1，轨道倾角均取45，升交点经度、近地点角距和平近点角均分别取30、40 和50.另外，由于空间碎片进出地影的时间很短，通常只有几秒钟，所以为了体现出进出地影边缘时不同地影模型的差别，在数值方法轨道预报中需要将步长调整到很小的水平，我们这里将步长调整为固定的1s.将步长固定为1s后会显著地影响轨道预报的速度，举例来说，一条7 0 0 km高度的LEO轨道，取动力学模型为30 30 阶次引力场、考虑大气阻力和日月第三体引力，在变步长情况下将其轨道向前预报一天需要1.8 s左右(不同的硬件或软件条件下会有一些差异)，而相同条件下的固定步长1s的预报需要6 8 s左右，时间消耗大大增加.幸运的是，从胡松华等：空间碎片轨道预报中动力学模型基准的选取4期我们的结果来看，绝大多数预报场景不需要考虑如此精细的地影模型差别，也就不必如此大程度地牺牲预报速度.图11展示了使用柱形地影与锥形地影两种情况下的轨道预报结果的差别，纵坐标是使用两种地影模型的预报结果的位置差异，或者说以锥形地影为基准模型，圆柱形地影产生的预报误差.我们可以看到，随着时间的推移，使用不同的地影模型的预报结果差异逐渐增大.位置误差在一天内可以积累到几个毫米，一些误差较大的轨道可以积累到1cm左右.我们注意到在4条轨道中，LEO3的位置偏差显著地大于其他3条轨道，这是因为LEO3轨道进出地影次数较多，所以地影模型对它的影响更显著.0.07LEO1.-.LEO20.060.050.040.030.020.0100Fig.1l The influence of Earth shadow model on the accuracy of the four orbits3.5潮汐形变摄动、后牛顿效应和大行星引力摄动对轨道预报精度的影响潮汐形变摄动、后牛顿效应和大行星引力摄动量级相对于前面的摄动都比较小，在潮汐形变摄动中，考虑太阳和月球使得地球发生的固体潮汐形变，且考虑到二阶近似 18 相对于日月引力，太.LEO3LE0410图11不同地影模型对4条轨道预报精度的影响W2030Time/h阳系内的其他行星的引力摄动较小，另外，当行星与地球的相对位置变化时，其摄动效果也会发生变化，因此，各个行星的引力摄动的量级随时间会发生变化.一般情况下，比较显著的是木星的引力，除此之外，水星、金星、火星、土星的引力摄动比木星的引力摄动小一到两个量级，天王星和海王星的4050607042-1164卷引力摄动比木星小3个量级或更多.这里研究大行星引力摄动对轨道预报的影响时，我们对比的是动力学模型是否包括木星、水星、金星、火星、土星的引力摄动的预报结果.关于后牛顿效应，我们考虑的是地球引力场的广义相对论效应修正.由于广义相对论下的动力学与牛顿体系下的动力学差别较小，所以考虑后牛顿效应相当于再增加一个摄动源，我们这里只考虑参数化的后牛顿效应的一阶项 11.我们计算的轨道包含2.2 节所述的各个类型，总数2 0 0 0 条.图12 展示了潮汐形变摄动、后牛顿效应和大行星引力摄动对轨道预报精度的影响，图中每个点代表一个空间目标，方形散点和圆形散点分别代表近圆轨道和大偏心率轨道，纵坐标表示一天内的位置偏差作为由地球产生的摄动，潮汐力作用和后牛顿效应的摄动作用都是随轨道的增加而减小的，而作为来源于地外的摄动，大行星引力摄动效果随轨道的升高而增大.潮汐形变摄动在一天内可以使低轨空间碎片积累百米以内的位置误差，而对于地天文学报球同步轨道，位置误差在1m以下.所以对于地球同步轨道的米级精度和中低轨的百米级精度，潮汐形变摄动是不可或缺的.后牛顿效应积累的误差比潮汐形变摄动小一到两个量级左右，对于中低轨空间碎片，一天会积累10 m以下的位置误差，对于地球同步轨道，误差则减小到0.1