CN114662294A - Picometer-level satellite orbit simulation method and system for deep space exploration - Google Patents

Picometer-level satellite orbit simulation method and system for deep space exploration Download PDF

Info

Publication number
CN114662294A
CN114662294A CN202210230755.5A CN202210230755A CN114662294A CN 114662294 A CN114662294 A CN 114662294A CN 202210230755 A CN202210230755 A CN 202210230755A CN 114662294 A CN114662294 A CN 114662294A
Authority
CN
China
Prior art keywords
satellite
solver
simulation
orbit
satellite orbit
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.)
Granted
Application number
CN202210230755.5A
Other languages
Chinese (zh)
Other versions
CN114662294B (en
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.)
National Space Science Center of CAS
Original Assignee
National Space Science Center of CAS
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 National Space Science Center of CAS filed Critical National Space Science Center of CAS
Priority to CN202210230755.5A priority Critical patent/CN114662294B/en
Publication of CN114662294A publication Critical patent/CN114662294A/en
Application granted granted Critical
Publication of CN114662294B publication Critical patent/CN114662294B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

The invention relates to the field of satellite orbit simulation, in particular to a pico-meter satellite orbit simulation method and a pico-meter satellite orbit simulation system for deep space exploration, which comprise the following steps: step 1) adding gravitational fields of a plurality of celestial bodies, and establishing a satellite formation dynamic model; step 2) setting a required effective digit, satellite orbit simulation starting time, an optimization parameter for satellite orbit initialization, a simulation step length and a parameter of the solver; and calculating the dynamic model of the satellite formation through the solver and a basic operation library developed based on an MPFR library to obtain a position and speed matrix of the satellite at the moment corresponding to the simulation step length. The invention can carry out floating point operation with any significant digit, the simulation precision is not limited by double-precision 16-bit significant digits any more, and the precision is superior to 0.03 picometer magnitude.

Description

