航天器轨迹规划的变分技术

1、1航天器轨迹规划的变分技术摘要：本文描述了一种基于变分法的轨迹规划算法，它可以解决六自由度的航天器对接和接近操作问题。目标函 数综合考虑燃料消耗、避障距离和到达时间。以非线性轨道动力学方程作为动力学约束条件，推导了欧拉-拉格 朗日方程，依据庞德里亚金原理推导了最优控制的输入，给出了推力器实际的开关曲线。采用对初始条件不敏感 的间接配点法解决相应的边值问题，解的开拓可以进一步改进算法的鲁棒性。文章讨论了对欧拉-拉格朗日方程 和横截条件的处理，以使其适合目前的配点法的应用形式。文章最后通过一个与翻滚卫星进行首尾相连的对接机 动对所研究的方法和结果进行了验证。关键词：航天器，轨道，空间探测，算法，航

2、天器任务规划1.1引言一个成功的自主航天器必须能够实现智能导航。智能导航能够找到并沿着无碰撞的路径以合 理的时间和燃料穿过预定空间。很多自主航天任务目前已经发射，还有一些其他任务，包括 DAPPA（美国国防高级研究计划局）的SUMO/FREND服务航天器，AFRL （美国空军研究实验室） 和Lockheed-Martin公司的ANGELS观察卫星，目前都在设计之中。尤其是，DARPA的轨道快 车任务，该任务在技术方面向我们展示了自主交会、抓捕、燃料加注和重构等技术，并在2007 年7月成功完成了其主要任务。轨道快车 项目的目的是研发新一代在轨服务航天器。 SUMO/FREND目的是通过交会对接

3、和机械臂的抓捕技术，为失效的客户卫星提供轨道定位服务。 所有的这些飞行器必须在限定的环境中运动，SUMO/FREND和轨道快车之后的服务航天器也需 要在预定的环境中与目标接触。为了能够完成满足需求的先进的自主任务，必须找到一种燃料消耗较低同时又能避免碰撞的 轨迹规划算法。这是一个很困难的问题，因为这些飞行器相对于它们的对接目标一般是从不同的 （初始）轨道开始的；因此，这些航天器的轨迹规划必须综合考虑远程轨道机动和近程的临近操 作，而远程机动伴随有非线性动力学和燃料消耗问题，避免碰撞对近程临近操作是极为重要的。轨道机动中的轨迹规划技术已经成功使用了很多年，但是这些技术通常并不能处理避障问 题。另

4、一方面，大部分处理避障的轨迹规划技术都在陆地应用，如轮式机器人或机械臂。遗憾的 是，这些算法不能很好地在空间中应用，因为它们所基于的假设在空间并不适用。特别是动力学 效应不能被忽略。而且，在燃料的节省和保存方面，空间应用比地面应用更重要。对于给定的轨 道机动问题，到达时间也是一个很重要的方面，燃料消耗最少的轨迹可能需要无穷大的转移时间。 因此，需要这样的轨迹规划技术，这些技术能够生成满足轨道和姿态动力学限制的路径，权衡燃 料消耗和避障距离的约束，满足宽变的终端范围限制，能应对和躲避移动障碍物，以及能优化完 成时间。本文提出了一套基于变分法的算法，能够满足以上要求。该算法以惩罚燃料消耗和障碍 物

5、间距离为目标函数，找到一条满足一系列边界约束条件且使目标函数最优的轨迹，这些边界约 束条件允许终端位置和时间固定或自由。1.2假设航天器轨迹规划的研究领域极度宽阔，不同的研究者在试图解决特定的问题时，有时采用不 同的假设。比如，不同的研究者，要么包括，要么不包括多引力体作为次要影响因素；有些考虑 地球非均匀重力场和牛顿均匀场的作用（称为J2影响）。类似的，有些研究者试图直接包含导 航误差，而另一些人则没有。一般的轨迹规划研究，主要来自移动机器人领域，它们的假设条件 通常不同于空间飞行器的轨迹规划条件。特别是，机器人的研究者们通常（尽管不是一直）忽略 所有系统动力学，他们假设机器人能够精确跟踪任

6、何给定轨迹，即使轨迹是分段线性的，这是由 于距离很近时动力学作用可以忽略。反而文献中主要关注在约束环境中找到无碰撞的路径。显然， 航天器在约束环境中的系统动力学是不能忽略的。本文将两个领域的问题进行结合。所以，为了使问题既模型简单又易于运算，这里采取一定 的假设。特别是，航天器动力学采用标准的二体模型，J2项摄动、N体动力学、大气阻力等都忽 略不计。这些作用是影响航天器轨迹规划的重要因素，近几年的航天器轨迹优化在这方面做了很 多工作。但是，本文的目的是提出一种算法，能够在一定约束条件下找到转移轨迹，并不是在处 理扰动方面开辟一个新领地，因此这里对摄动项不加以考虑。总之，尽管这些扰动会使数学模型

