CN111475767A - 考虑地球自转影响的最小能量弹道严格构造方法 - Google Patents

Abstract

本发明公开了一种考虑地球自转影响的最小能量弹道严格构造方法，采用了精确的地球物理模型(考虑了地球扁率和自转)，可在二体运动模型下准确构造两个三维大地点之间的最小能量弹道。另外，本发明通过深入挖掘地球自转对弹道构造的影响，从理论上找出最小能量弹道的严格解析解，并引入目标点和发射点地心距之比与1的差值ε和椭圆轨道近地点幅角偏置量Δω两个小量使得本方法具有普适性，公式简洁，为适用于任意两空间点间的最小能量自由弹道严格构造方法。

Description

考虑地球自转影响的最小能量弹道严格构造方法

技术领域

本发明涉及天体力学框架下解决弹道式导弹的弹道构造技术领域，具体而言涉及一种考虑地球自转影响的最小能量弹道严格构造方法。

背景技术

导弹的射程是指发射点与弹着点之间沿地球大圆测量的距离，是有效攻击目标的关键因素，也是最主要的战术指标。针对某一型号的导弹，最大射程取决于发动机性能、燃料性能、气动特性、弹道特性等多种因素，对其估算需要掌握一些难以确知的技术参数和物理条件，不能满足一般性的弹道研究和应用需求；当关机点位置确定时，自由飞行段的射程仅取决于导弹速度的大小和方向，基于椭圆弹道基本假设(张毅等.弹道导弹弹道学[M].1999)，射程与关机点能量参数和弹道倾角之间存在确切的函数关系。当关机点能量一定时，总能找到一个最佳速度方向使射程最大；反之，当射程一定时，也能找到一个最佳速度方向使得关机点能量最小(能量最小意味着速度最小，燃料最省)，由此确定的弹道称为最小能量弹道。远程导弹自由飞行段射程占全弹道射程约90％，发展自由飞行段最小能量弹道快速和有效的构造方法具有重要的应用价值，不仅可提供实验室中重要弹道的仿真环境，而且可辅助应用于弹道诸元设计。

当构造地球表面两点间的最小能量弹道时，地球自转引起目标点的惯性坐标时刻发生变化，脱离地表的导弹如何击中变化的目标点，是弹道构造的本质问题。文献(白鹤峰,吴瑞林.战术弹道导弹(TBM)弹道的构造方法[J].现代防御技术,1998(1):39-43)提出一种通过迭代飞行时间构造两地面点间最小能量弹道的方法，以下简称为算法1。算法1迭代求解的思想是解决目标点位置动态变化的合理思路，但每一步迭代过程中的最小能量弹道构造公式并未统筹考虑动力学运动规律和地球自转两个独立的过程，因此由算法1仅能得到近似解，而无法得出严格解，这一点已经通过专利号为201910941400.5的发明“一种指定发射仰角的自由弹道构造方法”中提及的遍历发射仰角的弹道构造计算得到验证，因此，算法1输出的弹道参数不是能量上的最优。此外，算法1的弹道构造公式以地球为均质正球体，以发射点和目标点均位于正球体的表面为基本前提，对应上述前提条件的弹道构造公式虽然简单，但会在几何学和动力学两个方面引起弹道的构造偏差，该偏差既随射程增加而增大，也随发射点和目标点地心距之差的增加而增大。综上可知，算法1既非最小能量弹道的严格构造方法，也缺乏普适性，不适用于发射点和目标点地心距不等的情况，而且当发射点位于两极或其上空时算法失效。

发明内容

本发明目的在于提供一种考虑地球自转影响的最小能量弹道严格构造方法，该方法采用了精确的地球物理模型(考虑了地球扁率和自转)，可在二体运动模型下准确构造两个三维大地点之间的最小能量弹道。

为达成上述目的，结合图1，本发明提出考虑地球自转影响的最小能量弹道严格构造方法，所述构造方法包括：

S1，对输入的已知量t₀、(L_A,B_A,H_A)和(L_B,B_B,H_B)进行预处理，计算得到发射点A和目标点B的地固直角坐标矢量

和

以及A点处B点方向与北极方向的水平夹角φ，发射点A和目标点B两点既不重合，也不位于地球直径的两端；

其中，t₀是发射时刻，(L_A,B_A,H_A)是发射点A点的大地经度、大地纬度和大地高，(L_B,B_B,H_B)是目标点B点的大地经度、大地纬度和大地高；

S2，设置迭代的初始状态：令飞行时间T等于零；

S3，计算发射点能量参数ν的初值ν⁽⁰⁾和变量ζ的初值ζ⁽⁰⁾，利用二体运动模型迭代求解发射点能量参数ν和变量ζ，直至求解得到的发射点能量参数ν和变量ζ满足以下条件：|ν^(k)-ν^(k-1)|<δ_ν，令ν＝ν^(k)，ζ＝ζ^(k)；

其中，ν是发射点能量参数，为发射点导弹动能的两倍与位能的比值，ζ＝cotθ，θ是惯性空间中发射速度与发射点地心向径方向的夹角，ν^(k)是第k次迭代求解得到的发射点能量取值，ζ^(k)是第k次迭代求解得到ζ的取值，δ_ν为第一指定小量；

S4，计算新的飞行时间T^*；

S5，判断|T-T^*|<δ_T是否成立，δ_T为第二指定小量，若成立，进入步骤S6，否则令T＝T^*，转入步骤S3；

