CN106599410B - 一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法 - Google Patents

一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法 Download PDF

Info

Publication number
CN106599410B
CN106599410B CN201611075769.5A CN201611075769A CN106599410B CN 106599410 B CN106599410 B CN 106599410B CN 201611075769 A CN201611075769 A CN 201611075769A CN 106599410 B CN106599410 B CN 106599410B
Authority
CN
China
Prior art keywords
msub
mrow
mtd
mtr
mtable
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201611075769.5A
Other languages
English (en)
Other versions
CN106599410A (zh
Inventor
韦常柱
常晓华
李源
丰海
何飞毅
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Heilongjiang Industrial Technology Research Institute Asset Management Co ltd
Original Assignee
Harbin Institute of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology filed Critical Harbin Institute of Technology
Priority to CN201611075769.5A priority Critical patent/CN106599410B/zh
Publication of CN106599410A publication Critical patent/CN106599410A/zh
Application granted granted Critical
Publication of CN106599410B publication Critical patent/CN106599410B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]

Abstract

本发明提供了一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法。该系统采用模块化思想构建，计算速度快精度高，能够适应不同任务、不同弹道形态的扰动引力影响特性分析要求，提高了分析结果的可信度。本发明的优点在于提供了良好的交互界面、集成了通用的导弹动力学模型库和扰动引力计算方法库，便于系统的操作和功能扩充。

Description

