一种地球扰动引力赋值模型适用性能评估方法

摘要

本发明公开了一种地球扰动引力赋值模型适用性能评估方法，步骤包括：设置飞行轨迹样本试验因素及其试验水平，完成飞行轨迹样本的试验设计；设计单项指标；根据所述飞行轨迹样本的试验设计，对待评估地球扰动引力赋值模型进行弹道式飞行器飞行轨迹仿真，并计算各单项指标的数据得到各数据样本库；利用贝叶斯估计方法完成各个单项指标的评估分析；对各个单项指标的评估分析结果进行无量纲化处理，得到各个单项指标的隶属度；综合计算待评估地球扰动引力赋值模型的适用性能评估结果。本发明能够完成任意给定形式地球扰动引力赋值模型适用性能的评估和分析以确定其在弹道式飞行器控制领域的适用性能，具有评估全面合理、快速可靠、实施便捷的优点。

说明书

技术领域

本发明涉及飞行器控制技术，具体涉及一种地球扰动引力赋值模型适用性能评估 方法。

背景技术

引起远程弹道式飞行器落点偏差的因素主要包括制导工具误差和制导方法误差。 制导工具误差是由于惯性平台、加速度表、陀螺仪和弹载计算机等制导设备性能的不完善 引起的落点偏差；制导方法误差主要是制导方法本身无法在根本上消除扰动引力误差、发 动机特性误差、飞行器结构误差和气动系数误差等干扰因素对制导精度的影响从而造成的 落点偏差，且随着制导硬件设备的不断改进，制导方法误差的影响日益突出。根据现有经验 和相关资料，发动机特性误差、飞行器结构误差和气动系数误差等可以作用于飞行器视加 速度的干扰因素均能通过惯性导航系统和显式制导方法进行消除。相比较而言，作为制导 方法误差主要误差源的扰动引力场误差并不能通过惯性导航系统敏感，只能通过扰动引力 场建模进而在制导模型中进行修正；因此，为有效消除扰动引力场对远程弹道式飞行器制 导精度的影响，必须进行扰动引力场精确模型建模、扰动引力赋值、扰动引力场影响制导补 偿等层面的研究工作。需要说明的是，虽然高精度扰动引力场建模研究已逐步形成模型逼 近与数值逼近两种主要方法；扰动引力赋值和扰动引力模型重构快速赋值也基于不同数值 计算思想完成了扰动引力不同精度的赋值要求；然而将上述两项研究的相关成果应用于扰 动引力场影响制导补偿理论方法和工程应用之前，均要求对地球扰动引力场赋值模型的适 用性能进行系统和全面的评估，进而明确该地球扰动引力赋值模型的适用程度，因此开展 地球扰动引力场赋值模型适用性能评估方法的研究将具有重要的意义。

发明内容

本发明要解决的技术问题：针对现有技术的上述问题，提供一种能够完成任意给 定形式地球扰动引力赋值模型适用性能的评估和分析以确定其在弹道式飞行器控制领域 的适用性能，评估全面合理、快速可靠、实施便捷的地球扰动引力赋值模型适用性能评估方 法。

为了解决上述技术问题，本发明采用的技术方案为：

一种地球扰动引力赋值模型适用性能评估方法，步骤包括：

1)设置飞行轨迹样本试验因素及其试验水平，利用正交设计试验方法完成用于地 球扰动引力赋值模型的适用性能评估的飞行轨迹样本的试验设计；

2)从赋值精度、赋值速度和赋值模型适应性三方面设计用于适用性能评估的单项 指标以建立地球扰动引力赋值模型适用性能评估的指标体系；

3)根据所述飞行轨迹样本的试验设计，针对待评估地球扰动引力赋值模型进行弹 道式飞行器飞行轨迹仿真，并计算各个单项指标的数据，得到各个单项指标的数据样本库；

4)假设各个单项指标的样本数据变化规律满足正态分布统计规律，基于所述假设 利用贝叶斯估计方法完成各个单项指标的评估分析；

5)采用预设的无量纲隶属度函数，对各个单项指标的评估分析结果进行无量纲化 处理，得到各个单项指标的隶属度；

6)根据各个单项指标的预设权重及各个单项指标的隶属度，综合计算待评估地球 扰动引力赋值模型的适用性能评估结果。

优选地，所述步骤1)中设置飞行轨迹样本试验因素及其试验水平时，设置的飞行 轨迹样本试验因素包括：弹道式飞行器的发射位置地形、轨迹控制方式、航程、推进装置类 型、飞行轨迹类型以及初始时刻飞行方向。

优选地，所述步骤1)中设置飞行轨迹样本试验因素及其试验水平时，弹道式飞行 器的发射位置地形的试验水平为4水平，轨迹控制方式的试验水平为3水平，航程的试验水 平为6水平，推进装置类型的试验水平为3水平，飞行轨迹类型的试验水平为4水平，初始时 刻飞行方向的试验水平为10水平。

优选地，所述步骤2)中选取用于适用性能评估的单项指标时，根据赋值精度选取 的单项指标包括：扰动引力赋值残差δ_t、赋值模型落点偏差ε_x和ε_z、赋值残差百分比根据 赋值速度选取的单项指标包括：扰动引力赋值时间T_c、射前计算时间T_SQ、机载计算机解算时 间T_DS；根据赋值模型适应性选取的单项指标包括：不同飞行轨迹模式适应性不同轨 迹控制方式适应性不同地形适应性

优选地，所述扰动引力赋值残差δ_t的计算式如式(1)所示，所述赋值模型落点偏差 ε_x和ε_z的计算式如式(2)所示，所述赋值残差百分比的计算式如式(3)所示；所述扰动引 力赋值时间T_c的计算式如式(4)所示，所述射前计算时间T_SQ的计算式如式(5)所示，所述机 载计算机解算时间T_DS的计算式如式(6)所示；所述不同飞行轨迹模式适应性的计算式 如式(7)所示，所述不同轨迹控制方式适应性的计算式如式(8)所示，所述不同地形适 应性的计算式如式(9)所示；

${\delta }_t=\left|\right|g_r^t-g_c^t\left|\right|---\left(1\right)$

式(1)中，δ_t表示任意t时刻弹道式飞行器所在空间位置的扰动引力赋值残差，表示任意t时刻弹道式飞行器所在空间位置的实际扰动引力值，表示待评估地球扰动引 力赋值模型计算所得扰动引力值，表示实际扰动引力值与待评估地球扰动引力 赋值模型计算所得扰动引力值的差值的欧几里得范数；

$\left(\begin{array}{c}{ϵ}_x=|x_c-x_n|\\ {ϵ}_z=|z_c-z_n|\end{array}\right)---\left(2\right)$

式(2)中，以瞄准点为原点O、纵向和横向分别作为OX轴和OZ轴、落点偏离用(X,Z) 表示，ε_x表示X方向的落点偏差；ε_z表示Z方向的落点偏差；x_c表示利用待评估地球扰动引力 赋值模型进行飞行轨迹仿真所得落点的X向分量，x_n表示理想目标点位置坐标的X向分量，z_n表示理想目标点位置坐标的Z向分量，z_c表示利用待评估地球扰动引力赋值模型进行飞行 轨迹仿真所得落点的Z向分量；

${\zeta }_{\delta }^t={\delta }_t/\left|\right|g_r^t\left|\right|---\left(3\right)$

式(3)中，表示任意t时刻弹道式飞行器所在空间位置的赋值残差百分比，δ_t表 示任意t时刻弹道式飞行器所在空间位置的扰动引力赋值残差，表示任意t时刻弹道式 飞行器所在空间位置的实际扰动引力值的欧几里得范数；

$T_c=\underset{i=1}{\overset{n}{\Sigma }}{t}_c^i---\left(4\right)$

式(4)中，T_c表示扰动引力赋值时间，n表示调用待评估地球扰动引力赋值模型进 行扰动引力赋值计算的次数，表示第i次调用待评估地球扰动引力赋值模型进行扰动引 力赋值计算所用的时间；

T_SQ＝T_CX+T_ZY(5)

式(5)中，T_SQ表示射前计算时间，T_CX表示弹道式飞行器初始位置垂线偏差获取时 间，T_ZY表示弹道式飞行器诸元计算时间；