7、 变复杂，但是不能从根本上改变算法。在有障碍物限制的轨迹规划问题中，寻找全局最优解的方法中有我们熟知的NP-hard问题， 它被列为已知的最困难的计算问题之一。因此，很多有障碍物限制的轨迹规划算法，并不会试图 解决全局问题，而是试图得到局部最优的结果。在很多情况下，尤其障碍物相对较少时，局部的 最优解是可以接受的。本文也考虑这种情况，我们并不会试图寻找全局最优解，因为这对现有的 动力学模型进行计算太困难了。曾用到过的寻找全局最小值的方法的综述在下一节给出。1.3研究背景解决特定的轨道机动问题有很多经典的方法，这些方法在轨道力学领域是大家所熟知的。转 移时间自由的轨道转移，如霍曼转移1和双椭圆转

8、移2已经使用很多年了，同样的情况还有时间 自由或受限制特殊点拦截3,4。对于更加复杂的轨道机动，如考虑多引力体情况，圆锥曲线拼接 法通常被认为是已知的首选方法。尽管这些算法所获得的结果并不是最精确的，但是它能够为更 复杂最优化方案得到一个比较合适的初始值。这些方法一般都能得到燃料最优解，而且考虑动力 学限制。虽然这些方法假设推力是脉冲形式且大小不受限，但是考虑到化学推力器的输出能力， 以及大多数轨道转移问题的时间尺度，这个假设通常是有效的。但是这些方法既不能处理避障问 题，也不能处理姿轨耦合的机动问题。有很多研究将线性或者非线性最优方法应用在航天器轨迹规划中。大体来说，线性方法仅仅 在近距离操

9、作中很有效，因为真正的轨道动力学，甚至是最简单的二体问题也还是有很强的非线 性。而当两个航天器在相对于轨道周期的短时间内近距离相对运动时，它们的相对运动可以被线 性化，这就是著名的Hill方程。通过这样的简化，线性化最优方法可以有效使用。例如，目前 有很多基于线性或者非线性的创新方法，解决单个或者多个合作航天器的机动问题。在这些方法 当中，对于航天器的轨迹规划来说，MILP（混合整数线性规划）法在解决有约束线性最优问题时 是一种很有效的方法5问。MILP方法已经应用在包括碰撞规避及羽流的接近操作和编队飞行中。 这种方法必须使用样条或者其它近似方法对解空间（即所有可能的轨道集合区间）进行近似。定

10、 义一个目标函数，采用MILP法找到满足给定约束条件使目标函数最小化的样条曲线。因此，这 种方法是将轨迹规划问题转化为限制性最优问题，进而用一种已有的优化算法解决。MILP法在 计算方面异常高效。遗憾的是，MILP法似乎不能应用于明显非线性的动力学问题，因此，它不 能使用在轨道注入或者轨道（平移）和姿态（旋转）耦合的问题中。单元分解法也已经被用在最优轨迹规划中。单元分解法通常是形成一个离散化的网格状的状 态空间，然后使用一些标准的搜索技术搜索得到的离散空间。单元分解法已经用于解决航天器的 轨迹规划问题，著名的有文献7，他使用基于A *的寻找方法，应用于航天器近距离操作，得到 一个六自由度的轨迹

11、。路线图方法也能够产生一个邻接图表示搜索空间，然后用一个标准的图像 搜索方法从图像中得到一个路径。但是，与单元分解法不同的是，图像法不是状态空间离散分解 的结果，相反图像法试图获得空间的“内在连续性”。Generalized Voronoi diagrams在航天器在 轨迹规划中的应用被广泛接受8,9 这些方法主要用来处理避障问题，它们并不是燃料最优，尽 管燃料消耗在处理要求相对小间距的问题时似乎很合理。它们不能解决终端时间或者终点位置限 制的问题，并且不能够使用实际的推力器模型。也不能解决实际的轨道动力学（附加扰动）和姿 态耦合的问题。但是，对于大部分近距离操作，相对轨道动力学特性可以忽略，