一种多赋值法的扰动引力场对不同形态弹道影响特性分析系 统及方法
技术领域
本发明涉及一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法。
背景技术
引起导弹落点偏差的因素主要包括制导工具误差和制导方法误差。制导工具误差是制导设备，如惯性平台、加速度表、陀螺仪以及计算机等性能不完善引起的落点偏差，约占总误差的70％～80％。制导方法误差是由于制导方法不完善引起的落点偏差。随着制导设备性能的改进，制导方法误差的影响日渐突出。在引起制导方法误差的因素中，扰动引力场是主要误差源，对制导工具误差也有一定的影响。对于射程10000km左右的远程弹道导弹，全弹道上作用的扰动引力造成的落点偏差最大可达千米量级。同时，对于长时间低空机动飞行导弹、多头分导导弹等，其受扰动引力影响时间长，误差累积耦合作用明显。因此，深入研究扰动引力场对不同形态弹道及落点偏差的影响特性，开发集成简易直观、操作便捷的仿真分析平台具有重要意义。
确定扰动引力场的关键在于求解关于扰动引力位的外部边值问题。扰动引力赋值方法大致归纳为两类：模型逼近和算法逼近。模型逼近主要有以Stokes理论为代表的大地水准面边值问题解(直接积分法)、广义延拓法和覆盖层法等；算法逼近主要有点质量法、最小二乘配置法、有限元逼近法和球谐函数法等。几种主要赋值方法的特点见下表：
正因为有众多扰动引力的赋值方法，而且很多方法在理论方面也已经发展成熟，理论体系严密，所以在计算扰动引力时，可供选择的方法很多。如何选取最合适、高效的方法来计算某种形态的弹道，将各种计算方法整合在一起并突出每种方法的功能优势的集成平台可以很好地解决上述问题。对于扰动引力影响特性的分析方法，除了传统的直接积分求差法外，利用状态空间摄动理论推导的摄动法及误差传播分析法可以弥补直接积分求差法的不足，为影响特性的分析带来方便。而将扰动引力计算与弹道分析结合在一起，可以更加多元地进行弹道分析，增加结果间的对比性，为相关研发单位提供技术依据和数据支持。
发明内容
本发明的目的是为了实现导弹诸元射击、快速机动发射的能力，提高导弹精确打击能力，并获得扰动引力场对不同形态弹道及落点偏差的影响特性，由此本发明提供了一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法。
上述的发明目的是通过以下技术方案实现的：
一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统，包括仿真平台主界面模块，
所述的仿真平台主界面模块包括多种扰动引力赋值算法选择子模块、弹道形态选择子模块、扰动引力对弹道影响特性分析子模块、弹道仿真子模块和数据处理与结果显示子模块；
所述的多种扰动引力赋值算法选择子模块确定引力计算底层的数学形式；
所述的弹道形态选择子模块提供多种弹道模式，通过外部输入程序角完成指定弹道的生成；
所述的扰动引力对弹道影响特性分析子模块是引力影响特性分析的基础，提供三种不同的引力影响特性分析方法，用户根据实际需求选择分析方法，该三种分析方法可单独使用或任意组合使用；
所述的弹道仿真子模块根据底层程序代码及外部输入的程序角用于弹道仿真，得出扰动引力作用下的弹道计算结果；
所述的数据处理与结果显示子模块对得出的弹道计算结果进行数据读取、处理、查询、比较以及分析。
进一步的，所述的仿真平台主界面模块是面向用户的，仿真平台主界面模块有专门的帮助子模块的选项菜单和按钮，方便用户随时切换。
进一步的，所述的多种扰动引力赋值算法选择子模块从有限元逼近法、广义延拓法及球谐函数法中任意选择一种方法，三种方法是等价的，均可以提供对应模型在所需坐标点处扰动引力三分量，用户依据不同的任务要求选择合适的扰动引力计算方法。
进一步的，所述的弹道形态选择子模块提供多种弹道模式，通过外部输入程序角完成指定弹道的生成，如果没有设计好的程序角，则直接使用程序内置的多种典型弹道的程序角数据。
进一步的，所述的扰动引力对弹道影响特性分析子模块对积分求差法、摄动法和误差传播分析方法的优缺点和适应性给出提示，用户根据提示及实际需求选择相应的分析方法。
进一步的，所述的弹道仿真子模块用于弹道积分计算，得出扰动引力作用下的弹道计算数据，用户根据需要选取坐标系，导弹的总体参数由用户更新相关设置文件，其中相关设置文件包括典型导弹动力学模型、仿真步长、关机条件、落点偏差显示类型，其中典型导弹动力学模型包括气动参数、发动机参数以及不同坐标系间的转换模型。
进一步的，所述的数据处理与结果显示子模块对得出的弹道计算结果进行数据读取、处理、查询、比较以及分析，并将结果存入数据库与已经存在的仿真数据进行对照，其中，数据的读取遵循数据交互管理模块的接口定义，通过数据传输机制获得数据信息；数据查询功能依照数据管理系统进行操作，结果分析中还包含了数据处理、比较以及分析模块。
本发明还提供一种多赋值法的扰动引力场对不同形态弹道影响特性分析方法，包括如下步骤：
步骤一、在多赋值法的扰动引力场对不同形态弹道影响特性分析系统运行之前，进行Teechart绘制曲线插件的注册；
步骤二、Teechart绘制曲线插件注册成功之后，在VC6.0以上版本环境中打开多赋值法的扰动引力场对不同形态弹道影响特性分析系统软件，编译.exe文件通过后进入仿真平台主界面模块；
步骤三、在仿真平台主界面模块中的多种扰动引力赋值算法选择子模块对引力赋值方法进行选择，并对相应方法进行设置；其中参数设置包括：格网数划分选择，拟合点个数选择，扰动引力截断阶次，引力场位系数模型选择；
步骤四、在弹道形态选择子模块中设置导弹发射方式、发射点坐标、发射方位角以及射程信息，如果不使用默认数据，则可以自行输入指定弹道程序俯仰角；
步骤五、在扰动引力对弹道影响特性分析子模块中对分析方法进行设置，如果选择弹道积分求差分析方法，则设置标准弹道正常引力模型截断阶次；如果选择摄动法，则设置弹道平均角速度和变换步长；如果选择误差传播分析法，则设置权系数表示方式；在弹道仿真子模块中对仿真步长进行设置；
步骤六、利用数据处理与结果显示子模块对弹道过程数据读取、处理、查询、比较以及分析，并且弹道轨迹仿真数据结果、仿真过程中数据实时更新；到此分析方法结束。
进一步的，步骤三中多种扰动引力赋值算法包括：
有限元逼近法：
有限元逼近法的基本思路如下：(1)将所考虑的空域Ω用研究选定的网格划分为若干个单元Ωe，Ωe可由点质量法计算出扰动引力的量级确定或指定，Ωe的大小要满足扰动引力的计算精度；(2)将单元Ωe内的扰动引力构造成节点的多项式逼近函数，节点处的扰动引力由点质量法计算并存储；(3)为计算Ω中任一点A的扰动引力δ，需先判断计算点A所在的单元Ωe，再由点A位置与该单元节点值，用逼近函数计算δ；
记所考虑的空域Ω为一球壳，大小可由地固球坐标的取值区间表示；考虑到导弹的飞行空域，将Ω取为
其中，下标“0”、“f”分别表示起始点与终值点；
将Ω用所选择的网格分割为具有形状和大小的单元Ωe，单元的分割有以下要求：
(1)两个相邻单元Ωe在公共边界上是相容的，其顶点不能作为相邻单元边上的内点，即保持连续性；
(2)网格尽量规则，以减少计算量，减轻管理上的难度；
(3)靠近发射点的低空区域，对弹道产生影响的扰动引力集中、变化大，网格要密些，对应的单元Ωe体积小；反之，对远离发射点的高空区域，网格可疏些，Ωe体积大；
根据以上规则，取单元Ωe的形状为由地固球坐标截得的六面体，各面分别为：半径为r1与r2的球面，经度为λ1与λ2的子午面，纬度为与的纬圈，设r1＜r2，λ1＜λ2，沿径向看，Ωe的厚度为Δr，低空的Δr取值小；在球面上看，球面元方块大小为Δλ、为球面元覆盖的经度、纬度，离发射点近的区域，Δλ、取值小，为规整起见，令每一单元显然，Ωe可由8个节点的坐标表示，为了方便表示及计算，引入局部坐标系表示单元节点；
单元Ωe的局部曲线坐标系由半径rP＝r1+Δr/2的球面、经度λP＝λ1+Δλ/2的子午面、纬度的纬圈的交线组成；原点P为三交线的交点，局部坐标ξ、η、ζ分别沿点P的径向、纬圈与子午线方向；
原点P的球坐标为局部坐标为P(0,0,0)，单元内变点的局部坐标A(ξ,η,ζ)为：
单元顶点的局部坐标Ai(ξi,ηi,ζi)为：
记扰动引力δ在天、东、北方向上的分量为简记为因此δα是球坐标的函数，当将转换成局部坐标(ξ,η,ζ)后，δα可表示成(ξ,η,ζ)的函数；
在单元Ωe内，由于只取了8个节点的信息，因此逼近函数最多可以有8个待定系数；取逼近函数为如下的8次多项式：
δα＝a1+a2ξ+a3η+a4ζ+a5ξη+a6ξζ+a7ηζ+a8ξηζ (3)
节点的δαi由球谐函数法计算得到，记
则待定系数以向量表示为：
a＝H-1Z (5)
引入型函数Ni，满足
其中，
则变点A扰动引力分量计算式为：
经推导，Ni(A,Ai)可用下面的公式计算
广义延拓法：
考虑如下三维数值逼近问题；已知函数u(x,y,z):R3→R，在空域上的一组离散数据
{ui|ui＝u(xi,yi,zi),(xi,yi,zi)∈Ω,i＝1,2,…,n}
在空域Ω上构造u的一个近似函数U:Ω→R满足U(xi)＝ui(i＝1,2,…,n)；
采用分块逼近的方法求解，首先将区域进行划分，得m个互不重叠的子区域：
Ω＝Ω1∪Ω2∪…∪Ωm
设第e单元Ωe中包含r个插值节点；将单元Ωe与临近几个子域结合起来形成延拓域Ωe′，其中含s个节点，且s＞r，于是有：
在单元Ωe上，取三元多项式类
{gj(x,y,z)}＝{1,x,y,z,x2,y2,z2,xy,xz,yz,x3,y3,z3,
x2y,x2z,xy2,xz2,yz2,y2z,xyz,…}
的前t项为插值基函数，且r＜t＜s，即令：
其中，a1，a2，…，at为待定系数，可由下述问题解出：
显然，若r＝t＝s，即为通常意义下的分片插值模型，该模型为上节有限元模型；若r＝0,t＜s即为通常意义下的分片拟合模型；故广义延拓逼近是集插值、拟合之长的高精度分片逼近方法。
将空间域作立方体剖分，每个立方体子域Ωe上的节点编号依次为0、1、2、3、4、5、6和7，其延拓域Ωe′上有s个节点，s＞8；对于给定的插值基函数，可以给出待定系数求解的具体表达式；取r＝8，s＝32，t＝20；
将式(10)用矩阵形式表示，令
G＝{gj(xi,yi,zi)}ij,i＝9,10,...,32；j＝1,2,...,t
GI＝{gj(xi,yi,zi)}ij,i＝1,2,...,8；j＝1,2,...,t
u＝[u9 u10 … u32]T
uI＝[u1 u2 … u8]T
a＝[a1 a2 … a20]T
则式可表示为：
引入拉格朗日乘子λ＝[λ1 λ2 … λ8]T，则：
L(a,λ)＝(Ga-u)T(Ga-u)+2(G1a-u1)λ (12)
根据优化原理，待定系数a可通过代数方程
求解，其中各矩阵由如下公式表示，
F0＝GTu,F1＝u1
A＝GTG,C＝G1
A为方阵且可逆，结合分块矩阵求逆公式，最后有：
其中，D11和D12分别是A和C的逆；
令P(x,y,z)＝{gj(x,y,z)|j＝1,2,…,t}，则扰动引力分量
其中，型函数为：
N(x,y,z)＝P(x,y,z)[D11GT D12] (16)
球谐函数法：
以分量形式表示的扰动引力：
其中，μ为地球引力常数，ae为正常地球赤道半径，r为计算点至地心的距离，为计算点地心纬度，λ为计算点经度，为规格化的勒让德伴随函数，N为球谐函数最高阶，为位系数，n为扰动引力的截断阶，m≠n时表示位系数为田谐系数，m＝n时表示位系数为扇谐系数，d为微分号。
本发明的目的在于提供了一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法，为用户提供选择引力赋值模型的功能，用户根据需要选择或添加不同的引力模型和影响特性分析模型进行弹道计算仿真，可以实现对不同形态弹道的数值仿真，采用模块化思想构建计算速度快精度高，能够适应不同发射任务、不同引力影响条件，提高了影响特性分析的应用范围和通用性。本发明提供了一个弹道数据库，便于多种弹道的比较分析及系统功能的扩充，系统将目前成熟的引力计算方法集成在一个计算方法库中，采用模块化构建、面向对象设计，方便不同专业水平的用户使用，另外仿真平台还拥有一个覆盖了航空航天领域绝大多数类型导弹的动力学模型库，如巡航导弹、防空导弹、弹道导弹、运载火箭等。底层算法采用VC++语言编写，具有人机交互界面，计算速度快，易于不同平台间的移植，便于各种用户的使用。
附图说明
图1是多赋值法的扰动引力场对不同形态弹道影响特性分析系统Teechart的插件注册图；
图2是多赋值法的扰动引力场对不同形态弹道影响特性分析系统仿真平台主界面模块图；
图3是多赋值法的扰动引力场对不同形态弹道影响特性分析系统扰动引力赋值方法选择模块图。
具体实施方式
下面将结合本发明实施例中的附图对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
参考图2，本发明提供一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统，包括仿真平台主界面模块，所述的仿真平台主界面模块包括多种扰动引力赋值算法选择子模块、弹道形态选择子模块、扰动引力对弹道影响特性分析子模块、弹道仿真子模块和数据处理与结果显示子模块；
所述的多种扰动引力赋值算法选择子模块确定引力计算底层的数学形式；
所述的弹道形态选择子模块提供多种弹道模式，通过外部输入程序角完成指定弹道的生成；
所述的扰动引力对弹道影响特性分析子模块是引力影响特性分析的基础，提供三种不同的引力影响特性分析方法，用户根据实际需求选择分析方法，该三种分析方法可单独使用或任意组合使用；
所述的弹道仿真子模块根据底层程序代码及外部输入的程序角用于弹道仿真，得出扰动引力作用下的弹道计算结果；
所述的数据处理与结果显示子模块对得出的弹道计算结果进行数据读取、处理、查询、比较以及分析。
其中多种扰动引力赋值算法选择子模块、弹道形态选择子模块、扰动引力对弹道影响特性分析子模块，以上三种模块可以用其他模块替换，上述不同模块之间既有各自的功能又互相影响；
进一步的，所述的仿真平台主界面模块是面向用户的，仿真平台主界面模块有专门的帮助子模块的选项菜单和按钮，方便用户随时切换。
所述的仿真平台主界面模块是面向用户的，包括多种扰动引力赋值算法选择子模块、弹道形态选择子模块、扰动引力对弹道影响特性分析子模块、弹道仿真子模块和数据处理与结果显示子模块；仿真平台主界面模块有专门的帮助子模块的选项菜单和按钮，方便用户随时切换；
参考图3，所述的多种扰动引力赋值算法选择子模块从有限元逼近法、广义延拓法及球谐函数法中任意选择一种方法，三种方法是等价的，均可以提供对应模型在所需坐标点处扰动引力三分量，用户依据不同的任务要求选择合适的扰动引力计算方法。
所述的多种扰动引力赋值算法选择子模块确定引力计算底层的数学形式。提供的算法有有限元逼近法、广义延拓法及球谐函数法。三种方法是等价的，均可以提供对应模型在所需坐标点处扰动引力三分量。所以依据不同的任务要求，可以选择合适的扰动引力计算方法以提高效率。
进一步的，所述的弹道形态选择子模块提供多种弹道模式，通过外部输入程序角完成指定弹道的生成，如果没有设计好的程序角，则直接使用程序内置的多种典型弹道的程序角数据。
进一步的，所述的扰动引力对弹道影响特性分析子模块对积分求差法、摄动法和误差传播分析方法的优缺点和适应性给出提示，用户根据提示及实际需求选择相应的分析方法。
所述的扰动引力对弹道影响特性分析子模块是引力影响特性分析的基础，提供三种不同的引力影响特性分析方法，使用积分求差法来获得精度最高的影响特性分析结果，使用摄动法可以大大提高分析效率，使用误差传播分析方法可以观测飞行过程中扰动引力影响程度的变化情况，一般与前两种分析方法结合使用。
进一步的，所述的弹道仿真子模块用于弹道积分计算，得出扰动引力作用下的弹道计算数据，用户根据需要选取坐标系，导弹的总体参数由用户更新相关设置文件，其中相关设置文件包括典型导弹动力学模型、仿真步长、关机条件、落点偏差显示类型，其中典型导弹动力学模型包括气动参数、发动机参数以及不同坐标系间的转换模型。
进一步的，所述的数据处理与结果显示子模块对得出的弹道计算结果进行数据读取、处理、查询、比较以及分析，并将结果存入数据库与已经存在的仿真数据进行对照，其中，数据的读取遵循数据交互管理模块的接口定义，通过数据传输机制获得数据信息；数据查询功能依照数据管理系统进行操作，结果分析中还包含了数据处理、比较以及分析模块。
本发明还提供了一种多赋值法的扰动引力场对不同形态弹道影响特性分析方法，多赋值法的扰动引力场对不同形态弹道影响特性分析方法包括如下步骤：
步骤一、在多赋值法的扰动引力场对不同形态弹道影响特性分析系统运行之前，进行Teechart绘制曲线插件的注册；如图1所示；
步骤二、Teechart绘制曲线插件注册成功之后，在VC6.0以上版本环境中打开多赋值法的扰动引力场对不同形态弹道影响特性分析系统软件，编译.exe文件通过后进入仿真平台主界面模块；
步骤三、在仿真平台主界面模块中的多种扰动引力赋值算法选择子模块对引力赋值方法进行选择，并对相应方法进行设置；其中参数设置包括：格网数划分选择，拟合点个数选择，扰动引力截断阶次，引力场位系数模型选择；
步骤四、在弹道形态选择子模块中设置导弹发射方式、发射点坐标、发射方位角以及射程信息，如果不使用默认数据，则可以自行输入指定弹道程序俯仰角；
步骤五、在扰动引力对弹道影响特性分析子模块中对分析方法进行设置，如果选择弹道积分求差分析方法，则设置标准弹道正常引力模型截断阶次；如果选择摄动法，则设置弹道平均角速度和变换步长；如果选择误差传播分析法，则设置权系数表示方式；在弹道仿真子模块中对仿真步长进行设置；
步骤六、利用数据处理与结果显示子模块对弹道过程数据读取、处理、查询、比较以及分析，并且弹道轨迹仿真数据结果、仿真过程中数据实时更新；到此分析方法结束。
进一步的，步骤三中多种扰动引力赋值算法包括：
有限元逼近法：
有限元逼近法的基本思路如下：(1)将所考虑的空域Ω用研究选定的网格划分为若干个单元Ωe，Ωe可由点质量法计算出扰动引力的量级确定或指定，Ωe的大小要满足扰动引力的计算精度；(2)将单元Ωe内的扰动引力构造成节点的多项式逼近函数，节点处的扰动引力由点质量法计算并存储；(3)为计算Ω中任一点A的扰动引力δ，需先判断计算点A所在的单元Ωe，再由点A位置与该单元节点值，用逼近函数计算δ；
记所考虑的空域Ω为一球壳，大小可由地固球坐标的取值区间表示；考虑到导弹的飞行空域，将Ω取为
其中，下标“0”、“f”分别表示起始点与终值点；
将Ω用所选择的网格分割为具有形状和大小的单元Ωe，单元的分割有以下要求：
(1)两个相邻单元Ωe在公共边界上是相容的，其顶点不能作为相邻单元边上的内点，即保持连续性；
(2)网格尽量规则，以减少计算量，减轻管理上的难度；
(3)靠近发射点的低空区域，对弹道产生影响的扰动引力集中、变化大，网格要密些，对应的单元Ωe体积小；反之，对远离发射点的高空区域，网格可疏些，Ωe体积大；
根据以上规则，取单元Ωe的形状为由地固球坐标截得的六面体，各面分别为：半径为r1与r2的球面，经度为λ1与λ2的子午面，纬度为与的纬圈，设r1＜r2，λ1＜λ2，沿径向看，Ωe的厚度为Δr，低空的Δr取值小；在球面上看，球面元方块大小为Δλ、为球面元覆盖的经度、纬度，离发射点近的区域，Δλ、取值小为规整起见，令每一单元显然，Ωe可由8个节点的坐标表示，为了方便表示及计算，引入局部坐标系表示单元节点；
单元Ωe的局部曲线坐标系由半径rP＝r1+Δr/2的球面、经度λP＝λ1+Δλ/2的子午面、纬度的纬圈的交线组成；原点P为三交线的交点，局部坐标ξ、η、ζ分别沿点P的径向、纬圈与子午线方向；
原点P的球坐标为局部坐标为P(0,0,0)，单元内变点的局部坐标A(ξ,η,ζ)为：
单元顶点的局部坐标Ai(ξi,ηi,ζi)为：
记扰动引力δ在天、东、北方向上的分量为简记为因此δα是球坐标的函数，当将转换成局部坐标(ξ,η,ζ)后，δα可表示成(ξ,η,ζ)的函数；
在单元Ωe内，由于只取了8个节点的信息，因此逼近函数最多可以有8个待定系数；取逼近函数为如下的8次多项式：
δα＝a1+a2ξ+a3η+a4ζ+a5ξη+a6ξζ+a7ηζ+a8ξηζ (3)
节点的δαi由球谐函数法计算得到，记
则待定系数以向量形式表示为：
a＝H-1Z (5)
引入型函数Ni，满足
其中，
则变点A扰动引力分量计算式为：
经推导，Ni(A,Ai)可用下面的公式计算
广义延拓法：
考虑如下三维数值逼近问题；已知函数u(x,y,z):R3→R，在空域上的一组离散数据
{ui|ui＝u(xi,yi,zi),(xi,yi,zi)∈Ω,i＝1,2,…,n}
在空域Ω上构造u的一个近似函数U:Ω→R满足U(xi)＝ui(i＝1,2,…,n)；
采用分块逼近的方法求解，首先将区域进行划分，得m个互不重叠的子区域：
Ω＝Ω1∪Ω2∪…∪Ωm
设第e单元Ωe中包含r个插值节点；将单元Ωe与临近几个子域结合起来形成延拓域Ωe′，其中含s个节点，且s＞r，于是有：
在单元Ωe上，取三元多项式类
{gj(x,y,z)}＝{1,x,y,z,x2,y2,z2,xy,xz,yz,x3,y3,z3,
x2y,x2z,xy2,xz2,yz2,y2z,xyz,…}
的前t项为插值基函数，且r＜t＜s，即令：
其中，a1，a2，…，at为待定系数，可由下述问题解出：
显然，若r＝t＝s，即为通常意义下的分片插值模型，该模型为上节有限元模型；若r＝0,t＜s即为通常意义下的分片拟合模型；故广义延拓逼近是集插值、拟合之长的高精度分片逼近方法。
将空间域作立方体剖分，每个立方体子域Ωe上的节点编号依次为0、1、2、3、4、5、6和7，其延拓域Ωe′上有s个节点，s＞8；对于给定的插值基函数，可以给出待定系数求解的具体表达式；取r＝8，s＝32，t＝20；
将式(10)用矩阵形式表示，令
G＝{gj(xi,yi,zi)}ij,i＝9,10,...,32；j＝1,2,...,t
GI＝{gj(xi,yi,zi)}ij,i＝1,2,...,8；j＝1,2,...,t
u＝[u9 u10 … u32]T
uI＝[u1 u2 … u8]T
a＝[a1 a2 … a20]T
则式可表示为：
引入拉格朗日乘子λ＝[λ1 λ2 … λ8]T，则：
L(a,λ)＝(Ga-u)T(Ga-u)+2(G1a-u1)λ (12)
根据优化原理，待定系数a可通过代数方程
求解，其中各矩阵由如下公式表示，
F0＝GTu,F1＝u1
A＝GTG,C＝G1
A为方阵且可逆，结合分块矩阵求逆公式，最后有：
其中，D11和D12分别是A和C的逆；
令P(x,y,z)＝{gj(x,y,z)|j＝1,2,…,t}，则扰动引力分量
其中，型函数为：
N(x,y,z)＝P(x,y,z)[D11GT D12] (16)
球谐函数法：
以分量形式表示的扰动引力：
其中，μ为地球引力常数，ae为正常地球赤道半径，r为计算点至地心的距离，为计算点地心纬度，λ为计算点经度，为规格化的勒让德伴随函数，N为球谐函数最高阶，为位系数，n为扰动引力的截断阶，m≠n时表示位系数为田谐系数，m＝n时表示位系数为扇谐系数，d为微分号。
本发明目的在于提供了一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法，为用户提供选择引力赋值模型的功能，用户根据需要选择或添加不同的引力模型和影响特性分析模型进行弹道计算仿真，可以实现对不同形态弹道的数值仿真，采用模块化思想构建计算速度快精度高，能够适应不同发射任务、不同引力影响条件，提高了影响特性分析的应用范围和通用性。本发明还提供了一个弹道数据库，便于多种弹道的比较分析及系统功能的扩充：系统将目前成熟的引力计算方法集成在一个计算方法库中，采用模块化构建、面向对象设计，方便不同专业水平的用户使用，另外仿真平台还拥有一个覆盖了航空航天领域绝大多数类型导弹的动力学模型库，如巡航导弹、防空导弹、弹道导弹、运载火箭等。底层算法采用VC++语言编写，具有人机交互界面，计算速度快，易于不同平台间的移植，便于各种用户的使用。
本发明提出的扰动引力对弹道影响特性分析子模块：提供三种分析方法供用户选择，并提示每种方法的优缺点和适用性，用户可根据实际需求自由、便捷地进行选择相应的分析方法或算法组合模式，既可以单独使用又可以组合使用，并且状态空间摄动法获得的落点偏差信息会附加求差法的结果作为比对基准，提高了方法的可信度；利用传播模型分析法虽然不会提供落点偏差信息，但是可以直观反映出导弹飞行过程中各分量扰动引力对弹道的影响程度，一般与前两种分析方法结合使用。
以上对本发明所提供的一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法，进行了详细介绍，本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制。