$T_{DS}=\underset{i=1}{\overset{n}{\Sigma }}{t}_{DS}^i---\left(6\right)$

式(4)中，T_DS表示机载计算机解算时间，n表示弹道式飞行器飞行轨迹仿真过程中 待评估地球扰动引力赋值模型被调用的次数，表示弹道式飞行器模拟机载计算机第i次 调用待评估地球扰动引力赋值模型进行扰动引力解算所需的时间；

${ϵ}_{traj}^{adapt}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}\frac{e_{traj,i}^{fit}}{C_{traj,i}^{all}}---\left(7\right)$

式(7)中，表示不同飞行轨迹模式适应性，表示针对第i种飞行轨迹模式 进行的仿真飞行轨迹总数目，表示总飞行轨迹仿真数目中满足终端要求的飞行轨迹数 目，n表示飞行轨迹模式数目；

${ϵ}_{ctl}^{adapt}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}\frac{C_{ctl,i}^{fit}}{C_{ctl,i}^{all}}---\left(8\right)$

式(8)中，表示不同轨迹控制方式适应性，表示针对第i种轨迹控制方式 进行的仿真飞行轨迹总数目，表示总飞行轨迹仿真数目中满足轨迹控制要求的飞行轨 迹数目，n表示轨迹控制方式的数目；

${ϵ}_{hypso}^{adapt}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}\frac{C_{hypso,i}^{fit}}{C_{hypso,i}^{all}}---\left(9\right)$

式(9)中，表示不同地形适应性，表示针对第i种地形进行的仿真飞行 轨迹总数目，表示总飞行轨迹仿真数目中满足飞行任务要求的飞行轨迹数目，n表示 地形类型的数目。

优选地，所述步骤3)的详细步骤包括：

3.1)根据预设的空间边界约束或飞行器的飞行区域约束确定待评估地球扰动引 力赋值模型应用的空间范围约束；

3.2)针对空间范围约束的空间边界约束，结合待评估地球扰动引力赋值模型规定 的空间剖分准则完成待评估地球扰动引力赋值模型指定空域的空域剖分；

3.3)利用待评估地球扰动引力赋值模型根据空域剖分形成的空间有限元节点进 行扰动引力赋值，得到待评估地球扰动引力赋值模型的有限元节点数值模型；

3.4)在所述有限元节点数值模型的基础上，完成空间中任意位置扰动引力赋值的 解算和用于弹道式飞行器的飞行轨迹仿真；

3.5)在弹道式飞行器飞行轨迹仿真的过程中计算各个单项指标的数据，得到各个 单项指标的数据样本库。

优选地，所述步骤4)中假设的正态分布统计规律为其中，I_s.j表 示第j个单项指标，μ_Is.j表示待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性 误差决定的第j个单项指标的待估参数数值，表示数据样本中第j个单项指标值 相对单项指标I_s.j的散步程度，D_Is.j表示第j个单项指标的方差，σ_Is.j表示第j个单项指标的 标准差，第j个单项指标的方差D_Is.j、标准差σ_Is.j均为待估参数；所述利用贝叶斯估计方法完 成各项单项指标的评估分析的步骤包括：

4.1)从所有单项指标中取出一项作为当前单项指标j；

4.2)在当前单项指标j的数据样本库中任意选取n_QF组仿真样本数据作 为验前子样，则当前单项指标j的待估参数μ_Is.j的验前概率密度为其中 表示仿真样本数据均值，仿真样本数据均值的计算式如式(10)所示；D_Is.j表示当前 单项指标j的方差，当前单项指标j的方差D_Is.j的验前概率密度为逆Gamma分布的密度函数g (D_x；α₀,β₀)，其中间参数α₀和β₀的计算式如式(11)所示，n_QF表示当前单项指标j的验前子样 的样本数目；

${\bar{x}}^{QF}=\frac{1}{n_{QF}}\underset{i=1}{\overset{n_{QF}}{\Sigma }}{x}_i^{QF}---\left(10\right)$

式(10)中，表示仿真样本数据均值，n_QF表示当前单项指标j的验前子样的样本 数目，表示当前单项指标j的验前子样中第i个样本的样本值；

$\left(\begin{array}{c}{\alpha }_0=\frac{1}{2}\underset{i=1}{\overset{n_{QF}}{\Sigma }}{(x_i^{QF}-{\bar{x}}^{QF})}^2\\ {\beta }_0=\frac{n_{QF}-1}{2}\end{array}\right)---\left(11\right)$

式(11)中，α₀和β₀分别表示逆Gamma分布的密度函数g(D_x；α₀,β₀)的中间参数，n_QF表 示当前单项指标j的验前子样的样本数目，表示当前单项指标j的验前子样中第i个样 本的样本值，表示当前单项指标j的仿真样本数据均值；

4.3)假定当前单项指标j的数据样本库中共有n个样本x₁,…,x_n，确定(μ_Is.j,D_Is.j) 的验后和验前是共轭的且服从正态—逆Gamma分布如式(12)所示；

$\pi ({\mu }_{Is.j},D_{Is.j}|X)\propto \frac{1}{\sqrt{D_{Is.j}}}{e}^{-\frac{1}{2{\eta }_1{D}_{Is.j}}{({\mu }_{Is.j}-{\mu }_1)}^2}{D_{Is.j}}^{-({\beta }_1+1)}{e}^{-\frac{{\alpha }_1}{D_{Is.j}}}---\left(12\right)$

式(12)中，μ_Is.j表示待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性 误差决定的当前单项指标j的待估参数数值，D_Is.j表示当前单项指标j的方差，X表示当前单 项指标j的数据样本库中n个样本x₁,…,x_n构成的数据集，中间参数η₁、α₁、β₁、μ₁如式(13)所 示；

$\left(\begin{array}{c}{\eta }_1=\frac{1}{n+n_{QF}}\\ {\alpha }_1={\alpha }_0+\frac{n}{2}·u+\frac{1}{2}·\frac{n{(\bar{x}-{\bar{x}}^{QF})}^2}{\frac{n}{n_{QF}}+1},u=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}{(x_i-\bar{x})}^2,\bar{x}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}{x}_i\\ {\beta }_1={\beta }_0+\frac{n}{2}\\ {\mu }_1=\frac{n\bar{x}+n_{QF}{\bar{x}}^{QF}}{n+n_{QF}}\end{array}\right)---\left(13\right)$

式(13)中，η₁、α₁、β₁、μ₁分别为中间参数，n表示当前单项指标j的样本数据总数目， n_QF表示当前单项指标j的验前子样的样本数目，α₀,β₀为中间参数，u表示当前单项指标j的 数据样本库中样本数据的方差，x_i表示当前单项指标j的数据样本库中的第i个样本值，表示当前单项指标j的数据样本库中所有样本数据的均值，表示当前单项指标j的仿真 样本数据均值；

4.4)确定待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性误差决定 的当前单项指标j的待估参数数值μ_Is.j如式(14)所示；

${\mu }_{Is.j}={\left(\frac{{\alpha }_1{\eta }_1}{{\beta }_1}\right)}^{\frac{1}{2}}Y+{\mu }_1,p_Y\left(y\right)\propto {(1+\frac{y^2}{2{\beta }_1})}^{-\frac{2{\beta }_1+1}{2}}---\left(14\right)$

式(14)中，μ_Is.j表示待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性 误差决定的当前单项指标j的待估参数数值，η₁、α₁、β₁、μ₁分别为式(13)所示的中间参数，Y 表示学生氏t分布，p_Y(y)表示Y的概率密度函数；

4.5)确定当前单项指标j的待估参数μ_Is.j的贝叶斯估计值及其方差Var(μ_Is.j) 如式(15)所示，当前单项指标j的方差D_Is.j的贝叶斯估计值及其方差Var(D_Is.j)如式 (16)所示，当前单项指标j的标准差σ_Is.j的贝叶斯估计值及其方差Var(σ_Is.j)如式(17) 所示；

$\left(\begin{array}{c}{\hat{\mu }}_{Is.j}={\mu }_1\\ Var\left({\mu }_{Is.j}\right)=\frac{{\alpha }_1{\eta }_1}{{\beta }_1-1}\end{array}\right)---\left(15\right)$

