厘米级以下碎片轨道演化计算方法与流程

1.本发明属于航天技术领域，具体涉及一种厘米级以下碎片轨道演化计算方法。


背景技术：

2.航天活动遗留在太空中的废弃的由于地球引力而在轨道上运行的碎片称为太空碎片。碎片包括失去控制的航天器的有效载荷、航天器运行期间散落的残骸和碎片、由于空间目标爆炸和撞击产生的碎片以及发动机喷出的铝片和油漆片、航天器覆膜脱离等。厘米级以下碎片由于面值比较大，受到的大气阻力、辐射压力洛伦兹力等的加速度都比同轨道上的大尺寸碎片和航天器受到的该类作用力加速度大得多。这就导致厘米级以下碎片的轨道演化计算误差比大尺寸碎片和航天器的轨道演化计算误差大得多。


技术实现要素：

3.本发明的目的在于提供一种厘米级以下碎片轨道演化计算方法，解决了现有厘米级以下碎片的轨道演化计算误差较大的问题。
4.本发明所采用的技术方案是：厘米级以下碎片轨道演化计算方法，包括以下步骤：
5.步骤1、计算碎片分布特征；
6.步骤2、计算碎片初始的位置速度；
7.步骤3、叠加碎片精确加速度模型，外推计算碎片的轨道演化。
8.本发明的特点还在于，
9.步骤1中计算碎片分布特征为，将碎片速度的分布密度函数取为：
[0010][0011]
式(1)中，v
col
是产生大量碎片的目标碰撞速度，vm是碰撞产生的碎片的最大相对母体速度，β＝8.69,δ＝7.2
×
10-3s·
m-1
。
[0012]
步骤2中通过公式(2)计算碎片初始的位置速度：
[0013][0014]
式(2)中，r和v为碎片参考位置和速度，rd和vd为碎片初始位置和速度，δv为单位方向，指向随机分布。
[0015]
步骤3具体为，通过计算碎片受到的地球引力、日月引力、大气阻力、太阳辐射压、洛伦兹力及波应亭-罗伯逊拖曳力来对碎片轨道演化进行计算，在惯性空间中表示的碎片的运动外推按照公式(3)积分计算：
[0016][0017]
式(3)中，r是相对于惯性空间从地球中心到碎片的位置矢量，fs和fm是
地球、太阳和月球对碎片产生的重力加速度，v是地球的引力势；fa和f
sr
分别是来自大气阻力和太阳辐射压力的加速度；f
l
是地球磁场的洛伦兹力产生的加速度；f
pr
是波应亭-罗伯逊拖曳力产生的加速度。
[0018]
步骤3中地球引力势按照公式(4)计算：
[0019][0020]
式(4)中，g是引力常数，me是地球质量，r＝|r|，re是地球赤道半径，λ和ψ碎片位置的地心经纬度，和是n阶m度的归一化球谐系数，是归一化的勒让得多项式，n是公式使用的最大阶次。
[0021]
步骤3中太阳和月球对碎片产生的重力加速度按照公式(5)和(6)计算：
[0022][0023][0024]
式(5)和(6)中，g是引力常数，ms和mm分别是太阳和月球质量，r
ds
和r
dm
分别是太阳和月球到碎片的位置矢量，rs和rm分别是地心到日心和月心的位置矢量，r
ds
和r
dm
是r
ds
和r
dm
的模，rs＝|rs|，rm＝|rm|。
[0025]
步骤3中太阳辐射压力的加速度按照公式(7)计算：
[0026][0027]
式(7)中，太阳常数为s0＝1.3608
×
106ergs
·
cm-2
·
s-1
，as是碎片的光压截面积，md是碎片质量，c是光速，r
s0
＝1au，r
sd
是太阳到碎片的距离，q
pr
是太阳光压系数，r
sr
是从太阳到碎片的方向单位矢量。
[0028]
步骤3中碎片受到的大气阻力加速度按照公式(8)计算：
[0029][0030]
式(8)中，大气阻力系数为cd＝2.2，a为碎片的大气阻力截面积，ρ是碎片轨道位置处的大气密度，v是惯性空间的碎片速度矢量，v＝|v|。
[0031]
步骤3中碎片受到的洛伦兹力加速度按照公式(9)计算：
[0032][0033]
式(9)中，介电常数为ε0＝8.854187817
×
10
12f·
m-1
，碎片表面在空间物理环境中充电后达到电荷平衡的电压为u＝ 5v，s是碎片平均半径，md是碎片质量，ω是地球自旋的角速度矢量；碎片所在位置的磁场强度按照公式(10)计算：
[0034][0035]
式(10)中，磁标量势能为v
mag
按照公式(11)计算：
[0036][0037]
式(11)中，地磁平均参考半径为a＝6371.2km，θ和分别是碎片所在位置的地心固连坐标系经度和纬度，t为时间，n和m分别是磁场强度计算的度和阶，和分别是地磁场的高斯系数，是归一化勒让得函数，n
gm
是使用的最大的磁场强度计算阶次。
[0038]
步骤3中碎片受到的波应亭-罗伯逊拖曳力加速度按照公式(12)计算：
[0039][0040]
式(12)中，太阳常数为s0＝1.3608
×
106ergs
·
cm-2
·
s-1
，a为碎片的大气阻力截面积，md是碎片质量，c是光速，r
s0
＝1au，r
sd
是太阳到碎片的距离，q
pr
是太阳光压系数，r
ds
是太阳到碎片的位置矢量。
[0041]
本发明的有益效果是：本发明的厘米级以下碎片轨道演化计算方法，可以有效降低现有厘米级以下碎片的轨道演化计算误差，有利于航天器日常监管过程中对厘米级以下碎片轨道长期预报，延长碎片对航天器可能的碰撞预报时间，增大碰撞预警的有效时间。
附图说明
[0042]
图1是碎片轨道某时刻分布在惯性系xz平面分布图；
[0043]
图2是碎片轨道某时刻分布在惯性系yz平面分布图。
具体实施方式
[0044]
下面结合附图以及具体实施方式对本发明进行详细说明。
[0045]
本发明提供了一种厘米级以下碎片轨道演化计算方法，包括以下步骤：
[0046]
步骤1、计算碎片分布特征，将碎片速度的分布密度函数取为：
[0047][0048]
式(1)中，v
col
是产生大量碎片的目标碰撞速度，vm是碰撞产生的碎片的最大相对母体速度，β＝8.69,δ＝7.2
×
10-3s·
m-1
。
[0049]
步骤2、通过公式(2)计算碎片初始的位置速度：
[0050][0051]
式(2)中，r和v为碎片参考位置和速度，rd和vd为碎片初始位置和速度，δv为单位方向，指向随机分布。
[0052]
步骤3、叠加碎片精确加速度模型，外推计算碎片的轨道演化；具体为，通过计算碎片受到的地球引力、日月引力、大气阻力、太阳辐射压、洛伦兹力及波应亭-罗伯逊拖曳力来对碎片轨道演化进行计算，在惯性空间中表示的碎片的运动外推按照公式(3)积分计算：
[0053]
[0054]
式(3)中，r是相对于惯性空间从地球中心到碎片的位置矢量，fs和fm是地球、太阳和月球对碎片产生的重力加速度，v是地球的引力势；fa和f
sr
分别是来自大气阻力和太阳辐射压力的加速度；f
l
是地球磁场的洛伦兹力产生的加速度；f
pr
是波应亭-罗伯逊拖曳力产生的加速度。
[0055]
其中地球引力势按照公式(4)计算：
[0056][0057]
式(4)中，g是引力常数，me是地球质量，r＝|r|，re是地球赤道半径，λ和ψ碎片位置的地心经纬度，和是n阶m度的归一化球谐系数，是归一化的勒让得多项式，n是公式使用的最大阶次。
[0058]
其中太阳和月球对碎片产生的重力加速度按照公式(5)和(6)计算：
[0059][0060][0061]
式(5)和(6)中，g是引力常数，ms和mm分别是太阳和月球质量，r
ds
和r
dm
分别是太阳和月球到碎片的位置矢量，rs和rm分别是地心到日心和月心的位置矢量，r
ds
和r
dm
是r
ds
和r
dm
的模，rs＝|rs|，rm＝|rm|。
[0062]
其中太阳辐射压力的加速度按照公式(7)计算：
[0063][0064]
式(7)中，太阳常数为s0＝1.3608
×
106ergs
·
cm-2
·
s-1
，as是碎片的光压截面积，md是碎片质量，c是光速，r
s0
＝1au，r
sd
是太阳到碎片的距离，q
pr
是太阳光压系数，r
sr
是从太阳到碎片的方向单位矢量。
[0065]
其中碎片受到的大气阻力加速度按照公式(8)计算：
[0066][0067]
式(8)中，大气阻力系数为cd＝2.2，a为碎片的大气阻力截面积，ρ是碎片轨道位置处的大气密度，v是惯性空间的碎片速度矢量，v＝|v|。
[0068]
其中碎片受到的洛伦兹力加速度按照公式(9)计算：
[0069][0070]
式(9)中，介电常数为ε0＝8.854187817
×
10
12f·
m-1
，碎片表面在空间物理环境中充电后达到电荷平衡的电压为u＝ 5v，s是碎片平均半径，md是碎片质量，ω是地球自旋的角速度矢量；碎片所在位置的磁场强度按照公式(10)计算：
[0071]
[0072]
式(10)中，磁标量势能为v
mag
按照公式(11)计算：
[0073][0074]
式(11)中，地磁平均参考半径为a＝6371.2km，θ和分别是碎片所在位置的地心固连坐标系经度和纬度，t为时间，n和m分别是磁场强度计算的度和阶，和分别是地磁场的高斯系数，是归一化勒让得函数，n
gm
是使用的最大的磁场强度计算阶次。
[0075]
其中碎片受到的波应亭-罗伯逊拖曳力加速度按照公式(12)计算：
[0076][0077]
式(12)中，太阳常数为s0＝1.3608
×
106ergs
·
cm-2
·
s-1
，a为碎片的大气阻力截面积，md是碎片质量，c是光速，r
s0
＝1au，r
sd
是太阳到碎片的距离，q
pr
是太阳光压系数，r
ds
是太阳到碎片的位置矢量。
[0078]
此外，碎片个数的尺寸分布还可以按照公式(13)估计：
[0079][0080]
式(13)中，n
debris
是尺寸大于d的碎片个数，m
tot
是碎片总重量，β＝1.32是修正系数，ρ＝2.7
×
103kg
·
m-3
，d是碎片尺寸。
[0081]
实施例
[0082]
解体碎片总质量设为m
tot
＝740kg，计算不同尺寸碎片个数为n
debris
(≥0.2m)＝14，n
debris
(≥0.1m)＝60，n
debris
(≥0.01m)＝6587，n
debris
(≥0.0018m)＝2.18
×
105，n
debris
(≥0.001m)＝7.22
×
105。计算碎片初始位置速度，外推计算碎片的轨道演化，碎片轨道某时刻分布在惯性系xz平面分布如图1所示，在惯性系yz平面分布如图2所示。