基于hp自适应伪谱法的组合动运载器发射窗口拓展能力分析

周宏宇，王小刚，崔乃刚，许河川

（1. 哈尔滨工业大学 航天学院，哈尔滨 150001；2. 哈尔滨建成机械厂，哈尔滨 150030）

发射窗口（Launch Window, LW）是指完成航天任务所允许的发射时间间隔[1]。传统运载火箭在发射后基本只能在射面内飞行，因而有效载荷在主动段结束后所进入的轨道平面主要取决于发射位置、发射时刻和发射方向，其中，发射位置和发射方向决定了轨道倾角，发射时刻决定了轨道升交点赤经。通常发射位置和发射方向很容易保证，但发射时刻会由于各种问题而造成延迟。一旦错过了预定发射时刻，载荷进入的实际轨道和预定轨道间将产生一定的夹角，即轨道偏差，此时若要消除入轨偏差，就需要消耗上面级载荷的燃料进行变轨，导致载荷运行寿命减低，或是重新制定发射计划，等待下一个发射窗口的到来。因此，为顺利完成任务，采用运载火箭的发射方式必须在规定的较小的时间区间内完成点火发射，因而发射窗口很窄，在灵活性上具有很大局限性。

基于组合动力的水平起飞飞行器（Combined Cycle Launch Vehicle, CCLV）具有巡航飞行的能力，其高机动性使其能够在巡航段进行横向机动[2]，改变飞行方向，进而改变入轨的轨道参数，即使错过了发射时刻也能通过横向机动进入预定轨道面，由此实现了对发射窗口的拓宽。

为分析CCLV对发射窗口的拓宽能力，首先要通过轨迹优化获知其最大横向机动能力。由于飞行过程中气动环境复杂、动力性能受飞行环境影响较大、运动耦合性强、约束条件严格、控制量可能存在尖点[3]，传统的打靶法、Gauss伪谱法和间接法等难以解决这类优化问题[4]。文献[5]通过简化运动模型，利用间接法实现了高超声速飞行器轨迹快速优化；文献[6]针对含有复杂约束条件的非线性最优控制问题，提出了改进的高斯伪谱法，但为提高精度必须增加节点，由此增大了计算量，甚至导致问题病态；文献[7]结合 hp自适应伪谱法和模式搜索法，解决了水平起飞组合动力飞行器的上升段轨迹优化问题。

针对传统发射方式发射窗口窄、发射时刻不灵活的问题，提出了采用组合动力飞行器作为运载器的新思路。建立了组合动力运载器的数学模型，包括动力学模型和组合动力发动机模型；给出了发射窗口的计算方法，推导了横向机动范围与发射窗口的关系；给出了hp自适应伪谱法的算法流程，并针对爬升段和巡航段分别选取控制参数和性能指标，建立了优化模型。最后通过仿真验证了优化算法的有效性，并根据优化结果分析了发射窗口能力。

1 组合动力运载器数学模型

组合动力运载器的动力系统由多种动力模式组成：首先采用涡轮基组合动力发动机（Turbo based Combined Cycle, TBCC）水平起飞、爬升和加速，当高度和速度达到一定值时，切换到火箭基组合动力发动机（Rocket Based Combined Cycle, RBCC），进一步加速和爬升。由于TBCC和RBCC发动机有着各自的优缺点和适用环境，组合动力可以将二者结合使用，从而实现单级入轨和重复使用。

1.1 组合动力运载器动力学模型

忽略地球扁率与旋转，飞行器质心运动方程如下：

式中：V为飞行速度大小，γ为飞行路径角，P为推力大小，α为攻角，g为重力加速度大小，σ为倾侧角，ψ为飞行航向角，r为飞行器到地心的距离，m为飞行器质量，m˙为发动机秒耗量，X和Y分别为气动阻力和升力。

式中：q为飞行动压，S为特征面积；xc和yc分别为阻力系数和升力系数，是攻角和马赫数的函数，可通过多项式拟合获得。

1.2 发动机模型

对于TBCC，推力和比冲是飞行高度和飞行马赫数的函数，为便于分析，通过多项式计算推力和比冲：

式中：P为发动机推力，sI为发动机比冲，H为飞行高度，Ma为飞行马赫数。