式(15)中，表示当前单项指标j的待估参数μ_Is.j的贝叶斯估计值，μ_Is.j表示待 评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性误差决定的当前单项指标j的待 估参数数值，Var(μ_Is.j)表示贝叶斯估计值的方差，参数η₁、α₁、β₁、μ₁分别表示式(13)所 示的中间参数；

$\left(\begin{array}{c}{\hat{D}}_{Is.j}=\frac{{\alpha }_1}{{\beta }_1-1}\\ Var\left(D_{Is·j}\right)=\frac{{\alpha }_1^2}{{({\beta }_1-1)}^2({\beta }_1-2)}\end{array}\right)---\left(16\right)$

式(16)中，表示当前单项指标j的方差D_Is.j的贝叶斯估计值，D_Is.j表示当前单 项指标j的方差，Var(D_Is.j)表示当前单项指标j的方差D_Is.j的贝叶斯估计值，参数η₁、α₁、β₁、 μ₁分别为式(13)所示的中间参数；

$\left(\begin{array}{c}{\hat{\sigma }}_{Is.j}=\sqrt{{\alpha }_1}\frac{\Gamma ({\beta }_1-1/2)}{\Gamma \left({\beta }_1\right)}\\ Var\left({\sigma }_{Is.j}\right)=\frac{{\alpha }_1}{{\beta }_1-1}-{\alpha }_1{\left(\frac{\Gamma ({\beta }_1-1/2)}{\Gamma \left({\beta }_1\right)}\right)}^2\end{array}\right)---\left(17\right)$

式(17)中，表示当前单项指标j的标准差σ_Is.j的贝叶斯估计值，σ_Is.j表示当前 单项指标j的标准差，Var(σ_Is.j)表示贝叶斯估计值的方差，参数η₁、α₁、β₁、μ₁分别为式 (13)所示的中间参数，Γ为伽玛函数；

4.6)根据当前单项指标j的仿真样本数据均值待估参数μ_Is.j的贝叶斯估计值 及其方差Var(μ_Is.j)、方差D_Is.j的贝叶斯估计值及其方差Var(D_Is.j)、标准差σ_Is.j的 贝叶斯估计值及其Var(σ_Is.j)，对当前单项指标j进行定量或者定性评估；

4.7)判断是否所有单项指标已经处理完毕，如果尚未处理完毕，则跳转执行步骤 4.1)，否则判断利用贝叶斯估计方法完成各项单项指标的评估分析完毕，跳转执行步骤5)。

优选地，所述步骤5)对各个单项指标的评估分析结果进行无量纲化处理时：

◆针对根据赋值精度选取的单项指标采用的无量纲隶属度函数如式(18)所示；

$\bar{R}=e^{-k(R-R_0)},R>R_0---\left(18\right)$

式(18)中，表示赋值精度的各单项指标的无量纲化隶属度函数值，k和R₀为与各 单项指标的数值和评估要求相关的未知参数，R表示赋值精度的各单项指标自变量；针对每 一个单项指标，设计两个已知点位来确定未知参数k和R₀，从而得到该单项指标的无量纲隶 属度函数；针对扰动引力赋值残差δ_t，设置当扰动引力赋值残差δ_t等于扰动引力赋值残差的 最大允许值δ_r.max时，其无量纲隶属度函数值S_δ等于赋值残差无量纲隶属度函数的预设最小 值S_δ.min，其中S_δ.min∈[0,1]；当扰动引力赋值残差δ_t等于扰动引力赋值残差的最小预设值 δ_r.min时，无量纲隶属度函数值S_δ等于赋值残差无量纲隶属度函数的预设最大值S_δ.max，其中 S_δ.max∈[0,1]；针对赋值模型落点偏差ε_x和ε_z，其X,Z两个方向完全相同且相对独立，设置当 赋值模型落点偏差ε_x等于赋值模型X方向落点偏差的最大允许值ε_xr.max时，其无量纲隶属度 函数值S_x等于赋值模型X方向落点偏差无量纲隶属度函数的预设最小值S_x.min，其中S_x.min∈ [0,1]；当赋值模型落点偏差ε_x等于赋值模型X方向落点偏差的最小预设值ε_xr.min时，其无量 纲隶属度函数值S_x等于赋值模型X方向落点偏差无量纲隶属度函数的预设最大值S_x.max，其 中S_x.max∈[0,1]；针对赋值残差百分比当赋值残差百分比等于赋值残差百分比的最 大允许值ζ_r.max时，其无量纲隶属度函数值S_ζ等于赋值残差百分比无量纲隶属度函数的最小 预设值S_ζ.min，其中S_ζ.min∈[0,1]；当赋值残差百分比等于赋值残差百分比的最小允许值 ζ_r.min时，其无量纲隶属度函数值S_ζ等于赋值残差百分比无量纲隶属度函数的最大预设值 S_ζ.max，其中S_ζ.max∈[0,1]；

◆针对根据赋值速度选取的单项指标采用的无量纲隶属度函数如式(19)所示；

$\bar{R}=\left(\begin{array}{cc}{e}^{-k(R-R_0)},& R>R_0\\ 1,& R\le R_0\end{array}\right)---\left(19\right)$

式(19)中，表示赋值速度的各单项指标的无量纲化隶属度函数值，k和R₀为与各 单项指标的数值和评估要求相关的未知参数，R表示赋值速度的各单项指标自变量；针对每 一个单项指标，设计两个已知点位来确定未知参数k和R₀，从而得到该单项指标的无量纲隶 属度函数；针对扰动引力赋值时间T_c，当扰动引力赋值时间T_c等于扰动引力赋值时间的最大 允许值T_cr.max时，其无量纲隶属度函数值S_Tc等于扰动引力赋值时间无量纲隶属度函数的最 小预设值S_Tc.min，其中S_Tc.min∈[0,1]；当扰动引力赋值时间T_c小于或等于扰动引力赋值时间 的最小预设值T_cr.min时，其无量纲隶属度函数值S_Tc等于1；针对射前计算时间T_SQ，当射前计 算时间T_SQ等于射前计算时间的最大允许值T_SQr.max时，其无量纲隶属度函数值S_TSQ等于射前 计算时间无量纲隶属度函数的最小预设值S_TSQ.min，其中S_TSQ.min∈[0,1]；当射前计算时间T_SQ小于或等于射前计算时间的最小预设值T_SQr.min时，其无量纲隶属度函数值S_TSQ等于1；针对 机载计算机解算时间T_DS，当机载计算机解算时间T_DS等于机载计算机解算时间的最大允许 值T_DSr.max时，其无量纲隶属度函数值S_TDS等于机载计算机解算时间无量纲隶属度函数的最 小预设值S_TDS.min，其中S_TDS.min∈[0,1]，当机载计算机解算时间T_DS小于或等于机载计算机解 算时间的最小预设值T_DSr.min时，其无量纲隶属度函数值S_TDS等于1；

◆针对根据赋值模型适应性选取的单项指标采用的无量纲隶属度函数如式(20) 所示；

$\bar{R}=1-e^{-k{(R-R_0)}^2}---\left(20\right)$

式(20)中，表示无量纲化隶属度函数值，k和R₀为与各单项指标的数值和评估要 求相关的未知参数，R表示赋值模型适应性各单项指标自变量；针对每一个单项指标，设计 两个已知点位来确定未知参数k和R₀，从而得到该单项指标的无量纲隶属度函数；针对不同 飞行轨迹模式适应性当不同飞行轨迹模式适应性等于不同飞行轨迹模式适应 性最大预设值时，其无量纲隶属度函数值S_ε.traj等于不同飞行轨迹模式适应性无量 纲隶属度函数的最大预设值S_ε.traj.max，其中S_ε.traj.max∈[0,1]；当不同飞行轨迹模式适应性 等于不同飞行轨迹模式适应性最小预设值时，其无量纲隶属度函数值S_ε.traj等 于0；针对不同轨迹控制方式适应性当不同轨迹控制方式适应性等于不同轨迹 控制方式适应性最大预设值时，其无量纲隶属度函数值S_ε.ctl等于不同轨迹控制方式 适应性无量纲隶属度函数的最大预设值S_ε.ctl.max，其中S_ε.ctl.max∈[0,1]；当不同轨迹控制方 式适应性等于不同轨迹控制方式适应性最小预设值时，其无量纲隶属度函数值 S_ε.ctl等于0；针对不同地形适应性当不同地形适应性等于不同地形适应性的最 大预设值时，其无量纲隶属度函数值S_ε.hypso等于不同地形适应性无量纲隶属度函 数的最大预设值S_ε.hypso.max，其中S_ε.hypso.max∈[0,1]；当不同地形适应性等于不同地形 适应性的最小预设值时，其无量纲隶属度函数值S_ε.hypso等于0。