S6，计算发射点轨道根数σ除半长径a以外的其它分量i,Ω,ξ,η,λ_A，以及导弹在地固坐标系和轨道坐标系中的发射速度

和

计算并输出在二体运动模型下构造的弹道设计参数v_p,v_r,h,

T^*,σ；

其中，v_p是地固坐标系下发射速度矢量的模，v_r是轨道坐标系下发射速度矢量的模，h是发射仰角，具体指导弹发射速度与当地水平面的夹角，T^*是迭代求出的新的飞行时间，

是发射速度相对于目标点方向在水平面上的偏离角，σ是发射点第一类无奇点根数。

当导弹射程超过1000公里，将主动段和再入段视为椭圆弹道的一部分是常见的处理手段，本发明将按照自由飞行段受力情况构造的全弹道称为自由弹道。为表述方便，用点A代表发射点，点B代表目标点，已知两点的三维大地坐标。本发明提供一种自A点至B点的最小能量自由弹道严格构造方法，该方法的实施涉及一系列的变量，主要变量的数学符号及含义列于表1，具体的实施流程如下：

表1 技术方案所涉及主要变量的符号及其含义

此外，本发明涉及到的地球物理常数和坐标转换矩阵有：

i.参考椭球体的赤道半径a_e＝6378136m；

ii.地心引力常数μ＝0.39860043770442×10¹⁵m³/s²；

iii.地球扁率f＝1/298.25781；

iv.地球子午圈偏心率：

v.地球自转速度

vi.地固坐标系至轨道坐标系的转换矩阵M，t时刻的M矩阵为R_Z[-θ_g(t)]，简写为M(t)，其中R_Z表示绕Z轴的旋转矩阵，θ_g为格林尼治恒星时角；

vii.地固坐标系至地平坐标系的转换矩阵W，其表达式为

其中L和B表示测站的大地经度和大地纬度，R_Y表示绕Y轴的旋转矩阵；

viii.地平坐标系至特定地平坐标系的旋转矩阵Q，其表达式为R_Z[φ]，其中φ为特定地平坐标系的X轴在地平坐标系中的方位角。

本发明提出的考虑地球自转影响的最小能量弹道严格构造方法包括以下步骤：

步骤一：数据的预处理，计算发射点A和目标点B的地固直角坐标矢量

和

以及A点处B点方向与北极方向的水平夹角φ；

步骤一由以下子步骤组成：

(1)根据公式组(1)将A、B两点的大地坐标由地理坐标形式(L,B,H)转换为空间直角坐标形式(X,Y,Z)，分别用矢量

和

表示。

公式组(1)中，N为卯酉圈曲率半径，

(2)根据公式(2)计算A点在轨道坐标系下的坐标矢量

(3)根据公式(3)计算A、B两点的地心纬度

和

(4)根据公式(4)计算B、A两点地心距之比与1的差值ε，并根据公式(5)计算ε₁：

(5)根据公式(6)计算地球自转对最小能量弹道构造施加影响项的小因子系数∈：

(6)根据公式组(7)-(8)计算A点处B点方向与北极方向的水平夹角φ：

公式组(7)中，q为A点天顶方向与B点天顶方向的夹角，其计算公式为：

特别地，当A点位于两极或其上空时，公式组(7)仍然成立，且退化为：

步骤二：设置迭代的初始状态：令飞行时间T等于零。

步骤三：利用二体运动模型迭代求解发射点能量参数ν和变量ζ。

步骤三由以下子步骤组成：

(1)根据公式(10)，结合飞行时间T计算B点在轨道坐标系下的坐标矢量

(2)根据公式(11)计算半射程角β。

由夹角公式得：

(3)根据公式组(12)计算发射点能量参数ν的初值ν⁽⁰⁾和变量ζ的初值ζ⁽⁰⁾：

其中，

式中，ε₁，η_τ和b_τ是为使公式(12)表达简洁引入的辅助量。

(4)结合A、B两点的大地坐标，飞行时间等量，利用公式(13)计算函数f(T,η_τ)：

自k＝1开始利用公式组(14)-(15)迭代计算ν^(k)和ζ^(k)，直至满足条件|ν^(k)-ν^(k-1)|<δ_ν时结束，δ_ν为一指定小量，然后令ν＝ν^(k)，ζ＝ζ^(k)，轨道偏心率e＝e^(k)进入步骤四：

ζ^(k)中各组成函数的表达式如下：

上式中：

步骤四：计算新的飞行时间T^*。

步骤四由以下子步骤组成：

(1)将步骤三得到的精确ζ代入公式(16)求近地点幅角偏置量Δω：

(2)根据公式组(17)计算A、B两点在椭圆轨道上的真近点角f_A和f_B：

(3)根据公式组(18)-(19)计算A、B两点在椭圆轨道上的平近点角M_A和M_B：

其中，

(4)根据公式组(20)-(21)计算导弹新的飞行时间T^*：

其中，

式中，a为轨道半长径，n为平均角速度。

步骤五：判断|T-T^*|<δ_T是否成立，δ_T为一指定小量，若成立，进入步骤六，否则令T＝T^*，转入步骤三。

步骤六：计算发射点轨道根数σ除半长径a以外的其它分量i,Ω,ξ,η,λ_A，以及导弹在地固坐标系和轨道坐标系中的发射速度

和

计算并输出在二体运动模型下构造的弹道设计参数v_p,v_r,h,