Claims (9)

1.一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统，包括仿真平台主界面模块，其特征在于：
所述的仿真平台主界面模块包括多种扰动引力赋值算法选择子模块、弹道形态选择子模块、扰动引力对弹道影响特性分析子模块、弹道仿真子模块和数据处理与结果显示子模块；
所述的多种扰动引力赋值算法选择子模块确定引力计算底层的数学形式；
所述的弹道形态选择子模块提供多种弹道模式，通过外部输入程序角完成指定弹道的生成；
所述的扰动引力对弹道影响特性分析子模块是引力影响特性分析的基础，提供三种不同的引力影响特性分析方法，用户根据实际需求选择分析方法，该三种分析方法可单独使用或任意组合使用；
所述的弹道仿真子模块根据底层程序代码及外部输入的程序角用于弹道仿真，得出扰动引力作用下的弹道计算结果；
所述的数据处理与结果显示子模块对得出的弹道计算结果进行数据读取、处理、查询、比较以及分析。
2.根据权利要求1所述的系统，其特征在于：
所述的仿真平台主界面模块是面向用户的，仿真平台主界面模块有专门的帮助子模块的选项菜单和按钮，方便用户随时切换。
3.根据权利要求2所述的系统，其特征在于：
所述的多种扰动引力赋值算法选择子模块从有限元逼近法、广义延拓法及球谐函数法中任意选择一种方法，三种方法是等价的，均可以提供对应模型在所需坐标点处扰动引力三分量，用户依据不同的任务要求选择合适的扰动引力计算方法。
4.根据权利要求3所述的系统，其特征在于：
所述的弹道形态选择子模块提供多种弹道模式，通过外部输入程序角完成指定弹道的生成，如果没有设计好的程序角，则直接使用程序内置的多种典型弹道的程序角数据。
5.根据权利要求4所述的系统，其特征在于：
所述的扰动引力对弹道影响特性分析子模块对积分求差法、摄动法和误差传播分析方法的优缺点和适应性给出提示，用户根据提示及实际需求选择相应的分析方法。
6.根据权利要求5所述的系统，其特征在于：
所述的弹道仿真子模块用于弹道积分计算，得出扰动引力作用下的弹道计算数据，用户根据需要选取坐标系，导弹的总体参数由用户更新相关设置文件，其中相关设置文件包括典型导弹动力学模型、仿真步长、关机条件、落点偏差显示类型，其中典型导弹动力学模型包括气动参数、发动机参数以及不同坐标系间的转换模型。
7.根据权利要求6所述的系统，其特征在于：
所述的数据处理与结果显示子模块对得出的弹道计算结果进行数据读取、处理、查询、比较以及分析，并将结果存入数据库与已经存在的仿真数据进行对照，其中，数据的读取遵循数据交互管理模块的接口定义，通过数据传输机制获得数据信息；数据查询功能依照数据管理系统进行操作，结果分析中还包含了数据处理、比较以及分析模块。
8.一种多赋值法的扰动引力场对不同形态弹道影响特性分析方法，其特征在于：多赋值法的扰动引力场对不同形态弹道影响特性分析方法包括如下步骤：
步骤一、在多赋值法的扰动引力场对不同形态弹道影响特性分析系统运行之前，进行Teechart绘制曲线插件的注册；
步骤二、Teechart绘制曲线插件注册成功之后，在VC6.0以上版本环境中打开多赋值法的扰动引力场对不同形态弹道影响特性分析系统软件，编译.exe文件通过后进入仿真平台主界面模块；
步骤三、在仿真平台主界面模块中的多种扰动引力赋值算法选择子模块对引力赋值方法进行选择，并对相应方法进行设置；其中参数设置包括：格网数划分选择，拟合点个数选择，扰动引力截断阶次，引力场位系数模型选择；
步骤四、在弹道形态选择子模块中设置导弹发射方式、发射点坐标、发射方位角以及射程信息，如果不使用默认数据，则可以自行输入指定弹道程序俯仰角；
步骤五、在扰动引力对弹道影响特性分析子模块中对分析方法进行设置，如果选择弹道积分求差分析方法，则设置标准弹道正常引力模型截断阶次；如果选择摄动法，则设置弹道平均角速度和变换步长；如果选择误差传播分析法，则设置权系数表示方式；在弹道仿真子模块中对仿真步长进行设置；
步骤六、利用数据处理与结果显示子模块对弹道过程数据读取、处理、查询、比较以及分析，并且弹道轨迹仿真数据结果、仿真过程中数据实时更新；到此分析方法结束。
9.根据权利要求8所述的方法，其特征在于：步骤三中多种扰动引力赋值算法包括：
有限元逼近法：
有限元逼近法的基本思路如下：(1)将所考虑的空域Ω用研究选定的网格划分为若干个单元Ωe，Ωe可由点质量法计算出扰动引力的量级确定或指定，Ωe的大小要满足扰动引力的计算精度；(2)将单元Ωe内的扰动引力构造成节点的多项式逼近函数，节点处的扰动引力由点质量法计算并存储；(3)为计算Ω中任一点A的扰动引力δ，需先判断计算点A所在的单元Ωe，再由点A位置与该单元节点值，用逼近函数计算δ；
记所考虑的空域Ω为一球壳，大小可由地固球坐标的取值区间表示；考虑到导弹的飞行空域，将Ω取为
其中，下标“0”、“f”分别表示起始点与终值点；
将Ω用所选择的网格分割为具有形状和大小的单元Ωe，单元的分割有以下要求：
(1)两个相邻单元Ωe在公共边界上是相容的，其顶点不能作为相邻单元边上的内点，即保持连续性；
(2)网格尽量规则，以减少计算量，减轻管理上的难度；
(3)靠近发射点的低空区域，对弹道产生影响的扰动引力集中、变化大，网格要密些，对应的单元Ωe体积小；反之，对远离发射点的高空区域，网格可疏些，Ωe体积大；
根据以上规则，取单元Ωe的形状为由地固球坐标截得的六面体，各面分别为：半径为r1与r2的球面，经度为λ1与λ2的子午面，纬度为与的纬圈，设r1＜r2，λ1＜λ2，沿径向看，Ωe的厚度为Δr，低空的Δr取值小；在球面上看，球面元方块大小为Δλ、为球面元覆盖的经度、纬度，离发射点近的区域，Δλ、取值小，为规整起见，令每一单元显然，Ωe可由8个节点的坐标表示，为了方便表示及计算，引入局部坐标系表示单元节点；
单元Ωe的局部曲线坐标系由半径rP＝r1+Δr/2的球面、经度λP＝λ1+Δλ/2的子午面、纬度的纬圈的交线组成；原点P为三交线的交点，局部坐标ξ、η、ζ分别沿点P的径向、纬圈与子午线方向；
原点P的球坐标为局部坐标为P(0,0,0)，单元内变点的局部坐标A(ξ,η,ζ)为：
单元顶点的局部坐标Ai(ξi,ηi,ζi)为：
记扰动引力δ在天、东、北方向上的分量为简记为因此δα是球坐标的函数，当将转换成局部坐标(ξ,η,ζ)后，δα可表示成(ξ,η,ζ)的函数；
在单元Ωe内，由于只取了8个节点的信息，因此逼近函数最多可以有8个待定系数；取逼近函数为如下的8次多项式：
δα＝a1+a2ξ+a3η+a4ζ+a5ξη+a6ξζ+a7ηζ+a8ξηζ (3)
节点的δαi由球谐函数法计算得到，记
<mrow> <mi>Z</mi> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>&amp;delta;</mi> <mrow> <mi>&amp;alpha;</mi> <mn>1</mn> </mrow> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>&amp;delta;</mi> <mrow> <mi>&amp;alpha;</mi> <mn>2</mn> </mrow> </msub> </mtd> </mtr> <mtr> <mtd> <mo>&amp;CenterDot;</mo> </mtd> </mtr> <mtr> <mtd> <mo>&amp;CenterDot;</mo> </mtd> </mtr> <mtr> <mtd> <mo>&amp;CenterDot;</mo> </mtd> </mtr> <mtr> <mtd> <msub> <mi>&amp;delta;</mi> <mrow> <mi>&amp;alpha;</mi> <mn>8</mn> </mrow> </msub> </mtd> </mtr> </mtable> </mfenced> <mo>,</mo> <mi>a</mi> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>a</mi> <mn>1</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>a</mi> <mn>2</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mo>&amp;CenterDot;</mo> </mtd> </mtr> <mtr> <mtd> <mo>&amp;CenterDot;</mo> </mtd> </mtr> <mtr> <mtd> <mo>&amp;CenterDot;</mo> </mtd> </mtr> <mtr> <mtd> <msub> <mi>a</mi> <mn>8</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mo>,</mo> <mi>H</mi> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <msub> <mi>&amp;xi;</mi> <mn>1</mn> </msub> </mtd> <mtd> <msub> <mi>&amp;eta;</mi> <mn>1</mn> </msub> </mtd> <mtd> <msub> <mi>&amp;zeta;</mi> <mn>1</mn> </msub> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>1</mn> </msub> <msub> <mi>&amp;eta;</mi> <mn>1</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>1</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>1</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;eta;</mi> <mn>1</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>1</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>1</mn> </msub> <msub> <mi>&amp;eta;</mi> <mn>1</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>1</mn> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <msub> <mi>&amp;xi;</mi> <mn>2</mn> </msub> </mtd> <mtd> <msub> <mi>&amp;eta;</mi> <mn>2</mn> </msub> </mtd> <mtd> <msub> <mi>&amp;zeta;</mi> <mn>2</mn> </msub> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>2</mn> </msub> <msub> <mi>&amp;eta;</mi> <mn>2</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>2</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>2</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;eta;</mi> <mn>2</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>2</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>2</mn> </msub> <msub> <mi>&amp;eta;</mi> <mn>2</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>2</mn> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> <mtd> <mo>...</mo> </mtd> </mtr> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <msub> <mi>&amp;xi;</mi> <mn>8</mn> </msub> </mtd> <mtd> <msub> <mi>&amp;eta;</mi> <mn>8</mn> </msub> </mtd> <mtd> <msub> <mi>&amp;zeta;</mi> <mn>8</mn> </msub> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>8</mn> </msub> <msub> <mi>&amp;eta;</mi> <mn>8</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>8</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>8</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;eta;</mi> <mn>8</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>8</mn> </msub> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;xi;</mi> <mn>8</mn> </msub> <msub> <mi>&amp;eta;</mi> <mn>8</mn> </msub> <msub> <mi>&amp;zeta;</mi> <mn>8</mn> </msub> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </mrow>
则待定系数以向量形式表示为：
a＝H-1Z (5)
引入型函数Ni，满足
<mrow> <mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mrow> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>8</mn> </munderover> <msub> <mi>N</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>A</mi> <mo>,</mo> <msub> <mi>A</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>=</mo> <mn>1</mn> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mi>N</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <msub> <mi>A</mi> <mi>j</mi> </msub> <mo>,</mo> <msub> <mi>A</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>=</mo> <msub> <mi>&amp;delta;</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>,</mo> <mi>i</mi> <mo>,</mo> <mi>j</mi> <mo>=</mo> <mn>1</mn> <mo>,</mo> <mn>2</mn> <mo>,</mo> <mo>...</mo> <mo>,</mo> <mn>8</mn> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>6</mn> <mo>)</mo> </mrow> </mrow>
其中，
<mrow> <msub> <mi>&amp;delta;</mi> <mrow> <mi>i</mi> <mi>j</mi> </mrow> </msub> <mo>=</mo> <mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <mrow> <mi>i</mi> <mo>&amp;NotEqual;</mo> <mi>j</mi> </mrow> </mtd> </mtr> <mtr> <mtd> <mn>1</mn> </mtd> <mtd> <mrow> <mi>i</mi> <mo>=</mo> <mi>j</mi> </mrow> </mtd> </mtr> </mtable> </mfenced> </mrow>
则变点A扰动引力分量计算式为：
<mrow> <msub> <mi>&amp;delta;</mi> <mi>&amp;alpha;</mi> </msub> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>8</mn> </munderover> <msub> <mi>N</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>A</mi> <mo>,</mo> <msub> <mi>A</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>&amp;CenterDot;</mo> <msub> <mi>&amp;delta;</mi> <mrow> <mi>&amp;alpha;</mi> <mi>i</mi> </mrow> </msub> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow>
经推导，Ni(A,Ai)可用下面的公式计算
<mrow> <msub> <mi>N</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>A</mi> <mo>,</mo> <msub> <mi>A</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mn>1</mn> <mn>8</mn> </mfrac> <mrow> <mo>(</mo> <mn>1</mn> <mo>+</mo> <mfrac> <mi>&amp;xi;</mi> <msub> <mi>&amp;xi;</mi> <mi>i</mi> </msub> </mfrac> <mo>)</mo> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>+</mo> <mfrac> <mi>&amp;eta;</mi> <msub> <mi>&amp;eta;</mi> <mi>i</mi> </msub> </mfrac> <mo>)</mo> </mrow> <mrow> <mo>(</mo> <mn>1</mn> <mo>+</mo> <mfrac> <mi>&amp;zeta;</mi> <msub> <mi>&amp;zeta;</mi> <mi>i</mi> </msub> </mfrac> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>8</mn> <mo>)</mo> </mrow> </mrow>
广义延拓法：
考虑如下三维数值逼近问题；已知函数u(x,y,z):R3→R，在空域上的一组离散数据
{ui|ui＝u(xi,yi,zi),(xi,yi,zi)∈Ω,i＝1,2,…,n}
在空域Ω上构造u的一个近似函数U:Ω→R满足U(xi)＝ui(i＝1,2,…,n)；
采用分块逼近的方法求解，首先将区域进行划分，得m个互不重叠的子区域：
Ω＝Ω1∪Ω2∪…∪Ωm
设第e单元Ωe中包含r个插值节点；将单元Ωe与临近几个子域结合起来形成延拓域Ω′e，其中含s个节点，且s＞r，于是有：
<mrow> <msubsup> <mi>Q</mi> <mi>e</mi> <mo>&amp;prime;</mo> </msubsup> <mo>&amp;Superset;</mo> <msub> <mi>&amp;Omega;</mi> <mi>e</mi> </msub> <mo>,</mo> <mi>e</mi> <mo>=</mo> <mn>1</mn> <mo>,</mo> <mn>2</mn> <mo>,</mo> <mo>...</mo> <mo>,</mo> <mi>m</mi> </mrow>
在单元Ωe上，取三元多项式类
{gj(x,y,z)}＝{1,x,y,z,x2,y2,z2,xy,xz,yz,x3,y3,z3,
x2y,x2z,xy2,xz2,yz2,y2z,xyz,…}
的前t项为插值基函数，且r＜t＜s，即令：
<mrow> <mtable> <mtr> <mtd> <mrow> <mi>U</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>t</mi> </munderover> <msub> <mi>a</mi> <mi>j</mi> </msub> <msub> <mi>g</mi> <mi>j</mi> </msub> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> <mo>&amp;Element;</mo> <msub> <mi>&amp;Omega;</mi> <mi>e</mi> </msub> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>9</mn> <mo>)</mo> </mrow> </mrow>
其中，a1，a2，…，at为待定系数，可由下述问题解出：
<mrow> <mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mrow> <mi>min</mi> <mi> </mi> <mi>I</mi> <mrow> <mo>(</mo> <msub> <mi>a</mi> <mn>1</mn> </msub> <mo>,</mo> <msub> <mi>a</mi> <mn>2</mn> </msub> <mo>,</mo> <mo>...</mo> <mo>,</mo> <msub> <mi>a</mi> <mi>t</mi> </msub> <mo>)</mo> </mrow> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mi>r</mi> <mo>+</mo> <mn>1</mn> </mrow> <mi>s</mi> </munderover> <msup> <mrow> <mo>(</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>t</mi> </munderover> <msub> <mi>a</mi> <mi>j</mi> </msub> <msub> <mi>g</mi> <mi>j</mi> </msub> <mo>(</mo> <msub> <mi>x</mi> <mi>i</mi> </msub> <mo>,</mo> <msub> <mi>y</mi> <mi>i</mi> </msub> <mo>,</mo> <msub> <mi>z</mi> <mi>i</mi> </msub> <mo>)</mo> <mo>-</mo> <msub> <mi>u</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mi>s</mi> <mo>.</mo> <mi>t</mi> <mo>.</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>t</mi> </munderover> <msub> <mi>a</mi> <mi>j</mi> </msub> <msub> <mi>g</mi> <mi>j</mi> </msub> <mrow> <mo>(</mo> <msub> <mi>x</mi> <mi>i</mi> </msub> <mo>,</mo> <msub> <mi>y</mi> <mi>i</mi> </msub> <mo>,</mo> <msub> <mi>z</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>=</mo> <msub> <mi>u</mi> <mi>i</mi> </msub> <mo>,</mo> <mi>i</mi> <mo>=</mo> <mn>1</mn> <mo>,</mo> <mn>2</mn> <mo>,</mo> <mo>...</mo> <mo>,</mo> <mi>r</mi> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>10</mn> <mo>)</mo> </mrow> </mrow>
显然，若r＝t＝s，即为通常意义下的分片插值模型，该模型为上节有限元模型；若r＝0,t＜s即为通常意义下的分片拟合模型；故广义延拓逼近是集插值、拟合之长的高精度分片逼近方法。
将空间域作立方体剖分，每个立方体子域Ωe上的节点编号依次为0、1、2、3、4、5、6和7，其延拓域Ω′e上有s个节点，s＞8；对于给定的插值基函数，可以给出待定系数求解的具体表达式；取r＝8，s＝32，t＝20；
将式(10)用矩阵形式表示，令
G＝{gj(xi,yi,zi)}ij,i＝9,10,...,32；j＝1,2,...,t
GI＝{gj(xi,yi,zi)}ij,i＝1,2,...,8；j＝1,2,...,t
u＝[u9 u10 … u32]T
uI＝[u1 u2 … u8]T
a＝[a1 a2 … a20]T
则式可表示为：
<mrow> <mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mrow> <mi>min</mi> <mi> </mi> <mi>I</mi> <mrow> <mo>(</mo> <mi>a</mi> <mo>)</mo> </mrow> <mo>=</mo> <msup> <mrow> <mo>(</mo> <mi>G</mi> <mi>a</mi> <mo>-</mo> <mi>u</mi> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mrow> <mo>(</mo> <mi>G</mi> <mi>a</mi> <mo>-</mo> <mi>u</mi> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mi>s</mi> <mo>.</mo> <mi>t</mi> <mo>.</mo> <msub> <mi>G</mi> <mn>1</mn> </msub> <mi>a</mi> <mo>-</mo> <msub> <mi>u</mi> <mn>1</mn> </msub> <mo>=</mo> <mn>0</mn> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>11</mn> <mo>)</mo> </mrow> </mrow>
引入拉格朗日乘子λ＝[λ1 λ2 … λ8]T，则：
L(a,λ)＝(Ga-u)T(Ga-u)+2(G1a-u1)λ (12)
根据优化原理，待定系数a可通过代数方程
<mrow> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mi>A</mi> </mtd> <mtd> <msup> <mi>C</mi> <mi>T</mi> </msup> </mtd> </mtr> <mtr> <mtd> <mi>C</mi> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> </mtable> </mfenced> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mi>a</mi> </mtd> </mtr> <mtr> <mtd> <mi>&amp;lambda;</mi> </mtd> </mtr> </mtable> </mfenced> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>F</mi> <mn>0</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>F</mi> <mn>1</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>13</mn> <mo>)</mo> </mrow> </mrow>
求解，其中各矩阵由如下公式表示，
F0＝GTu,F1＝u1
A＝GTG,C＝G1
A为方阵且可逆，结合分块矩阵求逆公式，最后有：
<mrow> <mi>a</mi> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>D</mi> <mn>11</mn> </msub> </mtd> <mtd> <msub> <mi>D</mi> <mn>12</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>F</mi> <mn>0</mn> </msub> </mtd> </mtr> <mtr> <mtd> <msub> <mi>F</mi> <mn>1</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mrow> <msub> <mi>D</mi> <mn>11</mn> </msub> <msup> <mi>G</mi> <mi>T</mi> </msup> </mrow> </mtd> <mtd> <msub> <mi>D</mi> <mn>12</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mi>u</mi> </mtd> </mtr> <mtr> <mtd> <msub> <mi>u</mi> <mn>1</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>14</mn> <mo>)</mo> </mrow> </mrow>
其中，D11和D12分别是A和C的逆；
令P(x,y,z)＝{gj(x,y,z)|j＝1,2,…,t}，则扰动引力分量
<mrow> <msub> <mi>&amp;delta;</mi> <mi>&amp;alpha;</mi> </msub> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>N</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>,</mo> <mi>y</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mi>u</mi> </mtd> </mtr> <mtr> <mtd> <msub> <mi>u</mi> <mn>1</mn> </msub> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>15</mn> <mo>)</mo> </mrow> </mrow>
其中，型函数为：
N(x,y,z)＝P(x,y,z)[D11GT D12] (16)
球谐函数法：
以分量形式表示的扰动引力：
其中，μ为地球引力常数，ae为正常地球赤道半径，r为计算点至地心的距离，为计算点地心纬度，λ为计算点经度，为规格化的勒让德伴随函数，N为球谐函数最高阶，为位系数，n为扰动引力的截断阶，m≠n时表示位系数为田谐系数，m＝n时表示位系数为扇谐系数，d为微分号。
CN201611075769.5A 2016-11-30 2016-11-30 一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法 Active CN106599410B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201611075769.5A CN106599410B (zh) 2016-11-30 2016-11-30 一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201611075769.5A CN106599410B (zh) 2016-11-30 2016-11-30 一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法