Picometer-level satellite orbit simulation method and system for deep space exploration
Technical Field
The invention relates to the field of satellite orbit simulation, in particular to a pico-meter satellite orbit simulation method and system for deep space exploration.
Background
The space laser interference gravitational wave antenna is an ultra-precise measuring system, the Taiji plan proposed by Chinese academy of sciences is to adopt a three-star formation group running on the solar centre orbit to make space gravitational wave detection, and the distance between the satellite formation group and the solar centre is about 1 astronomical unit (1AU is equal to 1.496 multiplied by 10)11m), each satellite contains two test qualities. The method comprises using the test mass as gravitational wave sensor, and measuring the picometer-level (1 × 10) between two free test masses (several million kilometers) caused by gravitational wave by using high-precision laser interference ranging system-12m) distance changes, and thereby inverting gravitational wave signals.
Due to extremely high detection precision requirements and task complexity, the ground simulation system for establishing the detection formation is necessary, and the on-orbit running state and the like are simulated, so that support is provided for optimizing a task design scheme and assisting key technology attack. However, for "taiji", it is necessary to simulate the change of picometer magnitude caused by gravitational wave in the orbit scale of 1AU satellite, and the accuracy requirement of long baseline measurement is high, which is far beyond the range of double accuracy. On one hand, the calculation precision brought by computer programming language is not high enough, and the error is too large due to the precision problem; on the other hand, when the orbit is solved by numerical integration, the accuracy of the ordinary differential equation solver is far insufficient to meet the requirement due to the insufficient double accuracy.
In recent years, aiming at high-precision rail simulation, in low-low satellite tracking (LL-SST) tasks, the precision of rail simulation is improved by adopting a mode of mixing double precision and quadruple precision, but the precision only reaches the nanometer level, and has a great difference with the picometer-level precision required by space gravitational wave detection task 'Taiji' plan. An extended high performance satellite dynamics simulator (XHPS) developed based On Grace Follow-On tasks, uses a Runge-Kutta and a multi-step integrator to realize different numerical integration schemes, but the precision level is limited by the selected integrator, so more bits are needed. In the orbit simulation of the geocentric task 'Tianqin' of space gravitational wave detection, 34 effective numbers are expanded through a Matlab quadruple precision algorithm, and a four-precision orbit numerical simulation program of the Tianqin is developed. However, the scale spanned by the above problems is not as wide as the detection task of the space gravitational wave of the centrobaric orbit, and the wide range is not considered in the cross-scale calculation.
MPFR C + + is an open-source header file library of a multi-precision efficient floating-point computing C language, floating-point computation (limited only by available memory) with any significant digit can be performed by including an mp. sin, sqrt, pow, log, etc. However, at present, no method for performing orbit numerical simulation by using an MPFR library exists.
In summary, the technical precision of the current orbit simulation is poor, the pico-meter precision simulation cannot be achieved, and no effective solution exists for the pico-meter precision simulation problem of the satellite orbit under the space scale of 1AU in the deep space exploration task such as the space gravitational wave exploration task of the solar central orbit.
Disclosure of Invention
The invention aims to solve the problem that the satellite orbit numerical simulation precision is insufficient in the current deep space exploration tasks, such as the heliocentric space gravitational wave exploration task, and provides a high-precision orbit numerical simulation method and system suitable for the heliocentric orbit spatial gravitational wave exploration task, which are used for simulating the orbit of a satellite under the scale of 1AU in picometer order and verifying the correctness and the effectiveness of the method.
In order to solve the problems, the invention designs a pico-meter satellite orbit simulation method and system for deep space exploration based on an MPFR library MPreal.h header file, solves the problem of insufficient orbit simulation precision of exploration formation through a high-precision ordinary differential equation solver capable of realizing floating point operation with any significant digit, and realizes the requirement that the orbit simulation precision reaches the pico-meter level or even higher under the spatial scale of 1 AU. According to the pico-meter magnitude satellite orbit simulation method for deep space exploration, provided by the invention, a satellite orbit is simulated through a solver and a basic operation library so as to obtain a position and speed matrix of a satellite; the method is characterized by comprising the following steps:
step 1) adding gravitational fields of a plurality of celestial bodies, and establishing a satellite formation dynamic model;
step 2) setting a required effective digit, satellite orbit simulation starting time, an optimization parameter for satellite orbit initialization, a simulation step length and a parameter of the solver; and solving the satellite formation kinetic model through the solver and a basic operation library developed based on an MPFR library to obtain a position and speed matrix of the satellite at the moment corresponding to the simulation step length.
As an improvement of the above method, the basic operation library is developed in the underlying framework by using an MPFR library, wherein the MPFR library includes: an mp. The basic operation library developed based on the MPFR library comprises: the device comprises a matrix operation module, a track element calculation module, a character string type conversion module, a time conversion module, a coordinate system conversion module, a data input module and a data output module.
As an improvement of the method, step 2-1) sets the required significant digit, the start time of satellite orbit simulation, the optimization parameters for satellite orbit initialization, the simulation step length and the parameters of the solver;
step 2-2) converting the set starting time of the satellite orbit simulation from Beijing time to julian day time through the time conversion module;
step 2-3) reading the optimization parameters for satellite orbit initialization through the data input module;
step 2-4), acquiring an initialized satellite position and speed matrix corresponding to the starting time of the satellite orbit based on the optimized parameters through the matrix operation module;
step 2-5), converting the initialized satellite position and speed matrix into a coordinate system required by the satellite formation kinetic model through the coordinate system conversion module;
step 2-6), performing ordinary differential equation solution on the satellite formation kinetic model by using the set parameters of the solver, the simulation step length, the effective digit and the initialized satellite position and speed matrix of the required coordinate system through the solver to obtain the satellite position and speed with the required effective digit at the moment corresponding to the simulation step length under the coordinate system required by the satellite formation kinetic model; and converting the satellite position and the satellite velocity into an inertial system through the coordinate system conversion module so as to obtain the satellite position and the satellite velocity with the significant digit at the moment corresponding to the simulation step length under the inertial system.
As an improvement of the above method, the solver is of the type comprising: 8-order variable step length Longge Kutta solver and 8-order variable step length DOP853 solver; the solver is further developed by using the basic operation library developed based on the MPFR library so as to perform floating point operation with any significant digit.
As an improvement of the above method, the solver adopts an 8-order variable step length DOP853 solver and adopts a 7-order dense output, wherein the parameters of the solver are set as: tolerance of error is 1 x 10-NWherein N is any positive integer.
As an improvement of the above method, the solver converts the solver coefficients in the form of a character string and the optimized parameters for satellite initialization into an mpreal data type using a stringtorereal function in the character string type conversion module, so that the solver coefficients and the optimized parameters for satellite initialization are converted into constants corresponding to the significands.
As an improvement of the above method, the satellite formation dynamics model is established by using the following formula:
Figure BDA0003538264100000031
wherein i is the serial number of the satellite,
Figure BDA0003538264100000032
an acceleration matrix of the ith satellite in the inertial system;
Figure BDA0003538264100000033
is a matrix of positions of the ith satellite in the inertial system,
Figure BDA0003538264100000034
is a position matrix of the kth celestial body in the inertial system; mkMass of kth celestial body; g is a universal gravitation constant; wherein the ordering of the plurality of celestial bodies is: sun, water star, venus venosus, earth and moon, mars, saturn, asterias, plundersiatus; when the celestial body is the earth and the moon, the mass of the celestial body is the sum of the earth mass and the moon mass.
As an improvement of the above method, the gravitational field of the plurality of celestial bodies comprises: the sun's central gravity, the water star gravity, the golden star gravity, the sum of the earth and moon gravities, the mars gravity, the wooden star gravity, the earth star gravity, the heaven star gravity, the sea king star gravity and the pluto gravity.
In order to achieve another object of the present invention, the present invention provides a pico-meter satellite orbit simulation system for deep space exploration, which is used for executing the pico-meter satellite orbit simulation method for deep space exploration, and includes: a solver and a basic computation library, wherein,
the basic operation library is developed by utilizing an MPFR library in a bottom layer framework, and comprises the following steps: the device comprises a matrix operation module, a track element calculation module, a character string type conversion module, a time conversion module, a coordinate system conversion module, a data input module and a data output module; wherein the content of the first and second substances,
the time conversion module is used for converting the initial time of the satellite orbit simulation from Beijing time to julian daily time;
the data input module is used for reading the optimization parameters for satellite orbit initialization;
the matrix operation module is used for acquiring an initialized satellite orbit position and speed matrix corresponding to the starting time of the satellite orbit based on the optimization parameters;
the coordinate system conversion module is used for converting the initialized satellite orbit position and speed matrix into a coordinate system required by the satellite formation kinetic model;
the solver is used for solving an ordinary differential equation of the satellite formation kinetic model by using the set parameters of the solver, the simulation step length, the effective digit and the initialized satellite orbit position and speed matrix in the required coordinate system so as to obtain the satellite orbit position and speed with the required effective digit of the satellite formation kinetic model in the required coordinate system at the moment corresponding to the simulation step length; and the satellite orbit position and the satellite orbit speed with the significant digit at the moment corresponding to the simulation step length under the inertial system are obtained by converting the satellite orbit position and the satellite orbit speed into the inertial system through the coordinate system conversion module.
As an improvement of the above system, the solver is of the type comprising: 8-order variable step length Longge Kutta solver and 8-order variable step length DOP853 solver; the solver is further developed by utilizing the development of the basic operation library based on the MPFR library to carry out floating point operation with any significant digit.
The pico-meter-level satellite orbit simulation method for deep space exploration provided by the invention has the following advantages:
1. the invention is about 1.496X 101110 caused by simulated gravitational wave under satellite orbit scale of m (1AU)-12The invention develops a basic operation library which can realize floating point operation with any significant digit and comprises matrix operation, data input, time conversion, coordinate system conversion, orbit element calculation, character string type conversion, model construction, solver and data output module based on an open source MPFR library, thereby realizing high-precision orbit numerical simulation calculation.
2. The basic operation library developed based on the MPFR library ensures that the bottom layer framework and the solver have wide applicability and can be also applied toThe orbit numerical value of other deep space exploration tasks is solved, and the required significant digit and solver parameters can be selected according to the task requirements so as to achieve different accuracies; through verification, when 34-bit significant digits are set, the integration time is 6 years, and the simulation step length is 86400s (the selection of the integration time period and the simulation step length does not cause great influence on the simulation precision level), the error tolerance of the 8-order variable-step Runge Kutta solver is set to be 1 x 10 for the simulation of the position and the speed of the satellite running in the orbit of the Sun's mind-27The maximum error of simulation for the satellite position is 8 x 10-18The maximum error of simulation for the satellite position is 1.5 multiplied by 10 at the precision level within m-24A level of accuracy within m/s; the error tolerance of the 8-order variable step DOP853 solver is set to be 1 multiplied by 10-29The maximum error of simulation for the satellite position is 2 x 10-17The maximum error of simulation for the satellite position is 4 multiplied by 10 at the precision level within m-24A level of accuracy within m/s. Furthermore, a solver developed by using a basic operation library developed based on an MPFR library has the characteristic of being suitable for floating point operation with any significant digit, and the solving precision is not limited by double-precision 16-bit significant digits any more, so that the precision is greatly improved.
3. The invention simulates the in-orbit running state of the Tai Chi planned satellite formation and verifies that the method can simulate the position of the satellite under the space scale of 1AU with the precision superior to 0.03 picometer by utilizing the problem of disomic solution, and the precision of the orbit can meet the requirement of subsequent processing, thereby providing support for the design of a detection task.
Drawings
FIG. 1 is a flowchart of a pico-meter satellite orbit simulation method for deep space exploration according to an embodiment of the present invention;
fig. 2(a) is a schematic diagram of a satellite velocity three-component difference between a numerical solution and an analytic solution when an 8-order fixed-step length longstotta solver is adopted in the embodiment of the present invention;
fig. 2(b) is a schematic diagram of a satellite position three-component difference between a numerical solution and an analytic solution when an 8-order fixed-step length longstotta solver is adopted in the embodiment of the present invention;
FIG. 3(a) is a schematic diagram of a satellite velocity three-component difference between a numerical solution and an analytic solution when an 8-step variable-step Runge Kutta solver is adopted in the embodiment of the present invention;
fig. 3(b) is a schematic diagram of a satellite position three-component difference between a numerical solution and an analytic solution when an 8-step variable-step longge stotta solver is adopted in the embodiment of the present invention;
FIG. 4(a) is a schematic diagram of the three-component difference of satellite velocity between a numerical solution and an analytic solution when an 8-step DOP853 solver is adopted in the embodiment of the present invention;
fig. 4(b) is a schematic diagram of the satellite position three-component difference between the numerical solution and the analytic solution when the 8-step DOP853 solver is adopted in the embodiment of the present invention.
Detailed Description
The invention designs a satellite orbit high-precision solving method suitable for a space gravitational wave detection task of a centrobaric orbit, and the correctness and precision of the solving of the method are evaluated through a two-body problem. The method mainly comprises four parts: the method comprises the following steps of (1) bottom layer framework (C + +), solver development, force model construction and track solving, wherein the method comprises the following steps:
step 1) adding gravitational fields of a plurality of celestial bodies, and establishing a satellite formation dynamic model;
step 2) setting a required effective digit, satellite orbit simulation starting time, an optimization parameter for satellite orbit initialization, a simulation step length and a parameter of the solver; and solving the satellite formation kinetic model through the solver and a basic operation library developed based on an MPFR library to obtain a position and speed matrix of the satellite at a moment corresponding to the simulation step length in the set starting time, wherein in the embodiment, because the position and the speed at each moment are both 1 x 3 matrixes, the position and the speed of the satellite at the moment corresponding to the simulation step length in the simulation starting time are respectively n x 3 matrixes.
Specifically, the technical solution and flow of the present invention are given below with reference to the accompanying drawings, and in addition, it is necessary to verify the above-mentioned high-precision orbit numerical solution method by a two-body problem with an analytic solution, so as to evaluate the correctness and precision of the orbit solution method.
The technical scheme provided by the invention is further illustrated by combining the following embodiments.
The technical scheme of the invention is shown in figure 1, and the functions of various parts are described below.
The underlying framework (C + +) first includes the mpreal. h header files in the open-source MPFR library, including simple arithmetic operations like "+" and "/", and mathematical functions: sin, sqrt, pow, log, etc. Based on the database, a basic operation library which can realize high-precision matrix operation, data input, time conversion, coordinate system conversion, orbit element calculation, character string type conversion and data output of floating point operation with any significant digit is developed for realizing high-precision orbit numerical value solution. The matrix operation module in the basic operation library is mainly used in the steps of development of a solver, construction of a force model, subsequent data output, coordinate system conversion, calculation of track elements and the like which need matrixes. The orbit element calculation module, the time conversion module and the coordinate system conversion module in the basic operation library are mainly used for the construction of a force model and the subsequent orbit solving process; the data input/output of the basic operation library is mainly used for parameter reading in the orbit solving process and outputting and storing data such as satellite position and speed obtained through calculation; in the process of constructing the force model, conversion between Beijing time and the Breviz day time is involved, conversion between an inertial system and other coordinate systems (required in satellite orbit solution) is involved in the process of orbit calculation, and operations such as matrix base and transposition in a basic operation base, matrix multiplication and division and the like are required in the process of a column motion equation, and the orbit element calculation module is required at multiple places, for example, in the process of initializing a position and speed matrix of a satellite, constructing the force model and subsequently solving the orbit. Orbit element calculations are common in the satellite field and will not be described in detail here. For a basic operation library, one-dimensional, two-dimensional and three-dimensional matrix libraries and related matrix operations which can carry out floating point operations with any significant digit are autonomously realized by using the mpreal data type; floating point data and data files with any significant digit can be read and output; the time conversion module, the coordinate system conversion module and the orbit element calculation module can also realize floating point operation with any significant digit so as to keep the uniformity of precision.
Secondly, in the development of a solver capable of realizing floating point operation of any significant digit, the solver is realized by applying an independently realized basic operation library based on the bottom framework, and the solution of any significant digit of an ordinary differential equation can be realized. In the realization process of the solver, the stringtoreall function of the character string type conversion module which is realized autonomously sets the coefficients of the solver to be constants with corresponding significant digits, so as to avoid precision loss caused by insufficient effective digits of double precision to accurately represent the coefficients (the stringtoreall function obtains accurate representation of any significant digit by converting the coefficients which are input in the form of character strings into the type of temporal data). And the solver step is unchanged, and the data type is converted into an mpreal data type and one-dimensional, two-dimensional and three-dimensional matrix types of a matrix library in the basic operation library.
And thirdly, in the construction of the force model, modeling a dynamic model of the satellite formation aiming at a space gravitational wave detection task 'Taiji' plan. In the force model, only the main celestial gravity field is added, including the sun's central gravity, the eight planets (earth and moon considered together), and the pluto. The force model is calculated as follows
Figure BDA0003538264100000071
Wherein
Figure BDA0003538264100000072
Indicating the position of the ith satellite,
Figure BDA0003538264100000073
the positions of the sun and nine planets are shown in the inertial system. MkThe masses of the sun and the nine planets are shown, i is the serial number of the satellite, and in this embodiment, the number of the satellites used is 3.
And fourthly, finally, carrying out orbit solution, selecting the significant digit, setting the initial time and the solution step length of the orbit solution, initializing the position and the speed of the orbit, setting the parameters of the solver in the second step according to the needed orbit precision, and carrying out high-precision solution on the orbit of the Taiji. And when the orbit solution is finally carried out, the selection of the significant digit and the solver is required to be carried out according to the task requirement.
For a space gravitational wave detection task 'tai chi' plan, simulation of pico-meter-level precision of a satellite orbit needs to be realized under the space scale of 1AU, and a solver and parameter setting thereof which pass through the two-body problem verification precision and meet the requirements are utilized to carry out high-precision solution on the position and the speed of a satellite formation of the 'tai chi' task, so that the precision of data can meet the requirement of subsequent processing. When 34 significant digits are adopted, 8-order step length DOP853 solver is selected, and the error tolerance parameter is set to be 1 multiplied by 10-26The simulation step size is set to 86400s, the integration time is six years, namely from 1 month and 1 day of 2030 to 1 month and 1 day of 2036, and the simulation of the satellite position can reach 3 multiplied by 10-14m is the precision level of 0.03 picometers; when 34-bit significant digits are adopted, an 8-order variable-step Runge Kutta solver is selected, and the error tolerance is set to be 1 multiplied by 10-23The simulated maximum error of the satellite position is 4 multiplied by 10 for the same integration time-14The accuracy level within m and the simulation maximum error of the satellite position are 6 multiplied by 10-21The precision level within m/s also reaches the precision level of 0.03 picometer, and the task requirements can be met.
The verification of the accuracy of the orbit by the two-body problem does not fall within the scope of the technical solution of the present invention, but the accuracy of the method of the present invention needs to be verified by using the same, so it is necessary to explain. A force model of the orbit solution method is constructed by using a bottom layer framework of the first step and a solver of the second step and through an undisturbed two-body problem with an analytic solution, and the accuracy and the correctness of the orbit solution method are evaluated. Here we use the "earth-sun" dyads problem for accuracy verification because it has the same spatial scale as the spatial gravitational wave detection task "tai chi" and its analytical solution can be easily calculated. The satellite states can be represented by kepler orbital elements (a, e, i, Ω, ω, M), and specific values are shown in table 1, which are semi-major axis, eccentricity, orbital inclination, ascension at ascending intersection, argument of near point and argument of near point, respectively. The analytic solution to the two-body problem, the position vector R of the satellite in the inertial reference frame is given by
Figure BDA0003538264100000081
Wherein R is1And R3Euler rotation matrices around the x-axis and z-axis respectively; e represents the off-proximal angle, which can be easily calculated by the kepler equation using the eccentricity and the flat-proximal angle:
E-esinE=M;
the mean-near-point angles at different times t can be solved by the following formula
M=M0+n(t-t0);
M0Is the mean anomaly angle at the initial time, n is the track angular velocity, where,
Figure BDA0003538264100000082
t0is the initial time, t is the desired time.
The velocity vector V of the satellite is given by:
Figure BDA0003538264100000083
in the precision verification module, 34-bit significant digits are adopted, and three solvers are used for respectively solving numerical orbits, namely an 8-order fixed-step-length Rungestota solver, an 8-order variable-step-length Rungestota solver and an 8-order variable-step-length DOP853 solver are used for solving the numerical orbits. The analytical solution is likewise solved using 34-bit significands. For the numerical solution and the solution of the analytic solution, the same initial conditions are set, and the initial keplerian orbit elements thereof are shown in table 1.
TABLE 1 initial Kepler orbital elements at earth 2030.01.01 null
Track element Value of orbital element of the earth
a(m) 1.495978852117977450000000000000000×1011
e 1.669880800000000000000000000000000×10-2
i(rad) -6.739879794784769382886985850980193×10-5
Ω(rad) 6.060138148349757504428920819197350
ω(rad) 2.021557407619410035649390804621808
M(rad) 6.234822626305905303072886103194597
For an 8-order fixed-step Runge Kutta solver, selecting an integral step length as one day, namely 86400 seconds; the integral step length of the 8-order variable-step Longge Kuta solver is automatically changed, and the error tolerance is set to be 1 multiplied by 10-23(ii) a For the 8-order variable step length DOP853 solver, 7-order dense output and fixed output interval of 86400s are adopted, and the error tolerance is set to be 1 × 10-26. In the field of the artIt should be understood by the person that the parameter settings of the solver are selected according to the requirements, and the error tolerance setting herein is set for the tai chi task to meet the requirements in the embodiment; the error tolerance can be set as desired at different significands, typically 1 × 10-NAnd N is a positive integer.
The simulation step length has almost no influence on the accuracy level of the result, can be set to be any positive integer seconds according to the task requirement, and is only limited by the size of the memory. Here, the setting of 86400s is only the setting adopted in one embodiment, and the setting of other positive integers has no influence. For two variable step length algorithms, the setting of the error tolerance is closely related to the precision of the result, and only the result of the finally adopted setting of the error tolerance is shown in the technical effect part of the application. The integration time is six years, from 1 month 1 day 2030 to 1 month 1 day 2036. The integral orbit (namely, a numerical solution) of the earth-sun two-body problem is subtracted from the analytic orbit (namely, the analytic solution) to obtain a difference value of the two, so as to evaluate the accuracy level of the orbit numerical integration.
It should be noted that the significand and the solver setting used in the above setting are both selected for meeting the task requirements of the "tai chi" plan, and we can select a proper significand and solver setting according to different task requirements, thereby achieving different precision levels.
As shown in fig. 2(a) -2 (b), when 34 significant digits, integration time of 6 years and simulation step length of 86400s are set, for the simulation of the position and speed of the satellite operating in the solar orbit, the simulation maximum error of the 8-step fixed-step Runge Tower solver on the satellite position is 2 × 10-3A level of accuracy within m, a simulated maximum error for satellite position of 4 x 10-10A level of accuracy within m/s; although the 8-order fixed-step Longge library tower solver does not reach the accuracy level of picometer magnitude, the accuracy is improved compared with the existing satellite position and speed simulation method.
As shown in fig. 3(a) -3 (b), the position and velocity for a satellite operating in the heliocentric orbitDegree simulation, the error tolerance of the 8-order variable-step length Longge Kuta solver is set to be 1 multiplied by 10-23The simulated maximum error for the satellite position is 4 x 10-14Accuracy level within m, simulated maximum error for satellite position of 6 x 10-21A level of accuracy within m/s; as shown in FIGS. 4(a) -4 (b), the error tolerance of the DOP853 solver with 8 th step is set to 1 × 10-26The simulated maximum error for the satellite position is 3 × 10-14The maximum error of simulation for the satellite position is 5 multiplied by 10 at the precision level within m- 21A level of accuracy within m/s. The time performance difference of the two algorithms is large, when the same precision level is realized, the time required by the 8-order variable-step Runge Kutta solver is about 12 times that of the 8-order variable-step DOP853 solver, however, both the 8-order variable-step Runge Kutta solver and the 8-order variable-step DOP853 solver can meet the precision level of 0.03 picometer, and when the two algorithms are actually applied, a proper solver can be selected according to task requirements.
From the above description of the invention it can be seen that:
1. the method can select different significands and solver parameters according to task requirements, achieve different precision levels, and meet different task requirements. For example, according to the requirement of the "tai chi" planning task, in this embodiment, when we use the quadruple precision, namely 34-bit significant digits, the evaluation through the dyadic problem shows that, as shown in fig. 4(a) -4 (b), the error tolerance of the use of the 8-step DOP853 solver is set to be 1 × 10-26The simulation of the satellite position can realize that the maximum error is 3 multiplied by 10 within 6 years-14Within m, namely the precision level of 0.03 picometer, the maximum error of 5 multiplied by 10 within 6 years can be realized for the speed simulation-21A level of accuracy within m/s. For the simulation of the position and the speed of the satellite running in the orbit of the Sun's heart, the error tolerance of the 8-step variable-step Runge Kutta solver is set to be 1 multiplied by 10-23The simulated maximum error for the satellite position is 4 x 10-14Accuracy level within m, simulated maximum error for satellite position of 6 x 10-21The precision level within m/s also reaches the precision level of 0.03 picometer。
2. The high-precision orbit numerical solving method can be suitable for other deep space exploration tasks with higher requirements on orbit precision.
3. The high-precision simulation of the satellite running state is carried out through the ground simulation verification platform, the complex tasks can be better verified and optimized, potential problems and risks are analyzed, and a foundation is laid for the subsequent launching of scientific satellites.
Finally, it should be noted that the above embodiments are only used for illustrating the technical solutions of the present invention and are not limited. Although the present invention has been described in detail with reference to the embodiments, it will be understood by those skilled in the art that various changes may be made and equivalents may be substituted without departing from the spirit and scope of the invention as defined in the appended claims.