对于RBCC，火箭发动机秒耗量越大，则冲压发动机所占比例越小，对应推力越大、比冲越小[8]。为便于分析，将推力和比冲视作火箭发动机所占比重（记为K）的函数：

式中： Tmax为最大推力， Imax为最大比冲， K ∈[0,1]称为相对流量。本文取 Tmax=360 kN， Imax=3000 s。

2 发射窗口建模

2.1 发射窗口计算方法

在此定义飞行器在主动段结束后所进入的轨道为初始轨道，期望进入的轨道为目标轨道。下面推导无横向机动情况下飞行器的发射窗口宽度计算方法及横向机动拓宽发射窗口的原理。

若飞行器始终在初始射面内飞行，则当发射方位角和发射点纬度一定时，可根据式(5)计算初始轨道的轨道倾角：

式中：i为轨道倾角，A为发射方位角即发射方向，φ为发射点纬度。

如果同时确定了发射时刻，则由式(6)可求得初始轨道的升交点赤经：

式中：sΩ为初始轨道升交点赤经，Lθ为发射时刻对应的恒星时。

当初始轨道和目标轨道的轨道倾角不同时，两轨道间的夹角为

式中： i1为初始轨道的倾角， i2为目标轨道的倾角，ΔΩs为初始轨道和目标轨道间的升交点赤经差。

因此，若要使最终的轨道面偏差不超过Δα，则发射窗口为

式中：eω地球自转角速度。

图1 机动拓宽发射窗口示意图Fig.1 Sketch for extending the launch window

假设预定的发射时刻为 Tc，则在无横向机动的情况下，飞行器必须在时间区间[Tc-Δ t / 2 ,Tc+Δt /2]内点火发射才能保证入轨精度。

由上文叙述可知，发射时刻上的偏差最终将导致升交点赤经上的偏差，即轨道面偏差，而通过横向机动可以调整轨道面，进而弥补发射时刻上偏差，实现拓宽发射窗口。如图1所示，初始轨道和目标轨道的轨道倾角相同，但升交点赤经由于发射时刻不准确而有所差异，两个轨道面间存在夹角ΔΩs= Ωs1- Ωs2；通过横向机动，飞行器从初始轨道转移到目标轨道，在不改变轨道倾角的情况下改变了升交点赤经。因此，由发射时刻上的提前或推迟所造成的轨道面偏差可以通过横向机动修正回来，进而获得更大的发射窗口和发射灵活性。

2.2 发射窗口与横向机动的关系建立

图2给出了飞行器横向机动的俯视图，其中，ab为初始发射纵平面，飞行器经横向机动后到达c点。近似认为飞行器的转弯轨迹为一个圆弧，则圆弧上 a点和c点处的切线相交于d点，且有ad = c d。β为机动后的轨道面和初始轨道面间的夹角，且有

式中：cL和dL分别为机动过程中飞行器飞过的横程和纵程，cL即为横向机动距离。

图2 β计算方法示意图Fig.2 The calculation method of β

若 Lc是在考虑实际运动状态和约束条件下所获得的最大横程，则β便是对初始轨道面的最大调整角度。因此，当给定所允许的最大轨道面偏差Δα时，具有横向机动能力的飞行器可允许的初始轨道面偏差为 Δα ′ = Δα + β 。若认为实际轨道和目标轨道间不存在轨道倾角偏差，则由式(8)可得对应的发射窗口：

3 轨迹优化方法

3.1 hp自适应伪谱法

伪谱法在配点上利用P阶多项式近似状态量，为了提高近似精度，常采用增加节点（h方法[9]）或增加多项式阶次（p方法[10]）的方式。h方法通过细化网格降低差值多项式的阶次，但网格过密会增大非线性规划（NLP）问题的规模，影响求解速度；p方法即增加多项式阶次，具有设计变量少、初值易选取等优点，但同样会影响求解速度。

hp自适应伪谱法通过估计各时间区间上的状态量曲率来判断是采用h方法还是p方法来提高求解精度，其中采用h法还能够自行确定新增节点数量及位置。hp自适应伪谱法结合了h方法和p方法的优点，避免了不必要的网格细化因而减小了NLP的规模，利于提高求解速度[11]，同时具备求解精度高、收敛快、初值易选择、鲁棒性强等特点。

3.2 算法流程

1） 时间区间转换

将时间区间 t ∈[t0,tf]分为K个区间，并将每个区间转换至 τ ∈[-1 ,1]：