T^*,σ。

步骤六由以下子步骤组成：

(1)根据公式组(22)计算椭圆轨道的倾角i和升交点赤经Ω：

(2)根据公式组(23)计算A、B两点在椭圆轨道上的真纬度角u_A和u_B：

(3)根据公式(24)计算椭圆轨道的近地点幅角ω：

ω＝ω₀+Δω 公式(24)

其中，

(4)根据公式组(25)计算发射点轨道根数的分量ξ、η和λ_A：

(5)根据公式组(26)-(28)计算导弹在轨道坐标系下的发射速度v_r和地固坐标系下发射速度v_p。

导弹在轨道坐标系和地固坐标系下的发射速度矢量

和

满足以下公式：

则有：

(6)根据公式(29)-(30)计算特定地平坐标系下的发射方位角A^*和发射仰角h。

特定地平坐标系下的发射速度矢量记为

可经过如下的坐标旋转得出：

(7)根据公式(31)计算发射速度方向相对于目标点方向在水平面上的偏离角

输出基于二体运动模型构造的自发射点至目标点最小能量自由弹道的设计参数v_p,v_r,h,

T^*,σ。

以上6个步骤形成流程图见图1。

以上本发明的技术方案，与现有相比，其显著的有益效果在于：

(1)当给定发射时刻，发射点和目标点大地坐标时，该方法能够在考虑地球自转情况下，采用二体运动模型严格构造从发射点到目标点的最小能量弹道。

(2)该方法的发射点和目标点可以灵活选择，发射点可以是常规意义上的发射点，也可以设定为入轨点；同理，目标点也可以设定为再入点。

(3)该方法排除了数学上的计算奇点，对发射点坐标没有特别的要求，当发射点选定为两极或其上空时，本方法仍然适用，弹道设计参数中的发射速度矢量有明确的意义。

(4)弹道构造完成后，可提供各种形式的弹道设计参数，不仅包含发射点处速度的大小(相对于轨道坐标系以及相对于地固坐标系)、方向(在发射点处水平面内速度的方位角和仰角)，导弹的飞行时间，而且还提供了导弹在发射时刻相对于轨道坐标系的六个轨道根数，对所构造的最小能量弹道特征进行全面和精确的描述。

(5)本发明深入挖掘地球自转对弹道构造的影响，从理论上找出最小能量弹道的严格解析解，并引入目标点和发射点地心距之比与1的差值ε和椭圆轨道近地点幅角偏置量Δω两个小量使得本方法具有普适性，公式简洁，为适用于任意两空间点间的最小能量自由弹道严格构造方法。

应当理解，前述构思以及在下面更加详细地描述的额外构思的所有组合只要在这样的构思不相互矛盾的情况下都可以被视为本公开的发明主题的一部分。另外，所要求保护的主题的所有组合都被视为本公开的发明主题的一部分。

结合附图从下面的描述中可以更加全面地理解本发明教导的前述和其他方面、实施例和特征。本发明的其他附加方面例如示例性实施方式的特征和/或有益效果将在下面的描述中显见，或通过根据本发明教导的具体实施方式的实践中得知。

附图说明

附图不意在按比例绘制。在附图中，在各个图中示出的每个相同或近似相同的组成部分可以用相同的标号表示。为了清晰起见，在每个图中，并非每个组成部分均被标记。现在，将通过例子并参考附图来描述本发明的各个方面的实施例，其中：

图1是本发明的考虑地球自转影响的最小能量弹道严格构造方法的流程图。

图2是本发明的准地固坐标系的定义示意图。

图3是本发明的t时刻轨道坐标系的定义示意图，其中Υ为真春分点，

为历元平春分点。

图4是本发明的地平坐标系XYZ与特定地平坐标系X*Y*Z*的定义示意图。

图5是本发明的地心天球坐标系示意图。

图6是实施例一中遍历发射仰角方法获得发射仰角与发射速度的关系示意图。

图7是实施例二中遍历发射仰角方法获得发射仰角与发射速度的关系示意图。

具体实施方式

为了更了解本发明的技术内容，特举具体实施例并配合所附图式说明如下。

为便于描述，用点A代表发射点，点B代表目标点。下面将结合实施例对技术方案部分进行清晰、完整地描述。实施例中既用到了技术方案中定义的变量及其符号，还通过给变量赋予具体的数值得出了相应的结果。首先，介绍本发明使用的坐标系的定义及其转换方法，再结合实施例详细阐述技术方案。

一、坐标系统及其转换矩阵

● 地固坐标系

地固坐标系是固连在地球上与地球一起旋转的坐标系，可方便地描述地球表面点的空间位置，根据Z轴指向的不同分为协议地固坐标系和准地固坐标系两种，二者之间的差别(由极移引起)对地面点坐标的影响极小。图2为准地固坐标系的定义，Z轴指向瞬时北极，参考平面为瞬时赤道面；协议地固坐标系的Z轴指向协议地极，参考平面为与地心和协议地极连线正交的平面。在本发明的研究背景下，不需要考虑协议地极和瞬时地极的区别，统称为地固坐标系(或大地坐标系)，有地理坐标和空间直角坐标两种形式。

● 轨道坐标系