Publications (2)

Publication Number Publication Date
CN106599410A CN106599410A (zh) 2017-04-26
CN106599410B true CN106599410B (zh) 2018-02-06

Family

ID=58593884

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201611075769.5A Active CN106599410B (zh) 2016-11-30 2016-11-30 一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法

Country Status (1)

Country Link
CN (1) CN106599410B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107977486B (zh) * 2017-11-06 2019-08-09 北京宇航系统工程研究所 一种地球扰动引力场球冠谐模型阶次扩展方法及系统
CN107832565A (zh) * 2017-12-06 2018-03-23 内蒙动力机械研究所 一种固体发动机一维内弹道建模与性能预示软件系统
CN110046439B (zh) * 2019-04-22 2020-05-19 中国人民解放军国防科技大学 考虑高阶扰动引力影响的弹道偏差解析预报算法
CN111475767B (zh) * 2020-03-18 2021-03-16 中国科学院紫金山天文台 考虑地球自转影响的最小能量弹道严格构造方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103487801A (zh) * 2013-09-02 2014-01-01 电子科技大学 雷达识别真假弹头的方法
CN103913991A (zh) * 2014-04-22 2014-07-09 西北工业大学 高速轴对称飞行器复合控制方法
CN104751012A (zh) * 2015-04-23 2015-07-01 中国人民解放军国防科学技术大学 沿飞行弹道的扰动引力快速逼近方法
CN104750983A (zh) * 2015-03-20 2015-07-01 中国人民解放军信息工程大学 一种空间分层网格扰动引力场模型构建与扰动引力快速确定方法
CN105184109A (zh) * 2015-10-27 2015-12-23 中国人民解放军国防科学技术大学 扰动引力作用下弹道助推段状态偏差解析计算方法
CN105447228A (zh) * 2015-11-11 2016-03-30 中国人民解放军国防科学技术大学 一种地球扰动引力赋值模型适用性能评估方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103487801A (zh) * 2013-09-02 2014-01-01 电子科技大学 雷达识别真假弹头的方法
CN103913991A (zh) * 2014-04-22 2014-07-09 西北工业大学 高速轴对称飞行器复合控制方法
CN104750983A (zh) * 2015-03-20 2015-07-01 中国人民解放军信息工程大学 一种空间分层网格扰动引力场模型构建与扰动引力快速确定方法
CN104751012A (zh) * 2015-04-23 2015-07-01 中国人民解放军国防科学技术大学 沿飞行弹道的扰动引力快速逼近方法
CN105184109A (zh) * 2015-10-27 2015-12-23 中国人民解放军国防科学技术大学 扰动引力作用下弹道助推段状态偏差解析计算方法
CN105447228A (zh) * 2015-11-11 2016-03-30 中国人民解放军国防科学技术大学 一种地球扰动引力赋值模型适用性能评估方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
弹道导弹上升段扰动引力影响分析及补偿方法研究;朱晨昊;《中国优秀硕士学位论文全文数据库》;20150228;全文 *
空间扰动引力场网格模型构建方法研究;范昊鹏;《中国优秀硕士学位论文全文数据库》;20160731;3-5、24、28-32 *