12、而且轨道和姿态 运动可以解耦，那么，这种方法可以得到很好的结果。人工势函数制导方法是一种最简单并且最流行的轨迹规划方法。在它的基本形式中，建立了 一个叫做势函数的数学函数，它在目标点上取得全局最小值，在每个障碍物的中心取得局部最大 值（这是一种紧约束）Ji。该路径的产生是通过一个基于势函数的梯度下降算法完成的，然后 将得到的轨迹投影到基本的欧式空间中。人工势函数制导方法在航天器轨迹规划中的应用也很广 泛1214。以上的这些方法都考虑是实际的推力器，固定的障碍物，自由的到达时间。它们并没有直接 最小化燃料消耗，但是他们使用的是Hill方程并且应用在近程接近操作中，所以燃料消耗很少。 另一个方面，

13、它们也不适用于远程轨道机动，不能处理自由终端位置和移动的障碍物问题。它们 没有考虑姿态规划，因此就不能解决与旋转目标对接的姿态与轨道耦合问题。最近，有人提出了用于解决多颗卫星星座重构的轨迹规划问题的一种两阶段的优化方法。在 这种方法中，可以利用快速扩展随机树法得到一个可行但不是最优的轨道。该轨道可以作为平滑 或者优化算法的初始估计。该算法或者是由针对线性系统动力学的线性优化方法Pl，或者是由 非线性优化方法如伪光谱法16组成。这种方法可以被成功应用在考虑碰撞和定向约束，为重构 卫星星座寻找合理的轨道上。与以上方法不同的是，该算法可以寻找到燃料最优轨迹，可以处理 姿态规划和移动目标的问题。但是该

14、方法仅仅能够应用在轨道动力学并不重要的深空星座当中。也有将非线性优化方法应用在航天器轨道优化问题中。一般的，非线性方法通常使用在单个 航天器在单一引力体的不同轨道间或者在不同引力体间机动多个轨道周期的情况。这些方法可以 分为直接法和间接法。基本上，直接法是将最优轨迹规划问题转化为NLP（非线性规划问题），通 过采用一种非线性规划算法求解。最常用的方法是直接打靶法和直接配点法。直接打靶法主要通 过参数化受设计器控制的变量（通常是控制输入）实现。系统动力学模型从初始时刻至终端时刻 计算，进而计算结果的优化值；然后反复调整参数对优化结果进行改进。直接配点法通常被认为 是另一种可供选择的方法，轨迹和控

15、制输入都通过样条函数来参数化，进而调整样条函数的系数 来改进优化结果。但是，Betts指出，配点法也可以看作是打靶方法的一个特例，就是将打靶法 的积分器采用Hermite-Simpson算法，并且参数化的参数数量与积分器时间步长数量相一致。间接法则是推导一组微分方程：Euler-Lagrange方程，方程的解即是最优轨迹规划问题的极 值曲线。简单的来说，必须解决边值问题去得到这些解。相比直接法是一个优化问题，间接法是 一个寻根问题。像直接法一样，有很多解决边值问题的数学方法，包括间接打靶法，间接多重打 靶法，间接配点法。直接法渐渐的成为轨迹优化问题中的首选方法是有原因的。首先，直接法比 间接法

16、更具有鲁棒性，有更广泛的收敛范围。此外，由于间接法以伴随变量的形式引入其他状态 变量，一定要保证对这些变量进行合理估计才会使迭代收敛。但是由于伴随变量并没有直接的物 理意义，所以提供一个好的初始估计往往很困难。并且使用间接法要求对每个问题推导 Euler-Lagrange方程，这就使得这种方法在完成一般的轨道优化中更具有挑战性。但是，在某些情况下，间接法比直接法更有优势。特别是它能直接（至少在理论上）包含几 乎所有假设约束条件，不用事先用数值解对约束条件进行修正。使用直接法，通过数值解直接满 足约束条件，因此，它的解要么很笼统，要么仅适合于特定的问题。特别是在处理较强约束条件 的问题上，例如使

17、用自由开关推力器时，就很难使用直接法，因为在解决这个问题之前需要知道 推力器点火的时刻，这样才能在整个过程中处理变化中的系统动力学。间接法在解决非常困难的 轨道优化问题中相当成功，比如小推力的多体转移问题17,18。与直接法类似，这些方法通常是基于牛顿法的变形形式。遗憾的是，这些方法的普遍特点是 其局部最优性；它们全部需要利用梯度信息，对梯度的模或指标函数值进行迭代改进。结果都是 趋于收敛到局部最小值，在不对目标函数形式作附加要求或者满足微分、代数约束时，结果不能 保证收敛到全局最小值。对与这些优化算法在航天器轨迹优化中应用的综述，见文献19。直接 法和间接法在“六自由度轨道机动”章节有详细讨