优选地，所述步骤6)的详细步骤包括：

6.1)针对各个单项指标的隶属度，根据基于各个单项指标作为底层，赋值精度、赋 值速度和赋值模型适应性三者作为中间层，待评估地球扰动引力赋值模型的适用性能评估 指标作为顶层的分层结构，基于各个单项指标的预设权重、各个中间层的预设权重采用分 层加权综合评估的方法最终得到待评估地球扰动引力赋值模型的适用性能评估指标A；

6.2)针对适用性能评估指标A进行模糊化定性处理，得到待评估地球扰动引力赋 值模型的适用性能评估结果。

优选地，所述步骤6.2)针对适用性能评估指标A进行模糊化定性处理具体是指采 用式(21)所示隶属度函数对适用性能评估指标A进行模糊化定性处理；

${\bar{F}}_R=1-e^{-k·A^2}---\left(21\right)$

式(21)中，A为待评估地球扰动引力赋值模型的适用性能评估指标，表示隶属 度函数值，k为模糊化定性处理参数。

本发明地球扰动引力赋值模型适用性能评估方法具有下述优点：本发明从赋值精 度、赋值速度和赋值模型适应性三方面设计用于适用性能评估的单项指标以建立地球扰动 引力赋值模型适用性能评估的指标体系；根据飞行轨迹样本的试验设计，针对待评估地球 扰动引力赋值模型进行弹道式飞行器飞行轨迹仿真，并计算各个单项指标的数据，得到各 个单项指标的数据样本库；假设各个单项指标的样本数据变化规律满足正态分布统计规 律，基于假设利用贝叶斯估计方法完成各个单项指标的评估分析；采用预设的无量纲隶属 度函数，对各个单项指标的评估分析结果进行无量纲化处理，得到各个单项指标的隶属度； 根据各个单项指标的预设权重及各个单项指标的隶属度，综合计算待评估地球扰动引力赋 值模型的适用性能评估结果，能够完成任意给定形式地球扰动引力赋值模型适用性能的评 估和分析以确定其在弹道式飞行器控制领域的适用性能，评估全面合理、快速可靠、实施便 捷的优点。

附图说明

图1为本发明实施例方法的基本流程示意图。

图2为本发明实施例中单项指标的分层结构示意图。

图3为本发明实施例步骤3)生成各个单项指标的数据样本库的流程示意图。

图4为本发明实施例针对单项指标的数据样本库的处理流程示意图。

具体实施方式

如图1所示，本实施例地球扰动引力赋值模型适用性能评估方法的步骤包括：

1)设置飞行轨迹样本试验因素及其试验水平，利用正交设计试验方法完成用于地 球扰动引力赋值模型的适用性能评估的飞行轨迹样本的试验设计；

2)从赋值精度、赋值速度和赋值模型适应性三方面设计用于适用性能评估的单项 指标以建立地球扰动引力赋值模型适用性能评估的指标体系；

3)根据飞行轨迹样本的试验设计，针对待评估地球扰动引力赋值模型进行弹道式 飞行器飞行轨迹仿真，并计算各个单项指标的数据，得到各个单项指标的数据样本库；

4)假设各个单项指标的样本数据变化规律满足正态分布统计规律，基于假设利用 贝叶斯估计方法完成各个单项指标的评估分析；

5)采用预设的无量纲隶属度函数，对各个单项指标的评估分析结果进行无量纲化 处理，得到各个单项指标的隶属度；

6)根据各个单项指标的预设权重及各个单项指标的隶属度，综合计算待评估地球 扰动引力赋值模型的适用性能评估结果。

地球扰动引力赋值模型适用性能是指该模型应用于弹道式飞行器控制系统设计 及仿真分析的适合程度；考虑到弹道式飞行器一般都会被赋予不同的飞行任务，相应的不 同的飞行任务也会对应不同的控制方式、不同的飞行轨迹模式和不同的飞行环境等，为充 分验证和分析地球扰动引力赋值模型在弹道式飞行器控制领域的适用性能，必须设计地球 扰动引力赋值模型适用性能评估需要的飞行轨迹样本。本实施例中，步骤1)中设置飞行轨 迹样本试验因素及其试验水平时，设置的飞行轨迹样本试验因素包括：弹道式飞行器的发 射位置地形、轨迹控制方式、航程、推进装置类型、飞行轨迹类型以及初始时刻飞行方向。本 实施例中，步骤1)中设置飞行轨迹样本试验因素及其试验水平时，弹道式飞行器的发射位 置地形的试验水平为4水平，轨迹控制方式的试验水平为3水平，航程的试验水平为6水平， 推进装置类型的试验水平为3水平，飞行轨迹类型的试验水平为4水平，初始时刻飞行方向 的试验水平为10水平。本实施例选取六项试验因素并利用正交设计试验方法来完成地球扰 动引力赋值模型适用性能评估飞行轨迹样本的试验设计，在明确飞行轨迹样本试验设计的 试验因素后，需要分别确定上述六项试验因素的试验水平才能利用正交设计试验方法来完 成飞行轨迹样本的试验设计；由于弹道式飞行器不同的飞行任务和多层面多侧重的研究目 的，本实施例对上述六项试验因素的试验水平不做具体的规定，使用者可根据研究需要自 行确定并利用正交设计试验方法完成地球扰动引力赋值模型适用性能评估飞行轨迹样本 的试验设计。

地球扰动引力赋值模型的赋值精度、赋值快速性和赋值模型适应性与弹道式飞行 器控制领域内的仿真分析有直接联系。赋值精度关系到弹道式飞行器的扰动引力受力准确 度，赋值快速性会影响弹道式飞行器机载计算机的在线计算能力，赋值模型适应性则直接 决定赋值模型的使用范围和对象。为此本发明中地球扰动引力赋值模型适用性能评估的指 标体系将围绕赋值精度、赋值快速性和赋值模型适应性三个方面建立。赋值精度指标要求 待评估的地球扰动引力赋值模型计算所得的空间任意位置扰动引力值与该位置的实际扰 动引力值的差异必须在可接受的范围内，只有赋值精度满足要求，才能保证扰动引力引起 的轨迹控制方法误差所造成的落点偏差在可接受的范围。由于弹道式飞行器准备时间和机 载计算机计算量的实际约束，必须保证待评估地球扰动引力赋值模型在可接受的时间内完 成计算，所以提出扰动引力赋值速度评估指标。地球扰动引力赋值模型适应性指标反映该 赋值模型可应用的范围和对象，只有适应性广泛的赋值模型才能被接受。如图2所示，本实 施例步骤2)中选取用于适用性能评估的单项指标时，根据赋值精度选取的单项指标包括： 扰动引力赋值残差δ_t、赋值模型落点偏差ε_x和ε_z、赋值残差百分比根据赋值速度选取的 单项指标包括：扰动引力赋值时间T_c、射前计算时间T_SQ、机载计算机解算时间T_DS；根据赋值 模型适应性选取的单项指标包括：不同飞行轨迹模式适应性不同轨迹控制方式适应 性不同地形适应性

本实施例中，扰动引力赋值残差δ_t的计算式如式(1)所示，赋值模型落点偏差ε_x和 ε_z的计算式如式(2)所示，赋值残差百分比的计算式如式(3)所示；扰动引力赋值时间T_c的计算式如式(4)所示，射前计算时间T_SQ的计算式如式(5)所示，机载计算机解算时间T_DS的 计算式如式(6)所示；不同飞行轨迹模式适应性的计算式如式(7)所示，不同轨迹控制 方式适应性的计算式如式(8)所示，不同地形适应性的计算式如式(9)所示；

${\delta }_t=\left|\right|g_r^t-g_c^t\left|\right|---\left(1\right)$

式(1)中，δ_t表示任意t时刻弹道式飞行器所在空间位置的扰动引力赋值残差，表示任意t时刻弹道式飞行器所在空间位置的实际扰动引力值，表示待评估地球扰动引 力赋值模型计算所得扰动引力值，表示实际扰动引力值与待评估地球扰动引力 赋值模型计算所得扰动引力值的差值的欧几里得范数；