Also Published As

Publication number Publication date
CN106599410A (zh) 2017-04-26

Similar Documents

Publication Publication Date Title
CN103557867B (zh) 一种基于稀疏a*搜索的三维多uav协同航迹规划方法
CN106599410B (zh) 一种多赋值法的扰动引力场对不同形态弹道影响特性分析系统及方法
Langener et al. Trajectory simulation and optimization of the LAPCAT MR2 hypersonic cruiser concept
CN102073755B (zh) 近空间高超声速飞行器运动控制仿真方法
CN103942401A (zh) 一种优化高精度自适应模块化的航天器弹道多约束轨迹工具包及方法
CN104007665A (zh) 一种固液动力飞行器飞行仿真测试系统
Qi et al. An improved heuristic algorithm for UAV path planning in 3D environment
CN103743402B (zh) 一种基于地形信息量的水下智能自适应地形匹配方法
CN103412493A (zh) 固液动力巡航飞行器飞行任务规划仿真系统
CN101256080A (zh) 卫星/惯性组合导航系统的空中对准方法
CN104751012A (zh) 沿飞行弹道的扰动引力快速逼近方法
CN112733409B (zh) 一种多源传感综合集成复合导航微系统协同设计平台
CN114280934B (zh) 一种可重复使用运载火箭全程轨迹规划方法
CN103942369B (zh) 一种面向临近空间的智能目标发生方法
CN103542853B (zh) 一种估计加速度计漂移的绝对导航滤波方法
CN105300383A (zh) 一种基于回溯搜索的无人机空中加油位姿估计方法
CN113268087A (zh) 多约束复杂环境下基于改进蚁群算法的多无人机协同工作的航迹规划方法
CN109211230A (zh) 一种基于牛顿迭代法的炮弹姿态和加速度计常值误差估计方法
CN105069311A (zh) 一种远程火箭发射初态误差传播估计方法
CN107122540A (zh) 一种基于斯托克斯积分法的地面测绘保障条件需求分析系统与分析方法
Zhang et al. Hybrid hierarchical trajectory planning for a fixed-wing UCAV performing air-to-surface multi-target attack
Zheng et al. Multidisciplinary integrated design of long-range ballistic missile using PSO algorithm
CN115342812A (zh) 一种基于改进蝴蝶优化算法的无人机三维航迹规划方法
CN105760573B (zh) 沿临近空间大范围机动弹道的扰动引力延拓逼近方法
CN106021835A (zh) 一种面向最优侦察的航迹设计方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
TR01 Transfer of patent right

Effective date of registration: 20210114

Address after: Building 9, accelerator, 14955 Zhongyuan Avenue, Songbei District, Harbin City, Heilongjiang Province

Patentee after: INDUSTRIAL TECHNOLOGY Research Institute OF HEILONGJIANG PROVINCE

Address before: 150006 No. 92 West straight street, Nangang District, Heilongjiang, Harbin

Patentee before: HARBIN INSTITUTE OF TECHNOLOGY

TR01 Transfer of patent right
TR01 Transfer of patent right

Effective date of registration: 20230330

Address after: 150027 Room 412, Unit 1, No. 14955, Zhongyuan Avenue, Building 9, Innovation and Entrepreneurship Plaza, Science and Technology Innovation City, Harbin Hi tech Industrial Development Zone, Heilongjiang Province

Patentee after: Heilongjiang Industrial Technology Research Institute Asset Management Co.,Ltd.

Address before: Building 9, accelerator, 14955 Zhongyuan Avenue, Songbei District, Harbin City, Heilongjiang Province

Patentee before: INDUSTRIAL TECHNOLOGY Research Institute OF HEILONGJIANG PROVINCE

TR01 Transfer of patent right