利用部件质量特性计算卫星腔体内万有引力的等效质点组方法

摘要

一种利用部件质量特性计算卫星腔体内万有引力的等效质点组方法，1，根据任务需求，设计确定卫星上的部件；2，利用质量特性综合测量仪，测量步骤1所确定的部件中的某一部件的质量特性并记录；3，根据步骤2的测量结果，建立等效质点组满足的数学方程；4，根据步骤3的方程，求解得到所需要的等效质点组；5，根据步骤4所求解的等效质点组，计算该部件的万有引力及其梯度；6，重复2-5，计算出所有部件的万有引力及其梯度，求和得到卫星在腔体内的万有引力及其梯度。本发明在计算卫星上质量分布复杂部件的万有引力时，能够利用实际可测的质量特性信息，以二阶精度实现万有引力的计算，且不受部件形状的影响，在工程上有较广的适用范围。

说明书

技术领域

本发明涉及航天动力学技术领域，尤其是涉及一种利用部件质量特性计算卫星腔体内万 有引力的等效质点组方法。

背景技术

一些以基础物理实验为任务目标的空间任务，如探测引力波和检验广义相对论的LISA和 ASTROD任务，需要位于卫星腔体内的参考质量沿着纯引力轨道飞行(参见期刊《经典与量子 引力》(Classical and Quantum Gravity)2003年第20卷的文章“LISA的集成模型(The LISA integrated model)”和期刊《原了核物理B》(Nuclear Physics B)2007年第166卷153-158 页文章“ASTROD(激光天文动力学)and ASTROD I”)。同样，利用参考质量沿着近地纯引力 轨道飞行，并获取参考质量的纯引力轨道，能够用于精确测量地球重力场(参见期刊《国际 宇航联会刊》(Acta Astronautica)2012年特刊文章“采用精密编队飞行技术获取纯引力轨 道(Acquirement of pure gravity orbit using precision formation flying technology)”)。 对于这些任务的科学目标而言，卫星作用在腔体内部参考质量上的万有引力是一个主要的干 扰力，影响纯引力轨道的性能水平(参见《经典与量子引力》(Classical and Quantum Gravity) 2004年第21卷第5期S653-S660页的文章“当前的LISA残余加速度误差估计(Current error estimates for LISA spurious accelerations)”)。精确计算纯引力轨道的万有引力干扰， 是克服该项干扰从而提高纯引力轨道性能的基础。

现有技术中，LISA模型团队建立了万有引力干扰的数值计算方法，采用卫星有限单元模 型提供的结点质量和位置，并将每个单元近似为质点计算其对参考质量的引力、力矩和梯度 作用，然后对所有单元求和得到整体量(参见《经典与量子引力》(Classical and Quantum Gravity)2005年第22卷第10期S395-S402页的文章“LISA自引力分析模型(Self-gravity modeling for LISA)”)。

但是，这种利用卫星模型计算万有引力干扰的方法，对卫星有限单元模型的质量分布精 度要求苛刻，使得材质非均匀部件的精确建模非常困难。在更广泛的纯引力飞行任务中，可 以采用本发明所提出的方法利用实际可测的部件质量特性精确计算万有引力干扰。

发明内容

本发明的目的在于设计一种新型的利用部件质量特性计算卫星腔体内万有引力的等效质 点组方法，解决上述问题。

为了实现上述目的，本发明采用的技术方案如下：

一种利用部件质量特性计算卫星腔体内万有引力的等效质点组方法，包括步骤如下：

步骤1，根据任务需求，设计确定卫星上的部件；

步骤2，利用质量特性综合测量仪，测量步骤1所确定的部件中的某一部件的质量特性 并记录；

步骤3，根据步骤2的测量结果，建立等效质点组满足的数学方程；

步骤4，根据步骤3的方程，求解得到所需要的等效质点组；

步骤5，根据步骤4所求解的等效质点组，计算该部件的万有引力及其梯度；

步骤6，重复步骤2-5，计算出所有部件的万有引力及其梯度，求和得到卫星在腔体内的 万有引力及其梯度。

优选的，所述步骤3中所示的根据步骤2的测量结果，建立等效质点组满足的数学方程， 具体包括：

