基于四面体质元划分的纯引力轨道万有引力干扰计算方法

摘要

本发明公开了一种涉及航天动力学和科学计算技术领域的基于四面体质元划分的纯引力轨道万有引力干扰计算方法，包括以下步骤：首先，选择四面体作为质元形状对航天器进行划分，并定义如下参数：尺度参数SR、长宽比参数AR、力几何因子FGP、梯度几何因子GGP，然后计算不同SR、AR下的FGP、GGP，得到四面体质元尺度、形状对计算精度的影响曲线，再根据所要求的计算精度及影响曲线确定四面体质元尺度，并根据该质元尺度对航天器模型采用四面体质元进行划分，最后对万有引力进行计算。本发明的方法能够按照要求的精度计算航天器对验证质量的万有引力。在给定精度的情况下，可以得到对应的四面体质元尺度要求，并依此要求划分质元，完成万有引力作用的计算，满足所要求的精度。

主权项

1.基于四面体质元划分的纯引力轨道万有引力干扰计算方法，包括以下步骤：第一步：选择四面体作为质元形状对航天器进行划分；第二步：定义如下参数：(1)表示四面体尺度与其距离外部一点的相对大小的尺度参数SR；(2)表示四面体形状变化的长宽比参数AR；(3)表示万有引力计算质点近似与精确表达式之间的相对误差的力几何因子FGP；(4)表示万有引力梯度计算质点近似与精确表达式之间的相对误差的梯度几何因子GGP；、其中，所述表示四面体尺度与其距离外部一点的相对大小的尺度参数SR定义如下：$\mathrm{SR}=\frac{r_{P/G}}{L_1}$式中，r_P/G是点P到四面体质心的距离，L₁为四面体最长边的长度；所述表示四面体形状变化的长宽比参数AR定义如下：$\mathrm{AR}=\frac{L_1}{L_2}$式中，L₁为四面体最长边的长度，L₂为四面体最短边的长度；表示万有引力计算质点近似与精确表达式之间的相对误差的力几何因子FGP定义如下：$\mathrm{FGP}=\frac{\left|F_{\mathrm{point}}^{P/C}\right|-\left|F^{P/C}\right|}{\left|F^{P/C}\right|}$式中，F^P/C为四面体C对位于其外部一点P的万有引力，为将四面体C近似为位于质心的质点时对点P的万有引力；所述的根据如下公式计算：$\left|F_{\mathrm{point}}^{P/C}\right|{\mathrm{Gm}}^P{\rho }^{C}V\frac{1}{r_{P/G}^2}$式中，G是万有引力常数、m^P为质点P的质量、ρ^C是四面体C的密度，V是四面体C的体积，点P到四面体质心的距离；所述的F^P/C通过如下方法计算：将四面体置于一个坐标系中，然后根据如下公式计算：$\left|F^{P/C}\right|=\sqrt{{\left(F_x^{P/C}\right)}^2+{\left(F_y^{P/C}\right)}^2+{\left(F_z^{P/C}\right)}^2}$四面体C对点P的万有引力在坐标轴x、y、z方向的分量；第三步：计算不同尺度参数SR、长宽比参数AR下的力几何因子FGP、梯度几何因子GGP，得到四面体质元尺度、形状对计算精度的影响曲线；所述第三步中，计算不同尺度参数SR、长宽比参数AR下的力几何因子FGP通过如下方法实现：将正四面体按照如下方法置于一个坐标系中：正四面体的一个顶点位于原点，该顶点所对的底面平行于x-y平面且位于上方，并且底面内的一条边平行于x轴；将所述正四面体做如下变形获得一系列不同长宽比参数AR的四面体：将正四面体拉伸，且与x-y平面平行的底面边长保持不变作为所述L₂，而其它三个边被拉长同等长度作为所述L₁，拉伸不同的长度即可得不同的长宽比参数AR，计算不同长宽比参数AR、尺度参数SR下的力几何因子FGP，所述的计算公式如下：$\frac{F_x^{P/C}}{{\mathrm{Gm}}^P}={\rho }^C{\int }_0^h{\int }_0^{\frac{\zeta }{h}a}{\int }_{\frac{\zeta }{h}b}^{\frac{b-c}{a}\xi +\frac{\zeta }{h}c}\frac{\xi -x}{{[{(\xi -x)}^2+{(\eta -y)}^2+{(\zeta -z)}^2]}^{3/2}}\mathrm{d\eta d\xi d\zeta }$$+{\rho }^C{\int }_0^h{\int }_{-\frac{\zeta }{h}a}^0{\int }_{\frac{\zeta }{h}b}^{-\frac{b-c}{a}\xi +\frac{\zeta }{h}c}\frac{\xi -x}{{[{(\xi -x)}^2+{(\eta -y)}^2+{(\zeta -z)}^2]}^{3/2}}\mathrm{d\eta d\xi d\zeta }$$\frac{F_y^{P/C}}{{\mathrm{Gm}}^P}={\rho }^C{\int }_0^h{\int }_0^{\frac{\zeta }{h}a}{\int }_{\frac{\zeta }{h}b}^{\frac{b-c}{a}\xi +\frac{\zeta }{h}c}\frac{\eta -y}{{[{(\xi -x)}^2+{(\eta -y)}^2+{(\zeta -z)}^2]}^{3/2}}\mathrm{d\eta d\xi d\zeta }$$+{\rho }^C{\int }_0^h{\int }_{-\frac{\zeta }{h}a}^0{\int }_{\frac{\zeta }{h}b}^{-\frac{b-c}{a}\xi +\frac{\zeta }{h}c}\frac{\eta -y}{{[{(\xi -x)}^2+{(\eta -y)}^2+{(\zeta -z)}^2]}^{3/2}}\mathrm{d\eta d\xi d\zeta }$$\frac{F_z^{P/C}}{{\mathrm{Gm}}^P}={\rho }^C{\int }_0^h{\int }_0^{\frac{\zeta }{h}a}{\int }_{\frac{\zeta }{h}b}^{\frac{b-c}{a}\xi +\frac{\zeta }{h}c}\frac{\zeta -z}{{[{(\xi -x)}^2+{(\eta -y)}^2+{(\zeta -z)}^2]}^{3/2}}\mathrm{d\eta d\xi d\zeta }$$+{\rho }^C{\int }_0^h{\int }_{-\frac{\zeta }{h}a}^0{\int }_{\frac{\zeta }{h}b}^{-\frac{b-c}{a}\xi +\frac{\zeta }{h}c}\frac{\zeta -z}{{[{(\xi -x)}^2+{(\eta -y)}^2+{(\zeta -z)}^2]}^{3/2}}\mathrm{d\eta d\xi d\zeta }$其中G是万有引力常数、m^P为质点P的质量、ρ^C是四面体C的密度，h为四面体高度，参数a为四面体C中与x-y平面平行的底面边长的一半，参数b为四面体C在x-y平面的投影三角形中平行于x轴的底边与y轴的交点到原点的距离，参数c为四面体C在x-y平面的投影三角形中平行于x轴的底边对应的顶点到原点的距离，(x、y、z)为质点P的坐标，(ξ、η、ζ)为四面体C内部一点的直角坐标；表示万有引力梯度计算质点近似与精确表达式之间的相对误差的梯度几何因子GGP定义如下：$\mathrm{GGP}=\frac{\left|{\lambda }_{\mathrm{point},p}\right|-\left|{\lambda }_p\right|}{\left|{\lambda }_p\right|}$式中，λ_p为四面体C对位于其外部一点P的引力梯度张量矩阵T^P/C的主特征值，λ_point，p为将四面体C近似为位于质心的质点时对点P的引力梯度张量矩阵的主特征值；所述第三步中，计算不同尺度参数SR、长宽比参数AR下的梯度几何因子GGP通过如下方法实现：将正四面体按照如下方法置于一个坐标系中：正四面体的一个顶点位于原点，该顶点所对的底面平行于x-y平面且位于上方，并且底面内的一条边平行于x轴；将所述正四面体做如下变形获得一系列不同长宽比参数AR的四面体：将正四面体拉伸，且与x-y平面平行的底面边长保持不变作为所述L₂，而其它三个边被拉长同等长度作为所述L₁，拉伸不同的长度即可得不同的长宽比参数AR，计算不同长宽比参数AR、尺度参数SR下的梯度几何因子GGP，所述的引力梯度张量矩阵T^P/C通过如下公式计算：$T^{P/C}=\left(\begin{array}{ccc}\frac{{\partial F}_x^{P/C}}{\partial x}& \frac{{\partial F}_x^{P/C}}{\partial y}& \frac{{\partial F}_x^{P/C}}{\partial z}\\ \frac{{\partial F}_y^{P/C}}{\partial x}& \frac{{\partial F}_y^{P/C}}{\partial y}& \frac{{\partial F}_y^{P/C}}{\partial z}\\ \frac{{\partial F}_z^{P/C}}{\partial x}& \frac{{\partial F}_z^{P/C}}{\partial y}& \frac{{\partial F}_z^{P/C}}{\partial z}\end{array}\right)$式中，四面体C对点P的万有引力在坐标轴x、y、z方向的分量；所述引力梯度张量矩阵通过如下公式计算：$T_{\mathrm{point}}^{P/C}=\left(\begin{array}{ccc}\frac{{\partial F}_{\mathrm{point},x}^{P/C}}{\partial x}& \frac{{\partial F}_{\mathrm{point},x}^{P/C}}{\partial y}& \frac{{\partial F}_{\mathrm{point},x}^{P/C}}{\partial z}\\ \frac{{\partial F}_{\mathrm{point},y}^{P/C}}{\partial x}& \frac{{\partial F}_{\mathrm{point},y}^{P/C}}{\partial y}& \frac{{\partial F}_{\mathrm{point},y}^{P/C}}{\partial z}\\ \frac{{\partial F}_{\mathrm{point},z}^{P/C}}{\partial x}& \frac{{\partial F}_{\mathrm{point},z}^{P/C}}{\partial y}& \frac{{\partial F}_{\mathrm{point},z}^{P/C}}{\partial z}\end{array}\right)$式中，$\frac{{\partial F}_{\mathrm{point},x}^{P/C}}{\partial x}={\mathrm{Gm}}^P{\rho }^{C}V\frac{{2x}^2-y^2-{(z_Q-z)}^2}{r_{P/Q}^5},$$\frac{{\partial F}_{\mathrm{point},y}^{P/C}}{\partial y}={\mathrm{Gm}}^P{\rho }^{C}V\frac{{-x}^2+2y^2-{(z_Q-z)}^2}{r_{P/Q}^5}$$\frac{{\partial F}_{\mathrm{point},z}^{P/C}}{\partial z}={\mathrm{Gm}}^P{\rho }^{C}V\frac{{-x}^2-y^2+2{(z_Q-z)}^2}{r_{P/Q}^5},$$\frac{{\partial F}_{\mathrm{point},x}^{P/C}}{\partial y}=\frac{{\partial F}_{\mathrm{point},y}^{P/C}}{\partial x}={\mathrm{Gm}}^P{\rho }^{C}V\frac{-3xy}{r_{P/Q}^5}$$\frac{{\partial F}_{\mathrm{point},x}^{P/C}}{\partial z}=\frac{{\partial F}_{\mathrm{point},z}^{P/C}}{\partial x}={\mathrm{Gm}}^P{\rho }^{C}V\frac{3x(z_Q-z)}{r_{P/Q}^5},$$\frac{{\partial F}_{\mathrm{point},y}^{P/C}}{\partial z}=\frac{{\partial F}_{\mathrm{point},z}^{P/C}}{\partial y}={\mathrm{Gm}}^P{\rho }^{C}V\frac{3y(z_Q-z)}{r_{P/Q}^5}$其中G是万有引力常数、m^P为质点P的质量、ρ^C是四面体C的密度，V是四面体C的体积，点P到四面体质心的距离，(x、y、z)为质点P的坐标，z_Q为四面体C在的质心在z轴上的坐标；第四步：根据所要求的万有引力干扰计算精度及所述的影响曲线确定相应的四面体质元尺度；第五步：根据第四步所确定的质元尺度对航天器模型采用四面体质元进行划分；第六步：将每一个四面体质元近似为点质量计算其对验证质量的万有引力及梯度；第七步：将所有四面体质元的计算结果求和，得到航天器对验证质量的万有引力干扰。