18、论。以上提到的方法中，没有哪个方法是普遍全局最优的。然而，有很多人试图解决全局最优问 题，至少是在统计的意义上。模拟退火法就是其中之一2。21。在模拟退火法中，随机扰动可以应 用于最优解最好的估计。扰动大小的变化是依据“退火过程”，扰动大小本质上是随着时间减小的一个函数。此外，扰动的大小还是随着当前目标函数的值而减小的函数。这样的话，这种算法 可以用统计学的方式搜索到很大部分的解空间。遗憾的是，很难严格保证轨迹规划问题的随机退 火的全局收敛性，而且选择合理的退火过程通常需要反复试验。寻找全局最优的另一种方法是遗传算法22对于不同类型的遗传算法通常使用不同的变量， 但是一般来说首先需要为可能的问

19、题解建立一个编码表。在轨迹规划中，编码表由一些基向量的 权重或神经网络的权重组成。一组候选解的产生，通常是从编码表中随机选出。这些候选解通过 一个适应性函数进行评估，选择一部分子集作为新的候选解。对这些留存的候选解随机的扰动或 交叉，会产生新的候选解集合。反复迭代，直到达到最大设定迭代次数或候选解达到满意的品质。 遗传方法很普遍，可以应用到任何其他优化方法可以应用到的地方。但是，遗憾的是，遗传方法 不能很好利用梯度信息，因此对问题的计算量很大。因此还没有被广泛应用到轨迹规划问题中。1.4数学基础1.4.1直接法形式。本文使用Lagrange表述，目标函数取成积分形式，并对末端时刻f的状态进行惩

20、罚:J L(),u(),= j七 L L(),u(),山o直接法解决最优问题开始是构造一个目标函数和一些约束条件。目标函数有很多不同的规范动态约束条件：z(z(),u；根据待解决的问题本身的特点，状态变量双)、控制量 u也需要代数约束：g - gL,uo,k为。最优问题转化为求解最优控制量u*()和相应 状态量x*()，通过满足一些约束条件使得目标函数最小。有很多特定的情况，其中一种典型的情况是状态变量和输入控制量都是二次的，微分约束条 件是线性的，没有代数约束条件。这样的优化问题可以求得解析解。对于不满足这些标准的问题， 可以通过数值计算得到数值解。事实上这是相当复杂的，除了相当少的情况下，

21、允许的控制输入 是一个无限维数的函数空间。直接法通过限制函数控制量输入搜索空间使得最优化问题容易计 算。这个近似函数通常是B样条。因此，轨道最优问题可重新表示为：u* w w*b(t)% L z(t), wtb(t), t dto切 g V g r x (t), wb(t), t ! g和j L满足约束：” 2 M(t), wt,u。因此，问题变为NLP(非线其中，b(t)是基函数向量，w是常数权重向量。寻找权重向量w*最小化下式: J r Z (), wtb(t), t = j性规划问题)。有很多解决非线性规划问题的方法，包括：直接打靶法、多次直接打靶法和直接配点法。1.4.2间接法Eule

22、r-Lagrange 方程本文采用基于间接配点法的变分方法。变分法的提出，是为了解决找到使特定约束最小化(或 最大化)的函数。变分问题描述如下：23,24J z(t),u(t),t= jlf Lz(t),u(t),tdt问题一：给定指标函数t0，找到最优轨迹z*(t)和输入控制量u*(t)，使得指标函数加最小(或最大)，并且满足动态约束条件)=f z,(t),t。变分法的轨迹规划允许对燃料消耗和避障进行折中。使用间接优化方法允许非线性动态约束 和开关控制，满足推力器阈值约束。此外，到达时刻和终端位置都是最优的。这里，数学基础要 求总结出一个可以完成6自由度轨迹规划的间接配点法方法。“间接”这个

23、词是用来区分这两种方法，其中一种试图通过直接在适当的函数近似空间中近 似最小指标函数J日，如顺序二次规划法(“直接”法)；另一种试图通过解一组微分方程来 最小化指标函数(间接法)。间接变分优化方法的目的是要推导一组微分方程，它的解可使J日 最小。首先，通过在指标函数中加入时变Lagrange乘子项来强制系统动态约束：J z(t),v(t),t=j Lz(t),v(t),t+Xt(t)f z(t),v(t),tldt jf L z(t),v(t),tdtt0其中，e = f 如,知扇a,at0Lagrange乘数* (t)叫做协态向量，它们可以以Euler-Lagrange方程的解的形式表示:X