设步骤二中测得某一部件B的质量为M_o，在参考直角坐标系OXYZ中质心位置为 o(X_o，Y_o，Z_o)，在该部件主轴坐标系oyxz中的转动惯量矩阵为$\left(\begin{array}{ccc}{I}_{\mathrm{xx},o}& 0& 0\\ 0& I_{\mathrm{yy},o}& 0\\ 0& 0& I_{\mathrm{zz},o}\end{array}\right);$

设与所述部件B质量特性相同的等效质点组包含N_PG个质点，质点的质量记为 m_PG，i，i＝1，2，…N_PG，在该部件主轴坐标系oyxz中的坐标依次为(x_PG，i，y_PG，i，z_PG，i)，i＝1，2，…N_PG， 则等效质点组满足的数学方程为：

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}=M_o$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}(y_{\mathrm{PG},i}^2+z_{\mathrm{PG},i}^2)=I_{\mathrm{xx},o},\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}(x_{\mathrm{PG},i}^2+z_{\mathrm{PG},i}^2)=I_{\mathrm{yy},o},\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{P,i}(x_{\mathrm{PG},i}^2+y_{\mathrm{PG},i}^2)=I_{\mathrm{zz},o}$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0.---\left(1\right)$

优选的，所述步骤5中所示的，根据步骤4所求解的等效质点组，计算该部件的万有引 力及其梯度，具体包括：

等效质点组中质点的坐标(x_PG，i，y_pG，i，z_PG，i)，i＝1，2，…N_PG是在部件B的主轴坐标系中给出 的，部件B的主轴坐标系到参考直角坐标系的旋转矩阵为R_ob，则等效质点组在参考直角坐 标系中的的坐标为

$(X_{\mathrm{PG},b,i},Y_{\mathrm{PG},b,i},Z_{\mathrm{PG},b,i})=(x_{\mathrm{PG},i},y_{\mathrm{PG},i},z_{\mathrm{PG},i})R_{\mathrm{ob}}^T+(X_o,Y_o,Z_o)---\left(2\right)$

记录等效质点组在参考直角坐标系中的以上坐标，则可以始终利用等效质点组计算部件 B在腔体内一点(X₀，Y₀，Z₀)的万有引力：

$a_{\mathrm{PG},x}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{X_{\mathrm{PG},b,i}-X_0}{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2{]}^{3/2}}$

$a_{\mathrm{PG},y}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{Y_{\mathrm{PG},b,i}-Y_0}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{3/2}}---\left(3\right)$

$a_{\mathrm{PG},z}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{Z_{\mathrm{PG},b,i}-Z_0}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{3/2}}$

同样，可以利用等效质点组计算部件B的万有引力梯度。

$T=\left(\begin{array}{ccc}{V}_{\mathrm{xx}}& V_{\mathrm{xy}}& V_{\mathrm{xz}}\\ V_{\mathrm{yx}}& V_{\mathrm{yy}}& V_{\mathrm{yz}}\\ V_{\mathrm{zx}}& V_{\mathrm{zy}}& V_{\mathrm{zz}}\end{array}\right)---\left(4\right)$

