编程模拟飞船加速变轨过程-物理基础篇（5） 摄动方程

编程模拟飞船加速变轨过程物理基础篇（5）

摄动方程

根据前面的内容，我们已经能够做到当已知某时刻位矢和速度时，求出轨道形状了。但是，这些都是在只受中心引力的前提下求得的。本篇在原作中是讲其他微小力（稀薄大气的空气阻力、其他星球引力、太阳辐射压力等）对人造卫星轨道的影响，在此，本篇仍按教材中的微小力来讲述，但我希望掌握这种方法，应用到飞船主动加速的过程中去，也即将微小力换为并不微小的推进力。

摄动方程·正文开始

由于摄动力存在，飞船的运动轨道在开普勒轨道附近，真实轨道并非Kepler轨道，所以引入密切轨道概念。每一时刻的位矢、速度能确定唯一的Kepler轨道，该轨道称为此时刻真实轨道的密切轨道。可以理解为，若摄动突然消失，则飞船将在密切轨道上运动。密切轨道随时间变化是时间的函数。因此轨道根数并非常数，研究出轨道根数随时间的变化也就研究出了飞船在受力时轨道的变化。

 

下面介绍如何用常数变易法建立摄动方程。

常数变易法

对于一阶线性齐次微分方程，它的通解为，c为任意常数

对于一阶线性非齐次微分方程，设它的解仍有这种形式，只是将常数c换为t的函数c(t)，即，那么求出该解的微分并代入中，得到，从而可求得通解x：

这种方法就是常数变易法

用常数变易法建立摄动方程

对Kepler方程,它的解形式为，

c1,c2这些积分常数由初始条件确定，可以是3个位置坐标和3个速度分量。由前面所讲，初始条件可确定轨道根数，因此在这里我们可理解为这些常数对应于轨道根数。

根据存在摄动加速度时轨道运动的微分方程：，

我们设该方程的解仍有上面的形式，只是c1,c2,c3,c4,c5,c6不再是常数。

对t求导：

因为每一瞬时的密切轨道是由真实轨道(r,v)确定的，所以密切轨道和真实轨道切点具有相同的加速度，即，于是，前两式可推导出，

 

将对时间求二阶导，得

结合Kepler轨道中，，所以

简单来说，我们通过上面得到了，该方程组可以导出轨道根数变化的微分方程

这就是我们所需要的摄动运动方程了，下面我们用矩阵的思想将它简化，以便在程序中运用数组计算。

运动方程的解，我们写作：，方程组我们写作

这里符号的意义如下：

（1式），（2式），注：

对（1）式左乘，（2）式右乘，然后作差，得到：

其中

注：L是一个反对称矩阵，即有

 

写出方程中各项的具体表达形式，在地心惯性坐标系上投影，取c为：

（在最后的求变化率时，我们会忽略i，Ω和ω）

矢量r在拱线坐标系上的分量矩阵为：

速度矢量v在拱线坐标系上的分量矩阵为：

r和v在地心惯性坐标系上的投影可以这样表示：

，，

其中

考虑到关系E-esinE=M和，综上，可从导出：

一般来说，f在轨道坐标系中表示{fr，fu，fh}起来较为方便，需要把它转换到惯性坐标系中：

将它展开：

接下来我们讨论轨道平面与赤道面重合，即fh=0，i，Ω和ω无意义时的情况

这时候，上面的关系式可化简为：

分开写就是：

θ是真近点角。

为了简化程序，在程序中，我们假设有：飞船只能向垂直径矢方向喷气加速，那么，上式就变成了：

至此，编程实现飞船变轨过程的物理基础基本结束，若对此方面觉得有兴趣，可以阅读北航赵育善老师和师鹏老师编写的教材，再次鸣谢！