轨道坐标系(见图3)是一种混合地心坐标系，其参考平面为瞬时真赤道，X轴位于参考平面内指向历元(本发明选择J2000.0)平春分点，该点实为瞬时真赤道上的一个假想点，在真春分点以东(μ+Δμ)处，其中μ为赤经总岁差，Δμ为赤经章动，它们的计算公式参见文献(吴连大.人造卫星与空间碎片的轨道和探测[J].紫台专著,2011)。轨道坐标系一方面不会引起地球引力场位函数展开式球谐系数的变化，另一方面坐标系附加摄动很小，处理一般精度要求的问题时可当作惯性系，因此轨道坐标系是方便于人卫工作分析方法的坐标系统。

● 地平坐标系与特定地平坐标系

地平坐标系(见图4)以导弹发射点为原点，以过原点的水平面为参考平面，X轴在参考平面内指向北极。特定地平坐标系X^*轴在参考平面内指向目标点B。本技术方案中的发射方位角A^*是定义在特定地平坐标系中的。在本发明的研究背景下，可忽略大地水准面和参考椭球体切平面的差别，用参考椭球体的切平面代替水平面。

● 地固坐标系至轨道坐标系的转换矩阵M

准地固坐标系与轨道坐标系(J2000.0历元)仅在X轴方向相差格林尼治恒星时角θ_g，因此坐标转换矩阵M为：

其中，θ_g的计算公式为：

θ_g＝280°.460618375+360°.9856122882d

d＝MJD(UT1)-51544.5，为2000年1月1日12h UT1起算至导弹位置时刻t的累积天数。

考虑到矩阵的正交性，轨道坐标系至地固坐标系的转换矩阵为M^T。

● 地固坐标系至地平坐标系的转换矩阵W

不考虑地固坐标系和地平坐标系原点位置不同的情况下，二者的转换可通过2次旋转完成，转换矩阵为地面点大地经纬度的函数。设地面点的大地经纬度分别为L和B，则地固坐标系至地平坐标系的转换矩阵W为：

上式中，R_Y和R_Z分别表示绕Y轴和Z轴的旋转矩阵。考虑到矩阵的正交性，地平坐标系至地固坐标系的转换矩阵为W^T。

● 地平坐标系至特定地平坐标系的转换矩阵Q

地平坐标系和特定地平坐标系的原点和参考平面均相同，仅X轴的指向不同。在地平坐标系中，设目标点B的方位角为φ，则将地平坐标系绕Z轴逆时针旋转φ角可得到特定地平坐标系，因此二者的转换矩阵Q为：

同理，特定地平坐标系至地平坐标系的转换矩阵为Q^T。

二、结合实施例阐述技术方案：

实施例1：任意给定导弹的发射时刻t₀，发射点A和目标点B的大地坐标(具体参数见表2)，在二体运动模型下构造最小能量自由弹道，并与遍历发射仰角的弹道构造方法以及算法1的弹道构造结果进行对比。

表2.实施例1(发射时间：北京时2013年5月25日18时45分20秒)

注：迭代的收敛条件δ_ν和δ_T分别取1E-8和10μs。为应用算法1，本例特别设置A、B两点的地心距相等，但该设置并不是本专利算法的要求。

实施例2：任意给定导弹的发射时刻t₀和目标点B的大地坐标，发射点A位于北极(具体参数见表3)，在二体运动模型下构造最小能量自由弹道，并与遍历发射仰角弹道构造方法的结果进行对比。

表3.实施例2(发射时间：北京时2013年5月25日18时45分20秒)

注：迭代的收敛条件δ_ν和δ_T分别取1E-8和10μs。本例也特别设置A、B两点的地心距相等。

除上述已知条件外，本发明还用到了以下地球物理常数：

参考椭球体的赤道半径a_e＝6378136m。

地心引力常数μ＝0.39860043770442×10¹⁵m³/s²。

地球扁率f＝1/298.25781。

地球子午圈偏心率：

地球自转速度

分析以上实施例，当导弹射程超过1000公里，将主动段和再入段视为椭圆弹道的一部分是常见的处理手段，本发明将按照自由飞行段受力情况构造的全弹道称为自由弹道。最小能量自由弹道由于最具经济性，是弹道构造最重要的关注指标。专利号为201910941400.5的发明“一种指定发射仰角的自由弹道构造方法”中，通过遍历发射仰角可严格构造出两大地点之间的所有弹道，其中发射速度最小的弹道即为所有弹道中的最小能量弹道，但受限于遍历步长，利用该方法找出的最小能量弹道为最优值的近似。本发明深入挖掘地球自转对弹道构造的影响，从理论上找出最小能量弹道的严格解析解，并引入目标点和发射点地心距之比与1的差值ε和椭圆轨道近地点幅角偏置量Δω两个小量使得本方法具有普适性，公式简洁，为适用于任意两空间点间的最小能量自由弹道严格构造方法。结合实施例的具体流程如下：

步骤一：数据的预处理，计算发射点A和目标点B的地固直角坐标矢量

和

以及A点处B点方向与北极方向的水平夹角φ；

步骤一由以下子步骤组成：

1)根据公式组(1)将A、B两点的大地坐标由地理坐标形式(L,B,H)转换为空间直角坐标形式(X,Y,Z)，分别用矢量

和

表示：

公式组(1)中，N为卯酉圈曲率半径，

实施例2中，A点位于北极，其地理经度的定义不明确，在本专利中，通过指定A点地理经度为介于0和360度之间的任意值，使得弹道构造过程得以顺利执行，但不影响弹道构造结果。