$\left(\begin{array}{c}{ϵ}_x=|x_c-x_n|\\ {ϵ}_z=|z_c-z_n|\end{array}\right)---\left(2\right)$

式(2)中，以瞄准点为原点O、纵向和横向分别作为OX轴和OZ轴、落点偏离用(X,Z) 表示，ε_x表示X方向的落点偏差；ε_z表示Z方向的落点偏差；x_c表示利用待评估地球扰动引力 赋值模型进行飞行轨迹仿真所得落点的X向分量，x_n表示理想目标点位置坐标的X向分量，z_n表示理想目标点位置坐标的Z向分量，z_c表示利用待评估地球扰动引力赋值模型进行飞行 轨迹仿真所得落点的Z向分量；

${\zeta }_{\delta }^t={\delta }_t/\left|\right|g_r^t\left|\right|---\left(3\right)$

式(3)中，表示任意t时刻弹道式飞行器所在空间位置的赋值残差百分比，δ_t表 示任意t时刻弹道式飞行器所在空间位置的扰动引力赋值残差，表示任意t时刻弹道式 飞行器所在空间位置的实际扰动引力值的欧几里得范数；

$T_c=\underset{i=1}{\overset{n}{\Sigma }}{t}_c^i---\left(4\right)$

式(4)中，T_c表示扰动引力赋值时间，n表示调用待评估地球扰动引力赋值模型进 行扰动引力赋值计算的次数，表示第i次调用待评估地球扰动引力赋值模型进行扰动引 力赋值计算所用的时间；

T_SQ＝T_CX+T_ZY(5)

式(5)中，T_SQ表示射前计算时间，T_CX表示弹道式飞行器初始位置垂线偏差获取时 间，T_ZY表示弹道式飞行器诸元计算时间；

$T_{DS}=\underset{i=1}{\overset{n}{\Sigma }}{t}_{DS}^i---\left(6\right)$

式(4)中，T_DS表示机载计算机解算时间，n表示弹道式飞行器飞行轨迹仿真过程中 待评估地球扰动引力赋值模型被调用的次数，表示弹道式飞行器模拟机载计算机第i次 调用待评估地球扰动引力赋值模型进行扰动引力解算所需的时间；

${ϵ}_{traj}^{adapt}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}\frac{C_{traj,i}^{fit}}{C_{traj,i}^{all}}---\left(7\right)$

式(7)中，表示不同飞行轨迹模式适应性，表示针对第i种飞行轨迹模式 进行的仿真飞行轨迹总数目，表示总飞行轨迹仿真数目中满足终端要求的飞行轨迹数 目，n表示飞行轨迹模式数目，由设计者根据实际需要或弹道式飞行器具体的飞行轨迹确 定；

${ϵ}_{ctl}^{adapt}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}\frac{C_{ctl,i}^{fit}}{C_{ctl,i}^{all}}---\left(8\right)$

式(8)中，表示不同轨迹控制方式适应性，表示针对第i种轨迹控制方式 进行的仿真飞行轨迹总数目，表示总飞行轨迹仿真数目中满足轨迹控制要求的飞行轨 迹数目，n表示轨迹控制方式的数目；

${ϵ}_{hypso}^{adapt}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}\frac{C_{hypso,i}^{fit}}{C_{hypso,i}^{all}}---\left(9\right)$

式(9)中，表示不同地形适应性，表示针对第i种地形进行的仿真飞行 轨迹总数目，表示总飞行轨迹仿真数目中满足飞行任务要求的飞行轨迹数目，n表示 地形类型的数目。

如图3所示，步骤3)的详细步骤包括：

3.1)根据预设的空间边界约束或飞行器的飞行区域约束确定待评估地球扰动引 力赋值模型应用的空间范围约束；

3.2)针对空间范围约束的空间边界约束，结合待评估地球扰动引力赋值模型规定 的空间剖分准则完成待评估地球扰动引力赋值模型指定空域的空域剖分；

3.3)利用待评估地球扰动引力赋值模型根据空域剖分形成的空间有限元节点进 行扰动引力赋值，得到待评估地球扰动引力赋值模型的有限元节点数值模型；

3.4)在有限元节点数值模型的基础上，完成空间中任意位置扰动引力赋值的解算 和用于弹道式飞行器的飞行轨迹仿真；

3.5)在弹道式飞行器飞行轨迹仿真的过程中计算各个单项指标的数据，得到各个 单项指标的数据样本库。

其中，空间边界约束为使用者根据具体情况设置的空间范围；飞行区域约束为弹 道式飞行器实际飞行能力所确定的空间飞行区域范围；扰动引力赋值模型空间范围约束为 地球扰动引力赋值模型的可使用空间范围；空间剖分准则为地球扰动引力赋值模型进行有 限元划分所利用的空间剖分方法；扰动引力赋值模型指定空域剖分为利用空间剖分方法对 地球扰动引力赋值模型的空间使用范围进行空域剖分，得到该空间范围剖分后的有限元单 元、单元节点及其空间位置；地球扰动引力常用解算模型为目前地球扰动引力常用解算模 型，可以根据需要采用点质量计算方法、Stokes积分方法、球谐函数方法和球谐函数换极方 法；剖分空域特征节点扰动引力赋值为利用地球扰动引力常用解算模型对剖分后的有限元 单元各节点的扰动引力进行解算；基于空域特征节点的任意位置扰动引力赋值为利用剖分 后的有限元单元各节点的扰动引力值对指定空间内任意位置的扰动引力进行数值逼近计 算；各单项评估指标的样本数据库为根据本发明设计的各单项评估指标的计算公式，并结 合待评估的地球扰动引力赋值模型计算各单项评估指标的数值，进而形成单项评估指标的 样本数据库。

在利用待评估地球扰动引力赋值模型进行弹道式飞行器飞行轨迹仿真时，由于仿 真条件、赋值模型方法误差和不确定干扰等因素的影响，使得适用性能评估各单项指标的 不同时刻的计算数值呈现浮动变化。本实施例中，步骤4)中假设的正态分布统计规律为 其中，I_s.j表示第j个单项指标，μ_Is.j表示待评估地球扰动引力赋值模型 和飞行轨迹仿真过程中系统性误差决定的第j个单项指标的待估参数数值，表示 数据样本中第j个单项指标值相对单项指标I_s.j的散步程度，D_Is.j表示第j个单项指标的方 差，j＝0,1,…8，σ_Is.j表示第j个单项指标的标准差，第j个单项指标的方差D_Is.j、标准差σ_Is.j均为待估参数；鉴于以上分析，本实施例将地球扰动引力赋值模型适用性能单项指标的评 估问题转化为参数估计问题，进而利用贝叶斯估计方法完成上述假设条件下的单项指标评 估分析。利用贝叶斯估计方法完成各项单项指标的评估分析的步骤包括：

4.1)从所有单项指标中取出一项作为当前单项指标j；

4.2)在当前单项指标j的数据样本库中任意选取n_QF组仿真样本数据作 为验前子样，则当前单项指标j的待估参数μ_Is.j的验前概率密度为其中 表示仿真样本数据均值，仿真样本数据均值的计算式如式(10)所示；D_Is.j表示当前 单项指标j的方差，当前单项指标j的方差D_Is.j的验前概率密度为逆Gamma分布的密度函数g (D_x；α₀,β₀)，其中间参数α₀和β₀的计算式如式(11)所示，n_QF表示当前单项指标j的验前子样 的样本数目；

${\bar{x}}^{QF}=\frac{1}{n_{QF}}\underset{i=1}{\overset{n_{QF}}{\Sigma }}{x}_i^{QF}---\left(10\right)$

式(10)中，表示仿真样本数据均值，n_QF表示当前单项指标j的验前子样的样本 数目，表示当前单项指标j的验前子样中第i个样本的样本值；

$\left(\begin{array}{c}{\alpha }_0=\frac{1}{2}\underset{i=1}{\overset{n_{QF}}{\Sigma }}{(x_i^{QF}-{\bar{x}}^{QF})}^2\\ {\beta }_0=\frac{n_{QF}-1}{2}\end{array}\right)---\left(11\right)$