Claims (10)

1. A pico-meter magnitude satellite orbit simulation method for deep space exploration is characterized in that a satellite orbit is simulated through a solver and a basic operation library to obtain a position and speed matrix of a satellite; the method is characterized by comprising the following steps:
step 1) adding gravitational fields of a plurality of celestial bodies, and establishing a satellite formation dynamic model;
step 2) setting a required effective digit, satellite orbit simulation starting time, an optimization parameter for satellite orbit initialization, a simulation step length and a parameter of the solver; and solving the satellite formation kinetic model through the solver and a basic operation library developed based on an MPFR library to obtain a position and speed matrix of the satellite at the moment corresponding to the simulation step length.
2. The pico-meter satellite orbit simulation method for deep space exploration according to claim 1, wherein said basic operation library is developed in an underlying framework by using an MPFR library, wherein said MPFR library comprises: an mp. The basic operation library developed based on the MPFR library comprises: the device comprises a matrix operation module, a track element calculation module, a character string type conversion module, a time conversion module, a coordinate system conversion module, a data input module and a data output module.
3. The pico-meter satellite orbit simulation method for deep space exploration according to claim 2, wherein the step 2) specifically comprises:
step 2-1) setting a required effective digit, satellite orbit simulation starting time, an optimization parameter for satellite orbit initialization, a simulation step length and a parameter of a solver;
step 2-2) converting the initial time of the satellite orbit simulation from Beijing time to julian day time through the time conversion module;
step 2-3) reading the optimization parameters for satellite orbit initialization through the data input module;
step 2-4), acquiring an initialized satellite position and speed matrix corresponding to the starting time of the satellite orbit based on the optimized parameters through the matrix operation module;
step 2-5), converting the initialized satellite position and speed matrix into a coordinate system required in the solving process of the satellite formation dynamic model through the coordinate system conversion module;
step 2-6), performing ordinary differential equation solution on the satellite formation kinetic model by using the set parameters of the solver, the simulation step length, the effective digit and the initialized satellite position and speed matrix of the required coordinate system through the solver to obtain the satellite position and speed with the required effective digit at the moment corresponding to the simulation step length under the coordinate system required by the satellite formation kinetic model; and converting the satellite position and the satellite velocity into an inertial system through the coordinate system conversion module so as to obtain the satellite position and the satellite velocity with the significant digit at the moment corresponding to the simulation step length under the inertial system.
4. The pico-meter scale satellite orbit simulation method for deep space exploration according to claim 1, wherein the type of solver comprises: 8-order variable step length Longge Kutta solver and 8-order variable step length DOP853 solver; the solver is further developed by utilizing the development of the basic operation library based on the MPFR library to carry out floating point operation with any significant digit.
5. The pico-meter scale satellite orbit simulation method for deep space exploration according to claim 1, wherein the solver adopts an 8-order variable step DOP853 solver and adopts a 7-order dense output, wherein the parameters of the solver are set as: tolerance of error is 1 x 10-NWherein N is any positive integer.
6. The pico-meter scale satellite orbit simulation method for deep space exploration according to claim 2, wherein said solver utilizes a stringtoreall function in said string type conversion module to convert solver coefficients in the form of strings and optimized parameters for satellite initialization into an mpreal data type, so that solver coefficients and optimized parameters for satellite initialization are converted into constants corresponding to said significant digits.
7. The method for simulating the pico-meter satellite orbit for deep space exploration according to claim 1, wherein the satellite formation dynamics model is established by using the following formula:
Figure FDA0003538264090000021
wherein i is the serial number of the satellite,
Figure FDA0003538264090000022
an acceleration matrix of the ith satellite in the inertial system;
Figure FDA0003538264090000023
is a matrix of positions of the ith satellite in the inertial system,
Figure FDA0003538264090000024
is a position matrix of the kth celestial body in the inertial system; mkMass of kth celestial body; g is a universal gravitation constant; wherein the ordering of the plurality of celestial bodies is: sun, water star, venus venosus, earth and moon, mars, saturn, asterias, plundersiatus; when the celestial body is the earth and the moon, the mass of the celestial body is the sum of the earth mass and the moon mass.
8. The method for simulating the pico-meter satellite orbit for deep space exploration according to claim 1, wherein the gravitational field of the plurality of celestial bodies comprises: the sun's central gravity, the water star gravity, the golden star gravity, the sum of the earth and moon gravities, the mars gravity, the wooden star gravity, the earth star gravity, the heaven star gravity, the sea king star gravity and the pluto gravity.
9. A pico-meter scale satellite orbit simulation system for deep space exploration, which is used for executing the pico-meter scale satellite orbit simulation method for deep space exploration according to claims 1-8, and comprises the following steps: a solver and a basic operation library, characterized in that,
the basic operation library is developed by utilizing an MPFR library in a bottom layer framework, and comprises the following steps: the device comprises a matrix operation module, a track element calculation module, a character string type conversion module, a time conversion module, a coordinate system conversion module, a data input module and a data output module; wherein the content of the first and second substances,
the time conversion module is used for converting the initial time of the satellite orbit simulation from Beijing time to julian daily time;
the data input module is used for reading the optimization parameters for satellite orbit initialization;
the matrix operation module is used for acquiring an initialized satellite position and speed matrix corresponding to the starting time of the satellite orbit based on the optimization parameters;
the coordinate system conversion module is used for converting the initialized satellite position and speed matrix into a coordinate system required in the solving process of the satellite formation dynamic model;
the solver is used for solving an ordinary differential equation of the satellite formation kinetic model by using the set parameters of the solver, the simulation step length, the effective digit and the initialized satellite position and speed matrix in the required coordinate system so as to obtain the satellite position and speed with the required effective digit of the satellite formation kinetic model in the required coordinate system at the moment corresponding to the simulation step length; and the satellite position and the satellite speed with the significant digit at the moment corresponding to the simulation step length in the inertial system are obtained by converting the satellite position and the satellite speed into the inertial system through the coordinate system conversion module.
10. The pico-meter scale satellite orbit simulation system for deep space exploration according to claim 9, wherein said solver is of the type comprising: 8-order variable-step length Rungestota solver and 8-order variable-step length DOP853 solver; the solver is further developed by using the basic operation library developed based on the MPFR library so as to perform floating point operation with any significant digit.
CN202210230755.5A 2022-03-09 2022-03-09 Picometer-level satellite orbit simulation method and system for deep space exploration Active CN114662294B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210230755.5A CN114662294B (en) 2022-03-09 2022-03-09 Picometer-level satellite orbit simulation method and system for deep space exploration

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210230755.5A CN114662294B (en) 2022-03-09 2022-03-09 Picometer-level satellite orbit simulation method and system for deep space exploration