2)根据公式(2)计算A点在轨道坐标系下的坐标矢量

3)根据公式(3)计算A、B两点的地心纬度

和

4)根据公式(4)计算B、A两点地心距之比与1的差值ε，并根据公式(5)计算ε₁：

A、B两点地心距通常不相等，且它们的差距远小于地球半径，因此ε为非零的一个小量。

5)根据公式(6)计算地球自转对最小能量弹道构造施加影响项的小因子系数∈：

6)根据公式组(7)-(9)计算A点处B点方向与北极方向的水平夹角φ：

公式组(7)中，q为A点天顶方向与B点天顶方向的夹角，其计算公式为：

特别地，当A点位于两极或其上空时，公式组(7)仍然成立，且退化为：

公式组(7)-(8)由图5中的球面三角形O-PZ_AZ_B，通过应用球面三角公式推导而出。忽略大地水准面和参考椭球的差别，边PZ_B为B点大地纬度的余角，等于

边PZ_A为A点大地纬度的余角，等于

∠Z_BPZ_A为过A、B两点的子午线的夹角，等于(L_A-L_B)；边Z_AZ_B为A点天顶方向与B点天顶方向的夹角，用q表示；∠Z_BZ_AP记为φ，为待求量。应用球面三角形正弦定理和五元素公式得到如下两式：

整理得φ的计算公式：

再通过球面三角形余弦定理得到q的表达式：

整理得：

q的取值限定为(0,π)，要求A、B两点既不能重合，也不能位于地球直径的两端。这是因为当目标点和发射点重合，不存在构造弹道的实际需求；当目标点和发射点位于地球直径两端，在最小发射能量的约束下，经过两点的椭圆弹道有无数条，不能给出唯一一组弹道设计参数。

步骤二：设置迭代的初始状态：令飞行时间T等于零。

步骤三：利用二体运动模型迭代求解发射点能量参数ν和变量ζ。

步骤三由以下子步骤组成：

1)根据公式组(10)，结合飞行时间T计算A、B两点在轨道坐标系下的坐标矢量

和

2)根据公式(11)计算半射程角β。

由夹角公式得：

3)根据公式组(12)计算发射点能量参数ν的初值ν⁽⁰⁾和变量ζ的初值ζ⁽⁰⁾：

其中，

式中，ε₁，η_τ和b_τ是为使公式(12)表达简洁引入的辅助量。

4)结合A、B两点的大地坐标，飞行时间等量，利用公式(13)计算函数f(T,η_τ)：

自k＝1开始利用公式组(14)-(15)迭代计算ν^(k)和ζ^(k)，直至满足收敛条件|ν^(k)-ν^(k-1)|<δ_ν时结束，δ_ν为一指定小量，然后令ν＝ν^(k)，ζ＝ζ^(k)，轨道偏心率e＝e^(k)进入步骤四。

ζ^(k)中各组成函数的表达式如下：

上式中，

公式组(12)-(15)的推导过程如下：

根据椭圆弹道理论(张毅等.弹道导弹弹道学[M].1999)，自由弹道射程与发射点能量参数、弹道倾角(发射点地心向径与发射速度夹角θ的余角)、发射点和目标点的地心距之间存在确切的函数关系，用本发明定义的变量ζ,χ,ν和ε表示为：

其中，

且有ζ＝cotθ；χ＝tanβ，θ和β分别为惯性空间中发射速度与发射点地心向径的夹角和半射程角。

根据椭圆弹道理论(张毅等.弹道导弹弹道学[M].1999)，自由弹道飞行时间T也与发射点参数存在确切的关系：

Eq.(2)中，

和ε为常量，

为ν和θ的函数(张毅等.弹道导弹弹道学[M].1999)，因此T也是发射点能量参数ν和θ的函数，用数学语言描述为T＝T(ν,θ)。

Eq.(1)和Eq.(2)基于动力学理论推导得出，并未与地球自转发生联系。从地球自转的角度考虑，当给定A、B两点的大地坐标，则半射程角β是地球自转速度ω_e和导弹飞行时间T的函数，可用如下的公式表述：

Eq.(3)是利用地心天球上的球面三角形余弦公式推导得出。导弹飞行T时长后，发射点和目标点在地心惯性系中的地心经纬度分别为(L_A，

)和(L_B+ω_eT，

)，它们与北极点P构成球面三角形，大圆弧

即为射程角2β，应用余弦定理得到Eq.(3)，推导方法与公式(8)相同。

综合Eq.(1)～Eq.(3)发现，半射程角β是飞行时间T的函数，而飞行时间T是ν和θ的函数，用数学语言描述为β＝β[T(ν,θ)]＝β(ν,θ)。令

在给定发射点和目标点三维大地坐标的情况下，ε为一常量，则Eq.(1)中的χ(等同于β)也是发射点能量参数ν和θ的函数，由Eq.(1)得到ν和θ隐函数关系：

F(ν,θ)＝0 Eq.(4)

最小能量弹道的定义为惯性空间中的发射速度最小，当发射点坐标确定时，等价于ν值达到最小，则最小能量弹道需满足以下条件：

Eq.(4)两端同时对θ求全导数有

结合Eq.(5)得到：

将

代入Eq.(6)整理得：

其中，

Eq.(8)中，第3式是4个导数的乘积，由公式χ＝tanβ，Eq.(3)，Eq.(2)、

和ζ＝cotθ求导或求偏导得出。

和