式中： t ∈ [ tk-1,tk]，k = 1 ,2,… ,K 。

2）计算离散误差和约束误差

将最优控制问题离散化，计算各区间内各状态量的微分方程离散误差 e(k)和过程约束误差 b(k)。

3）精度提升方法选取

如果第k个区间内的 e(k)和b(k)都小于设定阈值εd，则认为该区间能够精确近似原最优控制问题，优化结束；否则进一步优化网格。

如果 rk＜rkmax（rkmax为设定的正常数），则认为该区间内各状态量的曲率差别较大，一致性不足，需要增加网格内插值多项式的阶次来提高精度；否则增加网格内的配点数目。

4）多项式阶次的确定

设 D1、D0分别为更新前后插值多项式阶次，为 ε(k)=[e(k), b(k)]中的最大值，ceil()为正向取整函数，则：

式中：A为任意常整数。

5）新增节点数与节点位置的确定

更新后配点数的计算方法如下：

式中：Y为任意正常数。

新节点位置通过曲率密度函数来决定。曲率密度函数定义为

式中：c为设定常数，满足 F （ 1） = 1 。

3.3 优化模型建立

飞行器先后分别在 TBCC和 RBCC的作用下飞行，故将飞行分为TBCC段和RBCC段。考虑到不同飞行阶段的动力模式特点、任务需求和飞行环境，需针对各飞行阶段选取不同的优化变量、约束条件和性能指标。

TBCC段的主要任务是爬升和加速，在该阶段飞行器不进行横向机动。该阶段优化变量为攻角，性能指标为燃料最省，约束条件为α＜30°，终端约束为动力切换时的高度、速度和飞行路径角。

RBCC段的任务是继续加速爬升，并通过横向机动进行轨道调整。该阶段的优化变量为攻角、倾侧角和相对流量，性能指标为横程最大，约束条件为α＜8°、σ＜ 40°、0 ≤ K ≤ 1。该阶段对攻角的约束更加严格，这是由于RBCC发动机必须保证足够的进气量，过大的攻角会导致发动机熄火。终端约束为主动段结束时的飞行路径角、高度和飞行器质量，即在给定燃料质量下完成飞行。

因此，共有4个优化变量和2个性能指标，且各阶段的约束条件也有所不同。此外，在 TBCC段和RBCC段，还需对法向过载和动压进行约束。本文统一规定动压不超过90 kPa，法向过载不超过5.0。

4 仿真分析

4.1 仿真初始条件

首先采用优化方法进行轨迹优化。由式(3)可知，当高度大于25 km时，TBCC推力趋于平稳，推力值较小，且TBCC的工作马赫数一般不超过 5.0。故本文将动力切换时机设定为高度20 km，速度1100 m/s。飞行器起飞质量90 t，燃料质量80 t，入轨点高度为110 km。其它仿真条件见表1。

表1 各阶段初始及终端运动状态Tab.1 Initial and final states in different phases

此外，由式(10)可知，发射窗口与目标轨道倾角有关，因此需要分析同样机动能力在不同目标轨道倾角下对发射窗口的拓展能力。本文所选择的目标轨道倾角范围为 20°~80°，优化出最大横程后，利用 2.2节中的方法求出发射窗口。αΔ取1°，即允许1°的轨道面偏差。

4.2 仿真结果

TBCC段燃料消耗5 506.5 kg，飞行73.4 s；RBCC段最大横程达到441.1 km，对应纵程652.4 km，飞行200.1 s。仿真结果表明，本文所采用的优化算法能够满足终端约束及过程约束，并且仅在需要的地方增加了节点，由此减小了NLP问题的规模。

图3中攻角在动力切换时发生了突变，这是由于攻角作为控制量并未对其连续性进行约束，而法向过载是与攻角有关的量，故图7中法向过载也在同一时刻产生了突变。实际应用中，可以在TBCC段结束时先通过姿控将攻角调整至RBCC段的起始攻角，然后开启RBCC。

对于横程441.1 km、纵程652.4 km，各目标轨道倾角对应的发射窗口如表2所示。可以看出，目标轨道倾角越低，同样的横向机动距离获得的发射窗口更宽，这一点也可由式(10)推导得出；同时，组合动力飞行器的横向机动能力大大的拓宽了发射窗口。