Publications (2)

Publication Number Publication Date
CN114662294A true CN114662294A (en) 2022-06-24
CN114662294B CN114662294B (en) 2022-12-06

Family

ID=82028980

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210230755.5A Active CN114662294B (en) 2022-03-09 2022-03-09 Picometer-level satellite orbit simulation method and system for deep space exploration

Country Status (1)

Country Link
CN (1) CN114662294B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116502399A (en) * 2023-03-02 2023-07-28 北京理工大学 Satellite orbit generation method and generator based on STK and MATLAB joint simulation

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2016091826A1 (en) * 2014-12-10 2016-06-16 Universite De Perpignan Via Domitia Method of adjusting the precision of a computer program manipulating at least one floating point number
CN106092096A (en) * 2016-06-03 2016-11-09 上海航天控制技术研究所 In high-precision orbital emulation, the satellite position based on iterative approach method determines method
CN106156414A (en) * 2016-06-30 2016-11-23 北京润科通用技术有限公司 A kind of satellite trajectory simulation method and device
CN109255096A (en) * 2018-07-25 2019-01-22 西北工业大学 A kind of uncertain evolution method of the geostationary orbits based on differential algebra
CN111104091A (en) * 2019-12-12 2020-05-05 北京科技大学 Detection and conversion method for precision specific calculation in dynamic floating point error analysis
CN112389679A (en) * 2020-11-03 2021-02-23 西北工业大学 Moonlet constant thrust orbit recursion method considering multiple perturbation forces
CN113128749A (en) * 2021-03-05 2021-07-16 中国科学院国家空间科学中心 Centralized online planning method for satellite observation network

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2016091826A1 (en) * 2014-12-10 2016-06-16 Universite De Perpignan Via Domitia Method of adjusting the precision of a computer program manipulating at least one floating point number
CN106092096A (en) * 2016-06-03 2016-11-09 上海航天控制技术研究所 In high-precision orbital emulation, the satellite position based on iterative approach method determines method
CN106156414A (en) * 2016-06-30 2016-11-23 北京润科通用技术有限公司 A kind of satellite trajectory simulation method and device
CN109255096A (en) * 2018-07-25 2019-01-22 西北工业大学 A kind of uncertain evolution method of the geostationary orbits based on differential algebra
CN111104091A (en) * 2019-12-12 2020-05-05 北京科技大学 Detection and conversion method for precision specific calculation in dynamic floating point error analysis
CN112389679A (en) * 2020-11-03 2021-02-23 西北工业大学 Moonlet constant thrust orbit recursion method considering multiple perturbation forces
CN113128749A (en) * 2021-03-05 2021-07-16 中国科学院国家空间科学中心 Centralized online planning method for satellite observation network

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
WENLIN TANG ET AL.: ""Chang’e 3 lunar mission and upper limit on stochastic background of gravitational wave around the 0.01 Hz band"", 《ADVANCES IN SPACE RESEARCH》 *
夏少鹏: ""液态熔盐堆高精度燃耗算法及钍铀增殖研究"", 《中国博士学位论文全文数据库工程科技Ⅱ辑》 *
张扬: ""基于混合网格的飞行器大迎角DES类数值模拟方法研究"", 《中国博士学位论文全文数据库工程科技Ⅱ辑》 *
方晓松: ""卫星轨道建模与仿真技术研究"", 《中国优秀硕士学位论文全文数据库工程科技Ⅱ辑》 *
罗子人等: ""中国空间引力波探测"太极计划"及"太极1号"在轨测试"", 《深空探测学报》 *
罗子人等: ""空间激光干涉引力波探测"", 《力学进展》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116502399A (en) * 2023-03-02 2023-07-28 北京理工大学 Satellite orbit generation method and generator based on STK and MATLAB joint simulation
CN116502399B (en) * 2023-03-02 2024-01-23 北京理工大学 Satellite orbit generation method and generator based on STK and MATLAB joint simulation