两项是将动力学规律和地球自转规律有机结合的部分，其内在含义为导弹飞行T时长后的惯性坐标恰好等于B点随地球自转T时长后的惯性坐标。

项求导较为简单；Eq.(2)关系式复杂，求

时可将Eq.(2)分成4部分单独求导，最后再合并，此处不再赘述。

将Eq.(8)中所有的导数或偏导数代入Eq.(7)，再代入Eq.(6)，整理得：

其中，

函数g(ν,ζ)中，ν₁＝1-ν；ν₂＝1+2ε-(1+ε)ν；ν₃＝2-ν。

引入辅助量b_τ和η_τ后，Eq.(1)变形为：

①在不考虑地球自转的情况下，即满足ω_e＝0，此时∈＝0，则Eq.(9)退化为：

将Eq.(12)代入Eq.(11)得：

当ω_e＝0，Eq.(3)退化为仅由A、B两点已知大地坐标决定的量，精确的β可直接解算出来，即η_τ是一个已知量。此时Eq.(13)是关于发射点能量参数ν的一元二次方程，方程的系数为已知，考虑能量参数ν＞0，其真解为：

Eq.(14)得到的ν代入Eq.(12)即可得到ζ，再将ν和ζ同时代入Eq.(2)，解算得出飞行时间T，至此解算过程结束。

②考虑地球自转，

当发射点或目标点其中之一位于两极或其上空时，因

或

此时∈＝0，Eq.(9)退化为Eq.(12)，这种情况下仍能采用Eq.(11)和Eq.(12)联合求解，由于β值也是确定的，与不考虑地球自转的情况相同。

③若考虑地球自转，且发射点和目标点均不在两极或其上空，此时∈≠0，则必须采用Eq.(9)和Eq.(11)严格求解。此时因飞行时间T未知，无法由Eq.(3)解算β，所以η_τ暂时未知。解决这一困难的唯一方式是通过迭代求解，先假定一个T值，如T＝0，利用T可由Eq.(3)解算β，从而得到η_τ的值，联立Eq.(9)和Eq.(11)理论上存在唯一一组解ν和ζ，利用它们可解算得到新的飞行时间T，重复以上过程，当前后两次T的差值满足指定的精度要求时，迭代求解过程结束。

针对确切T值的每一次迭代求解，由于∈的存在，联立Eq.(9)和Eq.(11)虽然在理论上能够求解，但只能得到一个关于单变量ν或ζ的非常复杂的超越方程，因此不能简单求解。注意到∈是一个无量纲量，并估算其量级为∈～10^-3，因此Eq.(9)右端第二项对弹道构造的影响是一个微小的扰动。据此，可先略去该扰动，由Eq.(14)和Eq.(12)得到ν和ζ一组相当精确的近似值，记为ν⁽⁰⁾和ζ⁽⁰⁾，然后以如下方式进行迭代：

从k＝1开始迭代，当满足|ν^(k)-ν^(k-1)|<δ_ν条件时，认为迭代收敛，此时即得到精确的ν和ζ。

以上弹道构造涉及的3种情况，第1和第2种是第3种情况的特例，在弹道构造过程中可以统一表述为公式(14)-公式(15)，当∈＝0时，在执行公式(14)-公式(15)时仅需一次迭代即可收敛得到准确的ν和ζ。理论上可以证明，仅当∈＝0和

同时成立时，本发明的弹道构造方法等价于算法1，因此算法1并不是严格意义上的考虑了地球自转影响的最小能量弹道构造方法，而只是一个近似方法，且不能适用于一般的应用场景。

步骤四：计算新的飞行时间T^*。

步骤四由以下子步骤组成：

1)将步骤三得到的精确ζ代入公式(16)求近地点幅角偏置量Δω：

2)根据公式组(17)计算A、B两点在椭圆轨道上的真近点角f_A和f_B。

当引入Δω后，A、B两点的真近点角计算公式为：

3)根据公式组(18)-(19)计算A、B两点在椭圆轨道上的平近点角M_A和M_B。

根据开普勒方程有以下公式成立：

其中，

4)根据公式组(20)-(21)计算导弹新的飞行时间T^*。

其中，

式中，a为轨道半长径，n为平均角速度。

步骤五：判断|T-T^*|<δ_T是否成立，δ_T为一指定小量，若成立，进入步骤六，否则令T＝T^*，转入步骤三。

步骤六：计算发射点轨道根数σ除半长径a以外的其它分量i,Ω,ξ,η,λ_A，以及导弹在地固坐标系和轨道坐标系中的发射速度

和

计算并输出在二体运动模型下构造的弹道设计参数v_p,v_r,h,

T^*,σ。

步骤六由以下子步骤组成：

1)根据公式组(22)计算椭圆轨道的倾角i和升交点赤经Ω。

式中，i∈(0,π)Ω∈[0,2π)。

2)根据公式组(23)计算A、B两点在椭圆轨道上的真纬度角u_A和u_B。

3)根据公式(24)计算椭圆轨道的近地点幅角ω。

由椭圆轨道的性质可知，弹道的远地点位于地球外部空间，介于A、B两点之间。当A、B两点的地心距相等时，远地点真纬度角等于A、B两点真纬度角的中值，且远地点真近点角为π，远地点真纬度角减去其真近点角得到轨道近地点幅角ω。然而A、B两点的地心距通常不相等，则近地点幅角ω与通过上述方法求得的值存在偏差，利用已求得的ω的偏置量Δω，则有：