24、 (t)=牛z *(t), v* (t), X (t), t(2)J d!z (t ) = 7 .7 (t ) = z最小化 ，需满足边界约束条件z(0) z0,z(f) zf和动态约束条件 z = (t) = fz(t),v(t),t, 一,24,25。注意到，Euler-Lagrange方程是与边界条件耦合的二阶微分方程组。除了一些特殊情况，即当系统动态函数f (D是线性的，并且指标函数仅包含二次型项，Euler-Lagrange方程通常是非线性的。如果推力器不受限制，可以通过下式解得最优控制输入少(t) 250 = a z(t), V(t), tav如果推力器是有饱和的或者有其他限制(如

25、bang-bang推力器)，最优控制输入V*(t)根据庞德 里亚金极小值原理25 来选择：L z(t), X(t), v*(t), t L z(t), X(t), v(t), t 其中，许可控制输入v(t)是满足约束条件的任意输入值。另一种情况，通常可能会通过状态向量 和协态向量得到v*(t)的解v*(t)=nz *(t ),X*(t)。解决任何微分方程的问题，都需要有边界条件。对于变分问题，边界条件通常是分开的，部分边界条件在t=0时是已知的，部分边界条件.t=t在 f时是已知的。这就被称为BVP (边值问题)。确切的边界条件取决于问题本身。例如，如果到达时间、起始点位置和终点位置给定，则边

26、界条件可以通过以下形式表示：z (t0) = Z0z (t f ) = Zf如果终端时刻自由，终端位置取决于终端时刻，即z =0，则边界条件是I-1dOL rz(t ),v*(t ),X(t ),t -Xt(t )(t ) = 0Lf f f f f dt fz (t f ) =O (t f) z (t0) = z01.4.3配点法求解边值问题有很多方法。最流行的是打靶法和它的变体25。打靶法必须为初始条件的未 知部分确定一个初始估计，对微分方程从0到tf积分，用理想值与实际值的差，用Newton法或 者其他相关方法更新估值。打靶法相当精确，但是对初始条件估计的准确性很敏感。遗憾的是， 协态方

27、程没有一个明显的物理意义，因此确定一个精确的初始估计可能很困难。配点法虽然不能和打靶方法一样精确，但是它对初始解估计的准确性没有那么敏感。配点法 使用三次样条或者其他的逼近方法去逼近微分方程的解2629。z(t) w cTb(t)zX (t) w cTb(t)这里向量Cz和向量CX是权重矢量，b(t)为基函数向量。每个基函数都以离散的时间节点t t t ., 一 cc .(t0，t1，.，m-1)为中心。配点法选择权重矢量z和cX，使得插值解满足边界条件和每个节点上 的微分方程。cb(t) = fb(t), 丫 cTb(t), cTb(t), t , tz i L z i X i i i/J

28、- r-|cTb(t) = - 一 Ub(t), cTb(t), 丫 cTb(t), cTb(t), t , tX i6zz i X iz i X iiicTb（t0） = z0crb（t 1） = Zf其中，前两个方程来自配点法要求，后两个方程确保满足边界条件。边界值问题由此转化为求解 权重矢量Cz和土的一个AME（代数矩阵方程）。对于一个4 -m的代数方程，其中n是状态空 间的大小，m是节点的数量，它的结果是：0 =中（cz, zJ = W （提该方程必须解出系数&。为了完成转化，假设微分方程的右边是利普希茨条件。对本文来说， AME通过Levenberg-Marquardt方法进行求解，

29、该方法同时具有梯度下降法的鲁棒性和Newton 法的快速收敛性。为了能够进一步加强配点法的鲁棒性，AME求解器可以融合延拓算法。这种算法可以用来 针对高非线性问题改进配点法的收敛性2729。用延拓修正方程，再去求解，例如，可以通过在 微分方程中为非线性项乘以一个增益。该增益起初为零，成为一个线性问题，这样可以很容易求 解。在每一步的迭代中，增益不断增加，前一步的解作为当前步的初始估计使用。当增益达到一 致之后，这个方法就可以找到原问题的解。1.4.5姿态描述为了规划六自由度轨迹，必须在系统状态空间中增加姿态的描述。所有的最小的（三参数） 姿态描述都在取值范围内某处有奇点。这一事实部分解释了非最

30、小描述，如四元数法，普及应用 的原因，因为它没有奇点。但是遗憾的是，对于轨迹规划问题，使用非最小姿态描述的方法会带 来一个很严重的问题：配点法的AME雅克比行列式是奇异的，这样就会给基于梯度的如Levenberg-Marquardt的AME求解器带来问题。克服这种问题的一个解决办法就是在有利的地方 使用具有动力学奇异的最小姿态描述。一种描述方法就是修正的Rodrigues向量，一个相对新的 三参数描述方法，仅仅是在从原姿态旋转360时才会出现奇点30。修正Rodrigues向量的元素=+ 2气卞，用四元数q=m，河定义为:8 b =1 +门2b8 =相反的:1 + b Tb1 -b Tb1 +

31、 b T b修正Rodrigues向量与Euler轴k和角度之间的关系如下，,/中、b = k tan()4其中，定义Euler角的范围0 V中 . . . U其中，p代表在非旋转坐标系中的航天器位置，此坐标系的原点固连于对接目标航天器；v为相 关的平动速度向量；at为对接目标的加速度，b为修正Rodrigues向量，为在体坐标系中 的角速度向量；m为航天器的质量，H为转动惯量向量。矩阵*(b)为旋转矩阵，将体坐标系转 化为惯性坐标系。该矩阵可方便地用等价的四元数法表示为R(q)=(中-t8)1 + 2跖t 2qS()1.5六自由度轨道机动本节推导了有障碍物情况下姿轨耦合轨道机动问题的指标函数

32、。算例中将用到的指标函数具有以下形式J=七Lz,v, T + Lz + L+ 人tf (z,v, T)dtt controlobstacletime,，V、一 1=ftf |u(t)| + t(t)tRt(t) + aEg (v) h (p) + y +人 tf z,呼(z, X) I dt TOC o 1-5 h z t0 L 1i=1 (7)式中f ()表示航天器动力学方程(6) ； Lcontrol =U(t)1 + T(t)T Wt为控制力的罚项；L =ag( v)h (p)T =obStaCl(；i=1i i 为障碍间距罚项；time-为完成时间罚项。1.5.1燃料消耗当轨迹指标函数

33、具有以下形式时，轨迹燃料消耗最小25J = f tf |u(t )| + t1(8)然而，由该指标函数得到的控制律:u X -1/ mu * (t) = 0-1/ m X 1/ mI maxv ,i是不连续的，因此导致欧拉一拉格朗日方程右端也不连续。为克服该问题，需使用与方程(9 )近似的连续控制律(见图1)图1控制律10 max 1。和m 1_ mu (e -1)1 2( e 1)e + 2m 2u(u m2 一 1)e + um2 一 122m(um2 + 1)e + u m2 一 1 TOC o 1-5 h z maxmax._ (2u2 m2 一 u)e1 (um2 一 1)e + u

34、m2maxmax._ (um2 一 1)e 2 一 (3u m2 一 1)eI = maxmax22m2(um2 + 1)e + u m2 一 1maxmaxu (t)=umas X+i1 v,i1s X+i2 v ,i21 f2m 2v ,iX ,i - i2S2X ,i - i1umaxX -2m - u-2m-u X , -(1 E)m-u +e /m-(1 -词-u+e /m X . -1/m-1/ m X 1/ m1/ m X . (1 -e )m - u +e/m (1-e)m-u+e/m X 2m - u(10)当式中e等于1时，上述控制律与最小燃耗律完全相同。本文中大部分例子将

35、使用该连续方法，e取值范围从0开始增长到e约等于1。姿态执行机构最小能量指标函数如下:J = j tf p T(t)t t(t) +t0式中，p为正常数。由此可得控制律为TmaxT*(t) = (人tH-1) /2p i广，max(人tH-1) 2p T-2 p T (人tH-1) 2 p t(11)1.5.2障碍间隔作为燃料消耗惩罚项的补充，在指标函数中引入障碍距离惩罚项:(12), P t , 一 j . 一一 ，一式中，Oj为与第j个障碍固联的坐标系原点位置O ，b 一 ,，一j与Oj表示障碍坐标系的方向。d(.)为航天器到障碍的距离，是障碍物外形的函数；5j为常值标量增益。障碍罚函数g

36、()值随着到障碍物边缘的距离的减小而增大，且当与障碍物的距离为某值时，罚函数值很小或为0。为使欧拉一拉格朗日方程满足利普希茨光滑条件，英严七须满足利普希茨条件；例如当g 为C 2时，该条件成立。一种满足这样要求的方程由分段三次函数组成:d = d (p, p , c ) sn rOjOjk (a d 3 + a d 2 + a d + a )1 n 2 n 3 n 4式中，S通过对相对距离dr进行缩放决定了 g ()非零区间的大小；k决定了罚函数的最大值;数3，a4的选择标准是使g ( 40 k)，g (s) = 0岛(d )/dd (s) = 0a2g(d)/dd2(s) = 0胞?r rr

37、 r；见见图2。图2分段三次障碍罚函数，光滑化的速度标准化函数有两种策略可以使方程(12)最小。一种是通过使轨迹离开障碍物以减小障碍物间距。另一种 是通过增大接近速度以减少接近障碍物消耗的时间(事实上，只需令-0就能使消耗任意小)。显然，第二种方法是不适用的。作为对相对于障碍物速度增加的补偿，需要引入：0j=1(13)J()丁 +a jgdr (p，pq .)叫(p，pq ，v，气.)+j jj jj j,“，_-,V 式中，0.为第j个障碍坐标系的速度，0.为其角速度。航天器与障碍物的相对速度r可由相 对距离对时间求导得出，且速度标准化函数”随七的增大而增大。很容易看出，要对速度改变进行精确

38、校正，h(Vr )应为七；但是绝对值的偏导在V广0时 无意义，因此，欧拉一拉格朗日方程在七=0处也无意义。为解决此问题，可以在= 0附近对绝对值函数进行三次插值以使其平滑:h(v )=rb v 3 + b v 2 + b v + b1 r 2 r 3 r 4V 其它其中，系数b1,，b4需满足h( )邛dh(v )/dv () = 152h(v )Mv2() = 0，r r，r rJ且枷匕)/叩(0) = 0 (见图2)。由于仞牧,)/存在且有界，该方程满足利普希茨条件。1.5.3欧拉一拉格朗日方程与性能函数对应的欧拉一拉格朗日方程为:V = bl |3口 p + R(b)u(t)(.a.=扁

39、F入=X +E vp j=iq = G (q)wq、ppt入1+ E ah(v )j=iOh(v) Ova t.jOvr, jOg(d ) OdOh(v ) OvFTrg (d )OdOpOvOpr ,jW = H-1S(Hw)w + T(t)O哄wOq0S(Hw)wOwTH-1,t XG (q)人 q q(14)g,h,d的大小如前文定义。式中，U(t)由方程(10)计算，T(t)由方程(11)计算。1.6结果翻滚目标对接下面的例子描述了一个对自主卫星服务来说很现实的任务，但超出了现有航天器轨迹规划程 序的能力。一个与ISS同轨道且在其后方约63km处的卫星发生故障并开始翻滚。将卫星视作直

40、径3m的球体，且两侧各有一沿y轴伸出体外7m的太阳能电池板。电池板宽3m。卫星上对接 点位于两太阳能电池板正中的星体上。卫星相对于惯性系的俯仰角速率为0，即相对ISS的俯仰 角速率约为1.14X 10-3rad/s （每92min旋转一周）。其偏航角速率约为5.2X 10-2rad/s （每2min 旋转一周）。要求设计一条服务星轨迹，轨迹的起点位于服务星在ISS上对接位置，距ISS中心 的距离为10m，且面向地球，终点是服务星要与翻滚卫星对接。翻滚卫星如图3所示。图3翻滚卫星模型（卫星俯仰角速率为1.14X10-3rad/s,偏航角速率为5.2X 10-2rad/s。初始轨道速度方向沿z轴正

41、向）假设机器人服务航天器的质量为2000kg，转动惯量1=32。 ml S 2。沿各轴 的最大作用力为1000N，最大作用力矩为2000Nm。力矩罚函数增益是10，时间增益Y是 0。如前例，采用千米、分和服务卫星质量单位构成非标准单位制对问题进行描述。而在表示结 果时，须转化回标准单位制。结果如图4所示。该轨迹需要81611际的推力，66.5Nms的力矩。总机动时间为2728s （45.46 min ）。整个轨道拦截轨迹见图4。卫星最后逼近段的时序显示在图5中。近距离观 察显示服务卫星在对接前“起旋”并持续了数分钟，以应对目标卫星的偏航转动。图4服务卫星对接问题轨迹。每隔约30s显示一次服务卫

42、星体坐标轴图5服务卫星逼近轨迹，每2s显示一次服务卫星体坐标轴该算例证明此轨迹规划程序能够解决需要质心运动与姿态运动高度协调的问题。在前面的例 子中，也可以获得满足bang-coast-bang （开-关-开）推力器约束的轨迹。最后，这个例子也说明 了该轨迹规划器能解决航天器飞行足够远距离的问题，这种情况下不能使用Hill方程或其他线性 化相对运动轨道方程。因此该轨迹规划器能够解决那些现有的轨迹规划器所不能解决的问题。1.7结论由以上算例可以看出，文中论述的变分轨迹规划算法能够满足引言中提出的要求。该算法具 有下列特点：能够在存在静止或运动障碍的情况下获得自由路径；能够在有引力场和其它环境力作

43、用的情况下寻找得到最小燃料消耗轨迹；能够处理大范围端点约束与时间尺度问题；考虑了现实的航天器动力学模型；考虑了现实的推力器，包括bang-coast-bang(开-关-开)推力器。同时，本文使用的动力学模型并非其它大部分轨迹规划法所采用的线性化模型。因此，相比 其它轨迹规划法，该算法能获得更大距离范围的最小燃料消耗轨迹。另外，变分轨迹规划法不局 限于质心运动轨迹规划；如最后一个例子所示，该方法能够解决姿轨紧耦合问题。综合以上特点, 可以看出，本文介绍的变分轨迹规划法对当前的轨迹规划有着深远的意义。参考文献Hohmann, W.(1925) Die erreichbarkeit der himm

44、elskarper (the attainability of celestial bodies), Munich, Germany.Hoelker, R. F., and Silber, P. S.(1961) “The bi-elliptical transfer between coplanar circular orbits.” Proc., Fourth AFBMD/STL Symp., Vol. 3, Pergamon, New York, 164-175.Tempelman, W.(1961) “Minimum-energy intercepts originating from

45、 a circular orbit.” J. Aerosp. Sci., 28, 924-929.Gedeon, G. S. (1958). “Orbital mechanics of satellites. Proc., American Astronautical Society, American Astronautical Society.Richards, A., Schouwenaars, T., How, J. P., and Feron, E.(2002). “Spacecraft trajectory planning with avoidance constraints u

46、sing mixedinteger linear programming. J. Guid. Control Dyn., 25_4_.Breger, L. S., and How, J. P.(2006). “Safe trajectories for autonomous rendezvous of spacecraft. Proc., 2006AIAA Guidance, Navigation and Control Conf., American Institute of Aeronautics and Astronautics.Jackson, M. C. (1994). “A six

47、 degree of freedom, plume-fuel optimal trajectory planner for spacecraft proximity operations using an anode search. MS thesis, Massachusetts Institute of Technology, Cambridge, Mass.Choset, H., Fernandez, J. J., and Kortenkamp, D.(1999). “Path planning and control for aercam, a free-flyinginspectio

48、n robot in space. Proc., 1999 IEEE Int. Conf, on Robotics and Automation (ICRA99), Institute of Electrical and Electronics Engineers _IEEE_, Detroit, 1396-1403.Choset, H., La Civita, M., and Park, J. C. (1999). “Path planning between two points for a robot experiencing localization error in known an

49、d unknown environments. Proc., Int. Conf, on Field and Service Robotics, Springer, Pittsburgh.Khatib, O. (1980). “Commande dynamique dans lespace operationnel des robots manipulateurs en presence dobstacles. Ph.D. thesis, Ecole Nationale Superieure de lAeronautique et de lEspace, Toulouse, France.Kh

50、atib, O. (1986). “Real time obstacle avoidance for manipulators and mobile robots. Int. J. Robot. Res., 5_1_.McInnes, C. R. (1999). “Autonomous path planning for on-orbit servicing vehicles.” Reducing space mission cost, American Institute of Aeronautics and Astronautics.McQuade, F., and McInnes, C.

51、 R. (1997). “Autonomous control for onorbit assembly using potential function methods.” Aeronaut. J., 101_1006_, 255-262.Roger, A. B., and McInnes, C. R. (2000). “Safety constrained free-flyer path planning at the international space station” J. Guid. Control Dyn., 23(6): 971-979.Garcia, I. M., and

52、How, J. P. (2005). “Trajectory optimization for satellite reconfiguration maneuvers with position and attitude constraints. Proc., 2005 IEEE American Control Conf., Institute of Electrical and Electronics Engineers _IEEE_, 889-895.Aoude, G. S., Garcia, I. M., and How, J. P. (2007). “Two-stage path p

53、lanning approach for designing multiple spacecraft reconfiguration maneuvers. Proc., 20th Int. Symp. on Space Flight Dynamics.Nah, R. S., and Vadali, S. R. (2001). “Fuel-optimal, low-thrust, threedimensional Earth-Mars trajectories. J. Guid. Control Dyn., 24.Raniera, C., and Ocampa, C. (2005). “Optimization of roundtrip, timeconstrained, finite burn trajectories via an indirect method. J. Guid. Control Dyn., 28.Betts, J. T. (1998). “Survey of numerical methods for trajectory optimization. J. Guid. Control Dyn., 21.Kirkpatrick