Also Published As

Publication number Publication date
CN114662294B (en) 2022-12-06

Similar Documents

Publication Publication Date Title
Peláez et al. A special perturbation method in orbital dynamics
CN102354123A (en) Cross-platform extendible satellite dynamic simulation test system
CN104298647A (en) Low earth orbit satellite based on-satellite determination method for earth shadow moment forecast
CN114662294B (en) Picometer-level satellite orbit simulation method and system for deep space exploration
Kramer et al. Outgassing-induced acceleration of comet 67P/Churyumov-Gerasimenko
CN115314101B (en) Low-orbit communication satellite constellation rapid modeling method based on parallel computing
Turner An open-source, extensible spacecraft simulation and modeling environment framework
Ye et al. Impacts of platform’s position errors on geolocation for a Moon-based sensor
Ray et al. Drag coefficient model to track variations due to attitude and orbital motion
Petrova et al. Lunar-Based Measurements of the Moon’s Physical Libration: Methods and Accuracy Estimates
Merisio et al. Characterization of ballistic capture corridors aiming at autonomous ballistic capture at Mars
Arora et al. A fast, accurate, and smooth planetary ephemeris retrieval system
De Pascuale et al. Simulations of Van Allen Probes plasmaspheric electron density observations
Mauland et al. Sesame: A power spectrum emulator pipeline for beyond-ΛCDM models
Putney et al. Precision orbit determination at the NASA Goddard Space Flight Center
Kelly et al. Geostationary debris mitigation using minimum time solar sail trajectories with eclipse constraints
Sales Trajectory optimization for spacecraft collision avoidance
Wawrzyniak et al. An adaptive, receding-horizon guidance strategy for solar sail trajectories
Yuan et al. Multisensor Integrated Autonomous Navigation Based on Intelligent Information Fusion
Acton et al. Extending NASA's SPICE ancillary information system to meet future mission needs
Dargent Automatic minimum principle formulation for low thrust optimal control in orbit transfers using complex numbers
Huckfeldt Study of space-environmental effects on interplanetary trajectories
Arora High performance algorithms to improve the runtime computation of spacecraft trajectories
McArdle et al. Point Mascon Global Lunar Gravity Models
KR102039108B1 (en) Apparatus for transformation of coordinate used in satellite and operation method of the apparatus

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