ω＝ω₀+Δω ω∈[0,2π) 公式(24)

其中，

4)根据公式组(25)计算发射点轨道根数的分量ξ、η和λ_A：

5)根据公式组(26)-(28)计算导弹在轨道坐标系下的发射速度v_r和地固坐标系下发射速度v_p。

导弹在轨道坐标系和地固坐标系下的发射速度矢量

和

满足以下公式：

则有：

6)根据公式(29)-(30)计算特定地平坐标系下的发射方位角A^*和发射仰角h。

特定地平坐标系下的发射速度矢量记为

可经过如下的坐标旋转得出：

特别地，当A点位于两极，本技术方案赋予其经度为0～360度之间的任意值，纬度为±90度，利用其经纬度得到地固坐标系至地平坐标系的旋转矩阵W，W与利用公式(9)求得的地平坐标系至特定地平坐标系的旋转矩阵Q相乘得到仅与B点经度L_B相关的旋转矩阵，即理论上可证明QW与L_A无关，因此地固坐标系至特定地平坐标系存在唯一的转换关系。由此通过引入特定地平坐标系，达到了消除两极奇点的效果。

7)根据公式(31)计算发射速度相对于目标点方向在水平面上的偏离角

输出基于二体运动模型构造自发射点至目标点最小能量自由弹道的设计参数v_p,v_r,h,

T^*,σ。

三、弹道构造结果

利用本发明提出的技术方案，当给定发射点和目标点的三维大地坐标，可确定能量最优弹道对应的发射参数，表4为实施例1分别利用本发明技术方案和算法1构造最小能量弹道的输出参数，比较发现利用本发明构造弹道的发射速度小于算法1相应结果，表明存在发射速度更小的弹道，也说明利用算法1得不到能量最优解。与算法1相比，本发明构造弹道的惯性速度小约2.5米/秒，发射仰角低约1.3度，飞行时间短约50秒，二者的差是A、B两点的大地坐标的复杂函数，在各项条件不利的情况下，由算法1构造的最小能量弹道与最优解的差距更大。图6为基于实施例1的已知条件，采用专利号为201910941400.5的发明“一种指定发射仰角的自由弹道构造方法”中的技术方案获得的发射仰角与惯性系下的发射速度的关系，据图可知存在发射速度的最小值，受限于遍历步长最小值不能确知，但根据趋势线的走向可确定出最小发射速度小于6903.739米/秒，对应的发射仰角介于26度和27度之间，与本专利的计算结果相吻合，说明本技术方案构造两大地点间最小能量弹道的正确性。

表5为实施例2利用本发明技术方案构造最小能量弹道的输出参数，由于发射点位于北极，此时等价于地球自转速度为零的情况，理论上弹道构造结果应与算法1相同，然而由于算法1采用了北天东坐标系，当发射点位于两极或其邻近区域时，坐标系的X轴失去意义，无法提供有效的弹道设计参数，而本专利采用了特定地平坐标系，即便发射点位于两极或其上空，方法仍然适用，因此是无奇点的。图7为利用遍历发射仰角方法获得实施例2发射仰角与发射速度的关系，根据趋势线可判断速度极小值与本专利的结果吻合，进一步验证了当发射点位于极点或其上空时，本技术方案的有效性和正确性。

表4.实施例1最小能量弹道主要输出参数与算法1的对比

表5.实施例2最小能量弹道的全部输出参数

在本公开中参照附图来描述本发明的各方面，附图中示出了许多说明的实施例。本公开的实施例不必定义在包括本发明的所有方面。应当理解，上面介绍的多种构思和实施例，以及下面更加详细地描述的那些构思和实施方式可以以很多方式中任意一种来实施，这是因为本发明所公开的构思和实施例并不限于任何实施方式。另外，本发明公开的一些方面可以单独使用，或者与本发明公开的其他方面的任何适当组合来使用。

虽然本发明已以较佳实施例揭露如上，然其并非用以限定本发明。本发明所属技术领域中具有通常知识者，在不脱离本发明的精神和范围内，当可作各种的更动与润饰。因此，本发明的保护范围当视权利要求书所界定者为准。

Claims (8)

1.一种考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，所述构造方法包括：

S1，对输入的已知量t₀、(L_A，B_A，H_A)和(L_B，B_B，H_B)进行预处理，计算得到发射点A和目标点B的地固直角坐标矢量

和

以及A点处B点方向与北极方向的水平夹角φ，发射点A和目标点B两点既不重合，也不位于地球直径的两端；

其中，t₀是发射时刻，(L_A，B_A，H_A)是发射点A点的大地经度、大地纬度和大地高，(L_B，B_B，H_B)是目标点B点的大地经度、大地纬度和大地高；

S2，设置迭代的初始状态：令飞行时间T等于零；

S3，计算发射点能量参数ν的初值ν⁽⁰⁾和变量ζ的初值ζ⁽⁰⁾，利用二体运动模型迭代求解发射点能量参数ν和变量ζ，直至求解得到的发射点能量参数ν和变量ζ满足以下条件：|ν^(k)-ν^(k-1)|＜δ_ν，令ν＝ν^(k)，ζ＝ζ^(k)；

其中，ν是发射点能量参数，为发射点导弹动能的两倍与位能的比值，ζ＝cotθ，θ是惯性空间中发射速度与发射点地心向径方向的夹角，ν^(k)是第k次迭代求解得到的发射点能量取值，ζ^(k)是第k次迭代求解得到ζ的取值，δ_ν为第一指定小量；