式(11)中，α₀和β₀分别表示逆Gamma分布的密度函数g(D_x；α₀,β₀)的中间参数，n_QF表 示当前单项指标j的验前子样的样本数目，表示当前单项指标j的验前子样中第i个样 本的样本值，表示当前单项指标j的仿真样本数据均值；

4.3)假定当前单项指标j的数据样本库中共有n个样本x₁,…,x_n，确定(μ_Is.j,D_Is.j) 的验后和验前是共轭的且服从正态—逆Gamma分布如式(12)所示；

$\pi ({\mu }_{Is.j},D_{Is.j}|X)\propto \frac{1}{\sqrt{D_{Is.j}}}{e}^{-\frac{1}{2{\eta }_1{D}_{Is.j}}{({\mu }_{Is.j}-{\mu }_1)}^2}{D_{Is.j}}^{-({\beta }_1+1)}{e}^{-\frac{{\alpha }_1}{D_{Is.j}}}---\left(12\right)$

式(12)中，μ_Is.j表示待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性 误差决定的当前单项指标j的待估参数数值，D_Is.j表示当前单项指标j的方差，X表示当前单 项指标j的数据样本库中n个样本x₁,…,x_n构成的数据集，中间参数η₁、α₁、β₁、μ₁如式(13)所 示；

$\left(\begin{array}{c}{\eta }_1=\frac{1}{n+n_{QF}}\\ {\alpha }_1={\alpha }_0+\frac{n}{2}·u+\frac{1}{2}·\frac{n{(\bar{x}-{\bar{x}}^{QF})}^2}{\frac{n}{n_{QF}}+1},u=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}{(x_i-\bar{x})}^2,\bar{x}=\frac{1}{n}\underset{i=1}{\overset{n}{\Sigma }}{x}_i\\ {\beta }_1={\beta }_0+\frac{n}{2}\\ {\mu }_1=\frac{n\bar{x}+n_{QF}{\bar{x}}^{QF}}{n+n_{QF}}\end{array}\right)---\left(13\right)$

式(13)中，η₁、α₁、β₁、μ₁分别为中间参数，n表示当前单项指标j的样本数据总数目， n_QF表示当前单项指标j的验前子样的样本数目，α₀,β₀为中间参数，u表示当前单项指标j的 数据样本库中样本数据的方差，x_i表示当前单项指标j的数据样本库中的第i个样本值，表示当前单项指标j的数据样本库中所有样本数据的均值，表示当前单项指标j的仿真 样本数据均值；

4.4)确定待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性误差决定 的当前单项指标j的待估参数数值μ_Is.j如式(14)所示；

${\mu }_{Is.j}={\left(\frac{{\alpha }_1{\eta }_1}{{\beta }_1}\right)}^{\frac{1}{2}}Y+{\mu }_1,p_Y\left(y\right)\propto {(1+\frac{y^2}{2{\beta }_1})}^{-\frac{2{\beta }_1+1}{2}}---\left(14\right)$

式(14)中，μ_Is.j表示待评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性 误差决定的当前单项指标j的待估参数数值，η₁、α₁、β₁、μ₁分别为式(13)所示的中间参数，Y 表示学生氏t分布，p_Y(y)表示Y的概率密度函数；

4.5)确定当前单项指标j的待估参数μ_Is.j的贝叶斯估计值及其方差Var(μ_Is.j) 如式(15)所示，当前单项指标j的方差D_Is.j的贝叶斯估计值及其方差Var(D_Is.j)如式 (16)所示，当前单项指标j的标准差σ_Is.j的贝叶斯估计值及其方差Var(σ_Is.j)如式(17) 所示；

$\left(\begin{array}{c}{\hat{\mu }}_{Is.j}={\mu }_1\\ Var\left({\mu }_{Is.j}\right)=\frac{{\alpha }_1{\eta }_1}{{\beta }_1-1}\end{array}\right)---\left(15\right)$

式(15)中，表示当前单项指标j的待估参数μ_Is.j的贝叶斯估计值，μ_Is.j表示待 评估地球扰动引力赋值模型和飞行轨迹仿真过程中系统性误差决定的当前单项指标j的待 估参数数值，Var(μ_Is.j)表示贝叶斯估计值的方差，参数η₁、α₁、β₁、μ₁分别表示式(13)所 示的中间参数；

$\left(\begin{array}{c}{\hat{D}}_{Is.j}=\frac{{\alpha }_1}{{\beta }_1-1}\\ Var\left(D_{Is·j}\right)=\frac{{\alpha }_1^2}{{({\beta }_1-1)}^2({\beta }_1-2)}\end{array}\right)---\left(16\right)$

式(16)中，表示当前单项指标j的方差D_Is.j的贝叶斯估计值，D_Is.j表示当前单 项指标j的方差，Var(D_Is.j)表示当前单项指标j的方差D_Is.j的贝叶斯估计值，参数η₁、α₁、β₁、 μ₁分别为式(13)所示的中间参数；

$\{\begin{array}{c}{\hat{\sigma }}_{Is.j}=\sqrt{{\alpha }_1}\frac{\Gamma ({\beta }_1-1/2)}{\Gamma \left({\beta }_1\right)}\\ Var\left({\sigma }_{Is.j}\right)=\frac{{\alpha }_1}{{\beta }_1-1}-{\alpha }_1{\left(\frac{\Gamma ({\beta }_1-1/2)}{\Gamma \left({\beta }_1\right)}\right)}^2\end{array}---\left(17\right)$

式(17)中，表示当前单项指标j的标准差σ_Is.j的贝叶斯估计值，σ_Is.j表示当前 单项指标j的标准差，Var(σ_Is.j)表示贝叶斯估计值的方差，参数η₁、α₁、β₁、μ₁分别为式 (13)所示的中间参数，Γ为伽玛函数；

4.6)根据当前单项指标j的仿真样本数据均值待估参数μ_Is.j的贝叶斯估计值 及其方差Var(μ_Is.j)、方差D_Is.j的贝叶斯估计值及其方差Var(D_Is.j)、标准差σ_Is.j的 贝叶斯估计值及其Var(σ_Is.j)，对当前单项指标j进行定量或者定性评估；获取第j个单 项指标的均值、方差和标准差的估计值后，参照实际飞行任务中的定量指标即可对每个单 项指标进行评估分析；即便没有定量指标要求，也可根据经验或者专家意见得到各单项指 标的定性评估结论；

4.7)判断是否所有单项指标已经处理完毕，如果尚未处理完毕，则跳转执行步骤 4.1)，否则判断利用贝叶斯估计方法完成各项单项指标的评估分析完毕，跳转执行步骤5)。

考虑到每个单项评估指标的单项评估流程基本一致，这里只以“赋值模型落点偏 差”单项指标为例进行详细的计算分析，其他8个单项指标的贝叶斯估计方法与分析流程和 该单项指标类似，这里不做赘述。根据实际要求设置赋值模型落点偏差单项指标的硬性可 接受量化范围为200～500m，则根据赋值模型落点偏差单项指标数据样本对该单项指标进 行评估分析。首先从赋值模型落点偏差指标的数据样本中选取50组数据作为验前子样，如 表1所示。

表1：赋值模型落点偏差指标的验前子样数据

根据贝叶斯估计方法得到X,Z方向的赋值模型落点偏差的贝叶斯估计中间参数及 估计结果如表2所示。

表2：赋值模型落点偏差的贝叶斯估计计算结果。

根据贝叶斯估计结果，利用该地球扰动引力赋值模型进行飞行轨迹仿真，可使得 弹道式飞行器的X向和Z向的落点偏差均值分别为358.1142m和259.6198m，两个方向落点偏 差的标准差在百米量级，根据事先设置的落点偏差的硬性指标可知，该扰动引力赋值模型 的赋值模型落点偏差指标是满足要求的。

与赋值模型落点偏差单项评估指标相同，其它8个单项评估指标的单项评估均可 以利用相同的方法实施；根据其它8个单项评估指标的样本数据，采用相同的单项评估流 程，计算得到扰动引力赋值残差、赋值残差百分比、扰动引力赋值时间、射前计算时间、机载 计算机解算时间、不同弹道模式适应性、不同轨迹控制方式适应性和不同地形适应性8个单 项评估指标以及赋值模型落点偏差单项评估指标的贝叶斯评估结果如表3所示。