$V_{\mathrm{xx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{2{(X_{\mathrm{PG},b,i}-X_0)}^2-{(Y_{\mathrm{PG},b,i}-Y_0)}^2-{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{yy}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2-{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{zz}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-{(X_{\mathrm{PG},b,i}-X_0)}^2-{(Y_{\mathrm{PG},b,i}-Y_0)}^2+2{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{xy}}=V_{\mathrm{yx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(X_{\mathrm{PG},b,i}-X_0)(Y_{\mathrm{PG},b,i}-Y_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{xz}}=V_{\mathrm{zx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(X_{\mathrm{PG},b,i}-X_0)(Z_{\mathrm{PG},b,i}-Z_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{yz}}=V_{\mathrm{zy}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(Y_{\mathrm{PG},b,i}-Y_0)(Z_{\mathrm{PG},b,i}-Z_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}.$

本发明所谓的质量特性综合测量仪，是指在卫星工程中，测量卫星或某部件的质量、质 心位置、转动惯量的设备，

本发明的目的是在精确计算卫星上复杂质量分布部件的万有引力时，利用部件的质量特 性以二阶精度实现万有引力计算。

本发明利用部件质量特性计算卫星腔体内万有引力的等效质点组方法包括：卫星部件选 型设计、测量部件的质量特性、求解与部件质量特性相同的等效质点组、利用等效质点组计 算部件的万有引力、加和得到卫星在腔体内的万有引力。

所述的卫星部件选型设计，是根据卫星任务需要，确定要使用的部件。所述的测量部件 的质量特性，是利用能够测量质量、质心和转动惯量的综合测试设备测量得到部件的质量、 质心和转动惯量这些质量特性值。所述的求解与部件质量特性相同的等效质点组，是根据已 经测量得到的部件质量特性值，建立等效质点组必须满足的数学方程并求解得到需要的等效 质点组。所述的利用等效质点组计算部件的万有引力，是根据求解得到的等效质点组，计算 其中每个质点的万有引力及其梯度，加和得到等效质点组也就是部件的万有引力及其梯度。 所述的加和得到卫星在腔体内的万有引力，是对所有涉及到的部件，将它们的万有引力及其 梯度求和，得到卫星在腔体内的万有引力及其梯度。

本发明的有益效果可以总结如下：

1，本发明在计算卫星上质量分布复杂部件的万有引力时，能够利用实际可测的质量 特性信息，以二阶精度实现万有引力的计算，且不受部件形状的影响，在工程上 有较广的适用范围。

2，本发明方法简单、实施成本低廉，计算结果精确。

附图说明

图1.参考直角坐标系OXYZ和某一部件B的主轴坐标系oyxz，其中o为部件B在参考直 角坐标系中的质心位置，x，y，z向分别为部件B的三个惯量主轴。从oyxz到OXYZ的旋转矩阵为 R_ob。

具体实施方式

为了使本发明所解决的技术问题、技术方案及有益效果更加清楚明白，以下结合附图及 实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释 本发明，并不用于限定本发明。

一种利用部件质量特性计算卫星腔体内万有引力的等效质点组方法，包括步骤如下：

步骤1，根据任务需求，设计确定卫星上的部件；

步骤2，利用质量特性综合测量仪，测量步骤1所确定的部件中的某一部件的质量特性 并记录；

步骤3，根据步骤2的测量结果，建立等效质点组满足的数学方程；

步骤4，根据步骤3的方程，求解得到所需要的等效质点组；

步骤5，根据步骤4所求解的等效质点组，计算该部件的万有引力及其梯度；

步骤6，重复步骤2-5，计算出所有部件的万有引力及其梯度，求和得到卫星在腔体内的 万有引力及其梯度。

其中，所述步骤3中所示的根据步骤2的测量结果，建立等效质点组满足的数学方程， 具体包括：

设步骤二中测得某一部件B的质量为M_o，在参考直角坐标系OXYZ中质心位置为 o(X_o，Y_o，Z_o)，在该部件主轴坐标系oyxz中的转动惯量矩阵为$\left(\begin{array}{ccc}{I}_{\mathrm{xx},o}& 0& 0\\ 0& I_{\mathrm{yy},o}& 0\\ 0& 0& I_{\mathrm{zz},o}\end{array}\right);$

设与所述部件B质量特性相同的等效质点组包含N_PG个质点，质点的质量记为 m_PG，i，i＝1，2，…N_PG，在该部件主轴坐标系oyxz中的坐标依次为(x_PG，i，y_PG，i，z_PG，i)，i＝1，2，…N_PG， 则等效质点组满足的数学方程为：

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}=M_o$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}(y_{\mathrm{PG},i}^2+z_{\mathrm{PG},i}^2)=I_{\mathrm{xx},o},\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}(x_{\mathrm{PG},i}^2+z_{\mathrm{PG},i}^2)=I_{\mathrm{yy},o},\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{P,i}(x_{\mathrm{PG},i}^2+y_{\mathrm{PG},i}^2)=I_{\mathrm{zz},o}$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0.---\left(1\right)$

其中，所述步骤5中所示的，根据步骤4所求解的等效质点组，计算该部件的万有引力 及其梯度，具体包括：

等效质点组中质点的坐标(x_PG，i，y_PG，i，z_PG，i)，i＝1，2，…N_PG是在部件B的主轴坐标系中给出 的，部件B的主轴坐标系到参考直角坐标系的旋转矩阵为R_ob，则等效质点组在参考直角坐 标系中的的坐标为

$(X_{\mathrm{PG},b,i},Y_{\mathrm{PG},b,i},Z_{\mathrm{PG},b,i})=(x_{\mathrm{PG},i},y_{\mathrm{PG},i},z_{\mathrm{PG},i})R_{\mathrm{ob}}^T+(X_o,Y_o,Z_o)---\left(2\right)$

记录等效质点组在参考直角坐标系中的以上坐标，则可以始终利用等效质点组计算部件 B在腔体内一点(X₀，Y₀，Z₀)的万有引力：

$a_{\mathrm{PG},x}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{X_{\mathrm{PG},b,i}-X_0}{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2{]}^{3/2}}$

$a_{\mathrm{PG},y}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{Y_{\mathrm{PG},b,i}-Y_0}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{3/2}}---\left(3\right)$

$a_{\mathrm{PG},z}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{Z_{\mathrm{PG},b,i}-Z_0}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{3/2}}$

同样，可以利用等效质点组计算部件B的万有引力梯度。

$T=\left(\begin{array}{ccc}{V}_{\mathrm{xx}}& V_{\mathrm{xy}}& V_{\mathrm{xz}}\\ V_{\mathrm{yx}}& V_{\mathrm{yy}}& V_{\mathrm{yz}}\\ V_{\mathrm{zx}}& V_{\mathrm{zy}}& V_{\mathrm{zz}}\end{array}\right)---\left(4\right)$

$V_{\mathrm{xx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{2{(X_{\mathrm{PG},b,i}-X_0)}^2-{(Y_{\mathrm{PG},b,i}-Y_0)}^2-{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{yy}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-{(X_{\mathrm{PG},b,i}-X_0)}^2+{2(Y_{\mathrm{PG},b,i}-Y_0)}^2-{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{zz}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-{(X_{\mathrm{PG},b,i}-X_0)}^2-{(Y_{\mathrm{PG},b,i}-Y_0)}^2+2{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{xy}}=V_{\mathrm{yx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(X_{\mathrm{PG},b,i}-X_0)(Y_{\mathrm{PG},b,i}-Y_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{xz}}=V_{\mathrm{zx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(X_{\mathrm{PG},b,i}-X_0)(Z_{\mathrm{PG},b,i}-Z_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{yz}}=V_{\mathrm{zy}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(Y_{\mathrm{PG},b,i}-Y_0)(Z_{\mathrm{PG},b,i}-Z_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}.$

以下举例说明计算的具体过程：

部件的质量特性包括质量、质心和某一特性坐标系下的转动惯量矩阵。

利用综合测试设备，测量得到：某一部件B的质量为M_o，在参考直角坐标系OXYZ中质 心位置为o(X_o，Y_o，Z_o)，在该部件主轴坐标系oyxz中的转动惯量矩阵为$\left(\begin{array}{ccc}{I}_{\mathrm{xx},o}& 0& 0\\ 0& I_{\mathrm{yy},o}& 0\\ 0& 0& I_{\mathrm{zz},o}\end{array}\right);$

设与部件B质量特性相同的等效质点组包含N_PG个质点，质点的质量记为 m_PG，i，i＝1，2，…N_PG，在该部件主轴坐标系oyxz中的坐标依次为(x_PG，i，y_PG，i，z_PG，i)，i＝1，2，…N_PG， 则等效质点组满足以下数学方程

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}=M_o$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}(y_{\mathrm{PG},i}^2+z_{\mathrm{PG},i}^2)=I_{\mathrm{xx},o},\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}(x_{\mathrm{PG},i}^2+z_{\mathrm{PG},i}^2)=I_{\mathrm{yy},o},\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{P,i}(x_{\mathrm{PG},i}^2+y_{\mathrm{PG},i}^2)=I_{\mathrm{zz},o}$

$\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{x}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0,\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}{y}_{\mathrm{PG},i}{z}_{\mathrm{PG},i}=0---\left(1\right)$

证明可知，满足上述方程的等效质点组与部件B，绕任意过质心的轴具有相同的转动惯 量。则根据引力位函数的二阶展开可知，等效质点组与部件B具有相同的万有引力计算结果。 因此，可用等效质点组代替部件B以二阶精度计算其外部的万有引力，且在实际计算中非常 方便，只需要计算等效质点组所包含的每个质点产生的引力，然后求和就得到等效质点组产 生的引力，也就是所代替的部件产生的引力。

求解方程组(1)就可以得到所需要的等效质点组。等效质点组中质点的坐标 (x_PG，i，y_PG，i，z_PG，i)，i＝1，2，…N_PG是在部件B的主轴坐标系中给出的。部件B的主轴坐标系到参 考直角坐标系的旋转矩阵为R_ob，则等效质点组在参考直角坐标系中的的坐标为

$(X_{\mathrm{PG},b,i},Y_{\mathrm{PG},b,i},Z_{\mathrm{PG},b,i})=(x_{\mathrm{PG},i},y_{\mathrm{PG},i},z_{\mathrm{PG},i})R_{\mathrm{ob}}^T+(X_o,Y_o,Z_o)---\left(2\right)$

记录等效质点组在参考直角坐标系中的以上坐标，则可以始终利用等效质点组计算部件 B在腔体内一点(X₀，Y₀，Z₀)的万有引力

$a_{\mathrm{PG},x}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{X_{\mathrm{PG},b,i}-X_0}{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2{]}^{3/2}}$

$a_{\mathrm{PG},y}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{Y_{\mathrm{PG},b,i}-Y_0}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{3/2}}---\left(3\right)$

$a_{\mathrm{PG},z}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{Z_{\mathrm{PG},b,i}-Z_0}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{3/2}}$

同样，可以利用等效质点组计算部件B的万有引力梯度。

$T=\left(\begin{array}{ccc}{V}_{\mathrm{xx}}& V_{\mathrm{xy}}& V_{\mathrm{xz}}\\ V_{\mathrm{yx}}& V_{\mathrm{yy}}& V_{\mathrm{yz}}\\ V_{\mathrm{zx}}& V_{\mathrm{zy}}& V_{\mathrm{zz}}\end{array}\right)---\left(4\right)$

$V_{\mathrm{xx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{2{(X_{\mathrm{PG},b,i}-X_0)}^2-{(Y_{\mathrm{PG},b,i}-Y_0)}^2-{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{yy}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-{(X_{\mathrm{PG},b,i}-X_0)}^2+{2(Y_{\mathrm{PG},b,i}-Y_0)}^2-{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{zz}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-{(X_{\mathrm{PG},b,i}-X_0)}^2-{(Y_{\mathrm{PG},b,i}-Y_0)}^2+2{(Z_{\mathrm{PG},b,i}-Z_0)}^2}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{xy}}=V_{\mathrm{yx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(X_{\mathrm{PG},b,i}-X_0)(Y_{\mathrm{PG},b,i}-Y_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{xz}}=V_{\mathrm{zx}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(X_{\mathrm{PG},b,i}-X_0)(Z_{\mathrm{PG},b,i}-Z_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

$V_{\mathrm{yz}}=V_{\mathrm{zy}}=G\underset{i=1}{\overset{N_{\mathrm{PG}}}{\Sigma }}{m}_{\mathrm{PG},i}\frac{-3(Y_{\mathrm{PG},b,i}-Y_0)(Z_{\mathrm{PG},b,i}-Z_0)}{{[{(X_{\mathrm{PG},b,i}-X_0)}^2+{(Y_{\mathrm{PG},b,i}-Y_0)}^2+{(Z_{\mathrm{PG},b,i}-Z_0)}^2]}^{5/2}}$

逐个得到各部件的万有引力及其梯度之后，对所有部件求和，就可以得到卫星在腔体内 的万有引力及其梯度。

综上可见，本发明在计算卫星上质量分布复杂部件的万有引力时，能够利用实际可测的 质量特性信息，以二阶精度实现万有引力的计算，且不受部件形状的影响，在工程上有较广 的适用范围。

以上通过具体的和优选的实施例详细的描述了本发明，但本领域技术人员应该明白，本 发明并不局限于以上所述实施例，凡在本发明的精神和原则之内，所作的任何修改、等同替 换等，均应包含在本发明的保护范围之内。