S4，计算新的飞行时间T^*；

S5，判断|T-T^*|＜δ_T是否成立，δ_T为第二指定小量，若成立，进入步骤S6，否则令T＝T^*，转入步骤S3；

S6，计算发射点轨道根数σ除半长径a以外的其它分量i，Ω，ξ，η，λ_A，以及导弹在地固坐标系和轨道坐标系中的发射速度

和

计算并输出在二体运动模型下构造的弹道设计参数ν_p，ν_r，h，

T^*，σ；

其中，ν_p是地固坐标系下发射速度矢量的模，ν_r是轨道坐标系下发射速度矢量的模，h是发射仰角，具体指导弹发射速度与当地水平面的夹角，T^*是迭代求出的新的飞行时间，

是发射速度相对于目标点方向在水平面上的偏离角，σ是发射点第一类无奇点根数。

2.根据权利要求1所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，步骤S1中，所述计算得到发射点A和目标点B的地固直角坐标矢量

和

以及A点处B点方向与北极方向的水平夹角φ包括以下步骤：

S11，据公式组(1)将A、B两点的大地坐标由地理坐标形式(L，B，H)转换为空间直角坐标形式(X，Y，Z)，分别用矢量

和

表示；

其中，N为卯酉圈曲率半径，

S12，根据公式(2)计算A点在轨道坐标系下的坐标矢量

S13，根据公式(3)计算A、B两点的地心纬度

和

S14，根据公式(4)计算B、A两点地心距之比与1的差值ε，并根据公式(5)计算ε₁；

S15，根据公式(6)计算地球自转对最小能量弹道构造施加影响项的小因子系数∈；

S16，根据公式组(7)-(8)计算A点处B点方向与北极方向的水平夹角φ；

公式组(7)中，q为A点天顶方向与B点天顶方向的夹角，其计算公式为：

3.根据权利要求2所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，当A点位于两极或其上空时，公式组(7)仍然成立，且退化为：

4.根据权利要求1所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，步骤S3中，所述利用二体运动模型迭代求解发射点能量参数ν和变量ζ的过程包括以下步骤：

S31，根据公式(10)，结合飞行时间T计算B点在轨道坐标系下的坐标矢量

S32，根据公式(11)计算半射程角β；

由夹角公式得：

S33，根据公式组(12)计算发射点能量参数ν的初值ν⁽⁰⁾和变量ζ的初值ζ⁽⁰⁾；

其中，

式中，ε₁，η_τ和b_τ是为使公式(12)表达简洁引入的辅助量。

S34，结合A、B两点的大地坐标，飞行时间等量，利用公式(13)计算函数f(T，η_τ)：

S35，自k＝1开始利用公式组(14)-(15)迭代计算ν^(k)和ζ^(k)，直至满足条件|ν^(k)-ν^(k-1)|＜δ_ν时结束，δ_ν为第一指定小量，然后令ν＝ν^(k)，ζ＝ζ^(k)，轨道偏心率e＝e^(k)进入步骤S4；

ζ^(k)中各组成函数的表达式如下：

上式中：

5.根据权利要求1所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，步骤S4中，所述计算新的飞行时间T^*的过程包括以下步骤：

S41，将步骤S3中得到的精确ζ代入公式(16)求近地点幅角偏置量Δω：

S42，根据公式组(17)计算A、B两点在椭圆轨道上的真近点角f_A和f_B；

S43，根据公式组(18)-(19)计算A、B两点在椭圆轨道上的平近点角M_A和M_B；

其中，

S44，根据公式组(20)-(21)计算导弹新的飞行时间T^*；

其中，

式中，a为轨道半长径，n为平均角速度。

6.根据权利要求1所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，步骤S6中，所述计算并输出在二体运动模型下构造的弹道设计参数ν_p，ν_r，h，

T^*，σ的过程包括以下步骤：

S61，根据公式组(22)计算椭圆轨道的倾角i和升交点赤经Ω；

S62，根据公式组(23)计算A、B两点在椭圆轨道上的真纬度角u_A和u_B；

S63，根据公式(24)计算椭圆轨道的近地点幅角ω；

ω＝ω₀+Δω (24)

其中，

S64，根据公式组(25)计算发射点轨道根数的分量ξ、η和λ_A；

S65，根据公式组(26)-(28)计算导弹在轨道坐标系下的发射速度ν_r和地固坐标系下发射速度ν_p；

导弹在轨道坐标系和地固坐标系下的发射速度矢量

和

满足以下公式：

则有：

S66，根据公式(29)-(30)计算特定地平坐标系下的发射方位角A^*和发射仰角h；

特定地平坐标系下的发射速度矢量记为

可经过如下的坐标旋转得出：

S67，根据公式(31)计算发射速度方向相对于目标点方向在水平面上的偏离角

S68，输出基于二体运动模型构造的自发射点至目标点最小能量自由弹道的设计参数ν_p，ν_r，h，

T*.σ。

7.根据权利要求1所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，所述发射点包括常规发射点和入轨点；所述目标点包括常规目标点和再入点。

8.根据权利要求1所述的考虑地球自转影响的最小能量弹道严格构造方法，其特征在于，所述构造方法适用的坐标系包括地固坐标系、轨道坐标系、地平坐标系。

Images