表3：各单项评估指标的贝叶斯评估结果。

各单项评估指标的评估可根据表3中所示的各单项指标的贝叶斯估计结果结合各 单项评估指标的具体要求进行对比得出具体的评估结论。

衡量地球扰动引力赋值模型的各单项指标属性不一、性质各异且各自的量纲不统 一，直接采用统一的度量方式进行适用性能综合评估缺乏合理性；为此，本实施例中步骤5) 将采用模糊综合评定方法将各单项指标进行无量纲化处理。本实施例中，步骤5)对各个单 项指标的评估分析结果进行无量纲化处理时：

◆针对根据赋值精度选取的单项指标采用的无量纲隶属度函数如式(18)所示；

$\bar{R}=e^{-k(R-R_0)},R>R_0---\left(18\right)$

式(18)中，表示赋值精度的各单项指标的无量纲化隶属度函数值，k和R₀为与各 单项指标的数值和评估要求相关的未知参数，R表示赋值精度的各单项指标自变量；针对每 一个单项指标，设计两个已知点位来确定未知参数k和R₀，从而得到该单项指标的无量纲隶 属度函数；针对扰动引力赋值残差δ_t，设置当扰动引力赋值残差δ_t等于扰动引力赋值残差的 最大允许值δ_r.max时，其无量纲隶属度函数值S_δ等于赋值残差无量纲隶属度函数的预设最小 值S_δ.min，其中S_δ.min∈[0,1]；当扰动引力赋值残差δ_t等于扰动引力赋值残差的最小预设值 δ_r.min时，无量纲隶属度函数值S_δ等于赋值残差无量纲隶属度函数的预设最大值S_δ.max，其中 S_δ.max∈[0,1]；针对赋值模型落点偏差ε_x和ε_z，其X,Z两个方向完全相同且相对独立，设置当 赋值模型落点偏差ε_x等于赋值模型X方向落点偏差的最大允许值ε_xr.max时，其无量纲隶属度 函数值S_x等于赋值模型X方向落点偏差无量纲隶属度函数的预设最小值S_x.min，其中S_x.min∈ [0,1]；当赋值模型落点偏差ε_x等于赋值模型X方向落点偏差的最小预设值ε_xr.min时，其无量 纲隶属度函数值S_x等于赋值模型X方向落点偏差无量纲隶属度函数的预设最大值S_x.max，其 中S_x.max∈[0,1]；针对赋值残差百分比当赋值残差百分比等于赋值残差百分比的最 大允许值ζ_r.max时，其无量纲隶属度函数值S_ζ等于赋值残差百分比无量纲隶属度函数的最小 预设值S_ζ.min，其中S_ζ.min∈[0,1]；当赋值残差百分比等于赋值残差百分比的最小允许值 ζ_r.min时，其无量纲隶属度函数值S_ζ等于赋值残差百分比无量纲隶属度函数的最大预设值 S_ζ.max，其中S_ζ.max∈[0,1]；

◆针对根据赋值速度选取的单项指标采用的无量纲隶属度函数如式(19)所示；

$\bar{R}=\left(\begin{array}{cc}{e}^{-k(R-R_0)},& R>R_0\\ 1,& R\le R_0\end{array}\right)---\left(19\right)$

式(19)中，表示赋值速度的各单项指标的无量纲化隶属度函数值，k和R₀为与各 单项指标的数值和评估要求相关的未知参数，R表示赋值速度的各单项指标自变量；针对每 一个单项指标，设计两个已知点位来确定未知参数k和R₀，从而得到该单项指标的无量纲隶 属度函数；针对扰动引力赋值时间T_c，当扰动引力赋值时间T_c等于扰动引力赋值时间的最大 允许值T_cr.max时，其无量纲隶属度函数值S_Tc等于扰动引力赋值时间无量纲隶属度函数的最 小预设值S_Tc.min，其中S_Tc.min∈[0,1]；当扰动引力赋值时间T_c小于或等于扰动引力赋值时间 的最小预设值T_cr.min时，其无量纲隶属度函数值S_Tc等于1；针对射前计算时间T_SQ，当射前计 算时间T_SQ等于射前计算时间的最大允许值T_SQr.max时，其无量纲隶属度函数值S_TSQ等于射前 计算时间无量纲隶属度函数的最小预设值S_TSQ.min，其中S_TSQ.min∈[0,1]；当射前计算时间T_SQ小于或等于射前计算时间的最小预设值T_SQr.min时，其无量纲隶属度函数值S_TSQ等于1；针对 机载计算机解算时间T_DS，当机载计算机解算时间T_DS等于机载计算机解算时间的最大允许 值T_DSr.max时，其无量纲隶属度函数值S_TDS等于机载计算机解算时间无量纲隶属度函数的最 小预设值S_TDS.min，其中S_TDS.min∈[0,1]，当机载计算机解算时间T_DS小于或等于机载计算机解 算时间的最小预设值T_DSr.min时，其无量纲隶属度函数值S_TDS等于1；

◆针对根据赋值模型适应性选取的单项指标采用的无量纲隶属度函数如式(20) 所示；

$\bar{R}=1-e^{-k{(R-R_0)}^2}---\left(20\right)$

式(20)中，表示无量纲化隶属度函数值，k和R₀为与各单项指标的数值和评估要 求相关的未知参数，R表示赋值模型适应性各单项指标自变量；针对每一个单项指标，设计 两个已知点位来确定未知参数k和R₀，从而得到该单项指标的无量纲隶属度函数；针对不同 飞行轨迹模式适应性当不同飞行轨迹模式适应性等于不同飞行轨迹模式适应 性最大预设值时，其无量纲隶属度函数值S_ε.traj等于不同飞行轨迹模式适应性无量 纲隶属度函数的最大预设值S_ε.traj.max，其中S_ε.traj.max∈[0,1]；当不同飞行轨迹模式适应性 等于不同飞行轨迹模式适应性最小预设值时，其无量纲隶属度函数值S_ε.traj等 于0；针对不同轨迹控制方式适应性当不同轨迹控制方式适应性等于不同轨迹 控制方式适应性最大预设值时，其无量纲隶属度函数值S_ε.ctl等于不同轨迹控制方式 适应性无量纲隶属度函数的最大预设值S_ε.ctl.max，其中S_ε.ctl.max∈[0,1]；当不同轨迹控制方 式适应性等于不同轨迹控制方式适应性最小预设值时，其无量纲隶属度函数值 S_ε.ctl等于0；针对不同地形适应性当不同地形适应性等于不同地形适应性的最 大预设值时，其无量纲隶属度函数值S_ε.hypso等于不同地形适应性无量纲隶属度函 数的最大预设值S_ε.hypso.max，其中S_ε.hypso.max∈[0,1]；当不同地形适应性等于不同地形 适应性的最小预设值时，其无量纲隶属度函数值S_ε.hypso等于0。

本实施例中根据表3所示各单项评估指标的贝叶斯估计结果，结合各单项评估指 标的样本数据和本发明提供的各单项评估指标的模糊隶属度函数，对各单项评估指标进行 无量纲化处理，各单项指标的无量纲化隶属度函数具体设计如下：

◆扰动引力赋值残差。设置当扰动引力赋值残差δ_t＝30mgal时，S_δ＝0.05；当δ_t＝ 0mgal时，S_δ＝1，其无量纲隶属度函数为式中S_δ为扰动引力赋值残差单项指标 的无量纲模糊隶属度函数值，δ_t为扰动引力赋值残差。

◆赋值模型落点偏差。设置当赋值模型落点偏差ε_x＝500m,ε_z＝500m时，S_x＝0.01, S_z＝0.01；当ε_x＝50m,ε_z＝50m时，S_x＝1,S_z＝1，其无量纲隶属度函数为

式中S_x,S_z分别为X,Z方向赋值模型落点偏差单项指标的 无量纲模糊隶属度函数值，ε_x,ε_z分别为X,Z方向赋值模型落点偏差。

◆赋值残差百分比。设置当赋值残差百分比时，S_ζ＝0.02；当时，S_ζ＝1；其无量纲隶属度函数为S_ζ＝e^-0.1956ζ，式中S_ζ为赋值残差百分比单项指标的无量纲隶属 度函数值，为赋值残差百分比。