此外，由式(5)可知，传统发射方式的初始轨道倾角不能小于发射点纬度，要进入低倾角轨道必须消耗上面的燃料级进行变轨。如果采用本文的方案，则能够通过横向机动在上升段进入低倾角轨道，避免消耗上面级燃料。

图3 攻角变化情况Fig.3 Time history of attack angle

图5 RBCC段倾侧角变化情况Fig.5 Time history of leaning angle in RBCC

图4 飞行路径角变化情况Fig.4 Time history of flight path angle

图6 RBCC段相对流量变化情况Fig.6 Time history of relative flow in RBCC

图7 动压和法向过载变化情况Fig.7 Time history of dynamic pressure and load

图8 飞行剖面Fig.8 The flight profile

表2 发射窗口计算结果Tab.2 The result of launch window

5 结 论

本文针对载人航天任务未来的发展趋势，设计了一种组合动力模式，并采用优化方法进行了轨迹优化设计，分析了横向机动能力和发射窗口拓展能力间的关系，具体如下：

1）采用基于组合动力的水平起飞运载器，可以克服传统发射方式窗口窄的问题，发射方式和时机更加灵活，还能在上升段进入轨道倾角小于发射点纬度的轨道。

2）hp自适应伪谱法能根据优化进程自适应调整离散网格，自行选择采用h方法或是p方法，在提高精度的同时降低了不必要的网格加密，以免NLP问题规模过大，因此提高了收敛速度，能够有效求解轨迹优化问题。

3）目标轨道倾角越低，发射窗口越易拓宽。

（

）：

[1] Deerwester J M, Mclaughlin J F, Wolfe J F. Earthdeparture plane change and launch window considerations for interplanetary missions[J]. Journal of Spacecraft &Rockets, 2015, 3(2): 169-174.

[2] Pu Z Q, Yuan R Y, Tan X M, et al. Active robust control of uncertainty and flexibility suppression for air-breathing hypersonic vehicles[J]. Aerospace Science & Technology,2015, 42: 429-441.

[3] Bowcutt K G.. Multidisciplinary optimization of air-breathing hypersonic vehicles[J]. Journal of Propulsion &Power, 2015, 17(6): 1184-1190.

[4] Takashima N, Lewis M J. Optimization of waverider based hypersonic cruise vehicles with off-design considerations[J]. Journal of Aircraft, 2015, 36(1): 235-245.

[5] Murillo O, Lu P. Fast ascent trajectory optimization for hypersonic air-breathing vehicles[C]//AIAA Guidance,Navigation, and Control Conference. 2013: 1-24.

[6] 孙勇. 基于改进 Gauss伪谱法的高超声速飞行器轨迹优化与制导[D]. 哈尔滨工业大学, 2012: 41-59.Sun Y. Trajectory optimization and guidance of hypersonic vehicle based on improved gauss pseudospectral method[M]. Harbin Institute of Technology, 2012: 41-59 .

[7] 周宏宇, 王小刚, 崔乃刚, 等. 基于hp自适应伪谱法的组合动力可重复使用运载器轨迹优化[J]. 中国惯性技术学报, 2016, 24(6): 832-837.Zhou H Y, Wang X G, Cui N G, et al. Trajectory optimazation of reusable vehicle with combined power based on hp adaptive pseudospectral algorithm[J]. Journal of Chinese Inertial Technology, 2016, 24(6): 832-837.

[8] Wang Y J, Jiang L I, Qin F, et al. Research on thermal adjustment method of RBCC in ramjet mode[J]. Journal of Solid Rocket Technology, 2016, 39(2): 151-159.

[9] Gandhi M. Trajectory optimization algorithm studies[J].Mathematics, 2015, 15(2): 267-281.

[10] Ross I M, Fahroo F. Pseudospectral knotting methods for solving non-smooth optimal control problems[J]. Journal of Guidance Control & Dynamics, 2004, 27(3): 397-405.

[11] 夏红伟, 李秋实, 李莉, 等. 基于hp自适应伪谱法的飞行器再入轨迹优化与制导[J]. 中国惯性技术学报, 2015,23(6): 818-823.Xia H W, Li Q S, Li L, et al. Trajectory optimization and guidance for reentry craft based on hp-adaptive pseudospectral method[J]. Journal of Chinese Inertial Technology, 2015, 23(6): 818-823.