◆扰动引力赋值时间。设置当扰动引力赋值时间T_c＝50ms时，S_Tc＝0.1；当扰动引 力赋值时间T_c≤3ms时，S_Tc＝1；其无量纲隶属度函数为：

$S_{Tc}=\left(\begin{array}{c}{e}^{-0.049(Tc-3)},T_c>T_{c\min }\\ 1,T_c\le T_{c\min }\end{array}\right)$

式中S_Tc为扰动引力赋值时间单项指标的无量纲隶属度函数值，T_c为扰动引力赋值 时间。

◆射前计算时间。设置当射前计算时间T_SQ＝5分钟时，S_TSQ＝0.01；当射前计算时间 T_SQ＝0.5分钟时，S_TSQ＝1；其无量纲隶属度函数为：

$S_{T_{SQ}}=\left(\begin{array}{c}{e}^{-1.0234(T_{SQ}-0.5)},T_{SQ}>T_{SQ}^{\min }\\ 1,T_{SQ}\le T_{SQ}^{\min }\end{array}\right)$

式中S_TSQ为射前计算时间单项指标的无量纲隶属度函数值，T_SQ为射前计算时间。

◆机载计算机解算时间。设置当机载计算机解算时间T_DS＝10ms时，S_TDS＝0.05；当 机载计算机解算时间T_DS＝0.5ms时，S_TDS＝1；其无量纲隶属度函数为：

$S_{T_{DS}}=\left(\begin{array}{c}{e}^{-0.3153(T_{DS}-0.5)},T_{DS}>0.5ms\\ 1,T_{DS}\le 0.5ms\end{array}\right)$

式中S_TDS为机载计算机解算时间单项指标的无量纲隶属函数值，T_DS为机载计算机 解算时间。

◆不同飞行轨迹模式适应性。设置当不同飞行轨迹模式适应性时， S_ε.traj＝0；当不同飞行轨迹模式适应性时，S_ε.traj＝0.99；其无量纲隶属度函数 为：

$S_{ϵ.traj}=1-e^{-8.187\times {10}^{-4}}{({ϵ}_{traj}^{adapt}-20)}^2$

式中S_ε.traj为不同飞行轨迹模式适应性单项指标的无量纲隶属度函数值，为 不同飞行轨迹模式适应性。

◆不同轨迹控制方式适应性。设置当不同轨迹控制方式适应性时， S_ε.ctl＝0；当不同轨迹控制方式适应性时，S_ε.ctl＝0.99；其无量纲隶属度函数 为：

$S_{ϵ.ctl}=1-e^{-8.187\times {10}^{-4}}{({ϵ}_{ctl}^{adapt}-20)}^2$

式中S_ε.ctl为不同轨迹控制方式适应性单项指标的无量纲隶属度函数值，为 不同轨迹控制方式适应性。

◆不同地形适应性。设置当不同地形适应性时，S_ε.hypso＝0；当不同地 形适应性时，S_ε.hypso＝0.99；其无量纲隶属度函数为

$S_{ϵ.hypso}=1-e^{-8.187\times {10}^{-4}}{({ϵ}_{hypso}^{adapt}-20)}^2$

式中S_ε.hypso为不同地形适应性单项指标的无量纲隶属度函数值，为不同地形 适应性。

扰动引力赋值残差、赋值模型落点偏差、赋值残差百分比、扰动引力赋值时间、射 前计算时间、机载计算机解算时间、不同飞行轨迹模式适应性、不同轨迹控制方式适应性和 不同地形适应性九个单项指标的无量纲化处理结果如表4所示。

表4：各单项指标无量纲化结果

地球扰动引力赋值模型适用性能综合评估的核心思想是对各单项指标进行分层 加权，不同层级指标的权重可由各指标对总体指标的贡献值和专家打分方法获得，结合无 量纲化处理后的单项指标的量值可得到适用性能总体指标的量值，而后进行模糊化处理得 到待评估地球扰动引力赋值模型适用性能的定性评估结论；需要说明，本实施例只是将分 层加权评估方法创新性地应用于地球扰动引力赋值模型适用性能评估问题，并未对分层加 权评估方法本身进行完善和补充，地球扰动引力赋值模型适用性能综合评估方法的具体思 路如图4所示。地球扰动引力赋值模型适用性能综合评估整体采用两条并行流程：1)单项指 标权重、中间指标权重的获取，适用性能评估各层级评估指标的权重是在充分考虑各单项 指标的重要性和其对总体指标贡献值的基础上，通过专家打分和具体客观的评估要求完成 设计的；2)本实施例单项指标样本数据库的处理。如图4所示，步骤6)的详细步骤包括：

6.1)针对各个单项指标的隶属度，根据基于各个单项指标作为底层，赋值精度、赋 值速度和赋值模型适应性三者作为中间层，待评估地球扰动引力赋值模型的适用性能评估 指标作为顶层的分层结构(如图2所示)，基于各个单项指标的预设权重、各个中间层的预设 权重采用分层加权综合评估的方法最终得到待评估地球扰动引力赋值模型的适用性能评 估指标A；

6.2)针对适用性能评估指标A进行模糊化定性处理，得到待评估地球扰动引力赋 值模型的适用性能评估结果。

本实施例中，各个单项指标的预设权重、各个中间层的预设权重如表5所示。

表5：各个单项指标的预设权重、各个中间层的预设权重设置表。

本实施例中，步骤6.2)针对适用性能评估指标A进行模糊化定性处理具体是指采 用式(21)所示隶属度函数对适用性能评估指标A进行模糊化定性处理；

${\bar{F}}_R=1-e^{-k·A^2}---\left(21\right)$

式(21)中，A为待评估地球扰动引力赋值模型的适用性能评估指标，表示隶属 度函数值，k为模糊化定性处理参数。本实施例中，模糊化定性处理参数k＝2.5；设置当 时，评估结果定义为“极差”；设置当时，评估结果定义为“差”；当 时，评估结果定义为“可以接受”；当时，评估结果定义为“良好”； 当时，评估结果定义为“优秀”。结合表4和表5，根据各单项指标无量纲化结果 和各单项指标的设计权重，利用分层加权方法计算得到的沿飞行轨迹地球扰动引力赋值模 型的适用性能指标值为0.7993，根据适用性能模糊定性处理结果可知，该地球扰动引力赋 值模型应用于弹道式飞行器控制领域和飞行轨迹仿真的评估结果为“良好”。沿飞行轨迹地 球扰动引力赋值模型用于弹道式飞行器仿真分析是可以接受的，可以满足弹道式飞行器仿 真计算所要求的扰动引力赋值精度、赋值快速性和赋值模型适应性的基本要求。

综上所述，本实施例基于上述六个步骤即可进行任意给定形式地球扰动引力赋值 模型适用性能的评估和分析，使用者可根据需要对地球扰动引力赋值模型进行全面、合理、 快速和可靠评估与分析，进而确定其在弹道式飞行器控制领域的适用性能。总体来讲，本发 明提出的地球扰动引力赋值模型适用性能评估方法具有以下特点：1)首次提出以赋值精 度、赋值快速性和赋值模型适应性作为核心指标并进行地球扰动引力赋值模型的适用性能 评估，规范了地球扰动引力赋值模型适用性能评估的评估流程，设计了合理全面且实施便 捷的评估方案；2)设计地球扰动引力赋值模型适用性能评估的指标体系，给出各单项指标 的定义并确定各单项指标的计算模型；3)提出地球扰动引力赋值模型适用性能评估“六步 三层”实施方案，“六步”是指本发明给出的步骤1)～步骤6)规定的具体操作步骤，“三层”是 指单项指标样本数据生成、单项指标贝叶斯估计和适用性能综合评估。本发明针对地球扰 动引力赋值模型用于弹道式飞行器仿真分析的适用性能评估问题，设计合理且涵盖全面的 适用性能评估指标体系，采用单项评估指标和适用性能综合指标顺序评估分析的思路依次 完成单项指标的评估分析和适用性能的综合评估；评估方法结合试验设计方法，贝叶斯估 计方法和模糊隶属方法和分层加权方法提出一种思路清晰、实施便捷、结论严谨的地球扰 动引力赋值模型适用性能评估方法，为地球扰动引力赋值模型在弹道式飞行器仿真分析领 域的应用提供基本依据。

以上所述仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施 例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域 的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，这些改进和润饰也 应视为本发明的保护范围。