Manifold splicing-based earth-moon orbit optimization design method and system under four-pulsar calendar model

Abstract

The invention discloses an earth-moon orbit optimization design method and system under a four-pulsar calendar model based on manifold splicing, and relates to the technical field of detector transfer design between earth-moon special orbits. According to the technical key points, the method comprises the following steps: processing Halo orbit initial values corresponding to a sun-earth system and an earth-moon system based on an invariant manifold, and respectively obtaining Halo orbit optimal manifold splicing points corresponding to the sun-earth system and the earth-moon system; constructing a double-circle restrictive four-body model based on four pulses; and taking the Halo orbit optimal manifold splicing points corresponding to a sun-earth system and an earth-moon system as initial values, and solving the four-pulse-based double-circle restrictive four-body model by using a particle swarm optimization algorithm and a sequential quadratic programming algorithm to obtain an optimal splicing orbit. According to the method, the low-energy transfer orbit is designed by using two optimization models, iteration can be accelerated, the calculation efficiency and the algorithm convergence can be improved, the low-energy transfer orbit in a specific launching window can be quickly searched, and the requirement of earth-moon low-energy orbit transfer of a space station is met.

Description

Technical Field

[0001] The present invention relates to the technical field of design of probe transfer between special Earth-Moon orbits, and in particular to an Earth-Moon orbit optimization design method and system based on a four-pulse ephemeris model with manifold splicing. Background Art

[0002] The unique space environment observation position and low-energy transit station for interplanetary exploration make the application value of the libration point orbit in space missions gradually obvious. The libration point orbit can provide a foothold for the deployment of lunar exploration relay satellites and support the navigation mission of landers and rovers on the far side of the moon. Deploying a manned space station on the libration point orbit facilitates the transportation of materials and fuel replenishment, as well as the implementation of deep space exploration missions such as Mars exploration and asteroid exploration. At the same time, with the help of invariant manifolds, low-energy transfer of the spacecraft can be achieved, effectively improving the ability of the spacecraft to carry payloads. my country's manned lunar space station is an important part of the mission of realizing manned lunar landing and long-term lunar base after the three-step lunar exploration of "orbiting, landing, and returning".

[0003] For the manifold splicing problem, most current studies focus on splicing the invariant manifold under the Earth-Moon circular restricted three-body model with the invariant manifold under the Sun-Earth circular restricted three-body model. The transfer problem is based on the Halo orbit between the Sun-Earth libration point and the Earth-Moon libration point. However, when designing low-energy transfer orbits, the traditional manifold splicing method is difficult to target under the ephemeris model, sensitive to the initial value of the transfer orbit of the three-body model, difficult to solve multi-point targeting, and slow and difficult to solve when optimizing and correcting under the ephemeris model. These limitations have largely hindered the diversity of transfer orbits and made the design of transfer orbits difficult. Therefore, considering the needs of manned and material transportation, the transfer problem from low-Earth orbit to the Halo orbit where the space station is located needs to be solved urgently. Summary of the invention

[0004] In view of this, the present invention proposes a method and system for optimizing the Earth-Moon orbit under a four-pulse ephemeris model based on manifold splicing, so as to solve the technical problems that the traditional method has a long transfer time and cannot achieve low-energy transfer of the spacecraft.

[0005] According to one aspect of the present invention, a method for optimizing the Earth-Moon orbit based on a four-pulse ephemeris model based on manifold splicing is proposed, comprising the following steps:

[0006] Step 1: Based on the invariant manifold, the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are processed to obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; the optimal manifold splicing points of the Halo orbits include the optimal position and optimal speed of the spacecraft;

[0007] Step 2: Construct a two-circle restrictive four-body model based on four pulses; take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, and use the particle swarm optimization algorithm and sequential quadratic programming algorithm to solve the two-circle restrictive four-body model based on four pulses to obtain the optimal splicing orbit.

[0008] Furthermore, the process of obtaining the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system in step 1 includes:

[0009] Step 11: Calculate and obtain the initial value array of the Halo orbit of the Sun-Earth system or the Earth-Moon system; the initial value array includes the positions and velocities of multiple discrete points during the flight of the spacecraft;

[0010] Step 12: using a numerical method to solve the initial value array of the Halo orbit to obtain an invariant manifold of the Halo orbit; the invariant manifold includes an unstable manifold and a stable manifold;

[0011] Step 1: Use the point selection method to splice the invariant manifold of the Halo orbit and calculate the optimal manifold splicing point for the comprehensive position deviation and velocity deviation.

[0012] Furthermore, the specific process of step one includes:

[0013] Under the circular restricted three-body dynamics model, the Halo orbit starts from the xz plane, and its initial value is calculated by Richardson's third-order analytical solution;

[0014] The initial value is subjected to a cycle of dynamic integration to obtain the final state; the error between the final state and the target state is calculated. If the error is greater than the given error, the final state is substituted into the differential correction solution equation to obtain the initial value of the next iteration, i.e. the corrected initial value;

[0015] Repeat the above process until the error between the final state and the target state is less than the given error, and then terminate the iteration to obtain the final corrected initial value.

[0016] The final corrected initial value is subjected to a period of dynamic integration to obtain the Halo orbit initial value array.

[0017] Furthermore, the differential correction solution equation in step 1 is expressed as:

[0018]

[0019] In the formula, is the state transfer matrix, which represents the linear estimation relationship between the initial state and the final state of the system; x0 represents the position of the spacecraft on the X axis; represents the speed of the spacecraft on the Y axis; Indicates the velocity of the spacecraft on the X axis at the final state; represents the acceleration of the spacecraft on the X axis at the final state; y f Indicates the position of the spacecraft on the Y axis at the final state; Indicates the velocity of the spacecraft on the Y axis at the final state; Indicates the velocity of the spacecraft on the Z axis at the final state; It represents the acceleration of the spacecraft on the Z axis at the final state; δ represents the correction amount, and the variable after it represents the correction amount of each iteration of the variable.

[0020] Furthermore, the specific process of steps one and two includes:

[0021] The state transfer matrix corresponding to the final state obtained by the dynamic integration of the final corrected initial value after one cycle is expressed as a single-value matrix, and the single-value matrix is ​​solved to obtain the eigenvalues ​​and eigenvectors; the eigenvalues ​​include: λ1>1 and

[0022] The initial value point of the unstable manifold is calculated by the following formula: In the formula, p i represents the column data in the Halo orbit initial value array; ε represents the perturbation factor; D u represents the eigenvector corresponding to the eigenvalue λ1;

[0023] The initial value point of the integration of the stable manifold is obtained by the following calculation: Where D s represents the eigenvector corresponding to the eigenvalue λ2;

[0024] After obtaining the initial value points of the integration of the two manifolds, the circular restricted three-body dynamics model is used for recursion to obtain the invariant manifold array of the Halo orbit; the invariant manifold array includes an unstable manifold array and a stable manifold array.

[0025] Furthermore, the specific process of step 1-3 includes:

[0026] Calculate the minimum position deviation Δr between the position in the unstable manifold array and the position in the stable manifold array i ;

[0027] Calculate the minimum speed deviation Δv between the speed in the unstable manifold array and the speed in the stable manifold array i ;

[0028] For the minimum position deviation Δr i Deviation from minimum speed Δv i Normalize and perform weighted summation to obtain the comprehensive deviation

[0029] right Sort and select the point corresponding to the smallest comprehensive deviation, which is the optimal manifold splicing point of the comprehensive position deviation and velocity deviation.

[0030] Furthermore, the four-pulse dual-circle restrictive four-body model in step 4 includes an upper optimization model and a lower optimization model; wherein the upper optimization model is expressed as:

[0031]

[0032] The lower layer optimization model is expressed as:

[0033]

[0034] In the formula, min means minimization; and They represent the optimal speeds at the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; ε s represents the disturbance factor corresponding to the Sun-Earth system, ε m represents the disturbance factor corresponding to the Earth-Moon system; Indicates the preset disturbance factor threshold; and They respectively represent the optimal positions of the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; t1 represents the time it takes for the spacecraft in the Halo orbit of the Sun-Earth system to run from the final corrected initial value to the optimal manifold splicing point, and t2 represents the time it takes for the spacecraft in the Halo orbit of the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

[0035] Furthermore, in step 2, the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are used as initial values, and the particle swarm optimization algorithm and the sequential quadratic programming algorithm are used to solve the four-pulse based double-circle restricted four-body model to obtain the optimal splicing orbit, which includes:

[0036] Step 21: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as the initial value, use the particle swarm optimization algorithm to solve the upper optimization model, and output the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; wherein, in the iterative optimization process of the particle swarm optimization algorithm, set a perturbation factor threshold, and at the end of each iteration, determine whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, punish the optimization result, that is, increase its performance index;

[0037] Step 22: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system, as well as the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, use the sequential quadratic programming algorithm to solve the lower-level optimization model, and output the time it takes for the spacecraft in the Halo orbits of the Sun-Earth system and the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

[0038] Furthermore, the method also includes step three: after obtaining the optimal stitched orbit, correcting the position and state of the data points in the optimal stitched orbit based on the four-pulse ephemeris model to obtain the optimal stitched orbit under the four-pulse ephemeris model.

[0039] According to another aspect of the present invention, a system for optimizing the Earth-Moon orbit based on a four-pulse ephemeris model based on manifold splicing is proposed, comprising:

[0040] The optimal manifold splicing point calculation module is configured to process the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system based on the invariant manifold, and obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; the optimal manifold splicing points of the Halo orbits include the optimal position and the optimal speed of the spacecraft; the specific process includes: calculating and obtaining the initial value array of the Halo orbits of the Sun-Earth system or the Earth-Moon system; the initial value array includes the position and speed of multiple discrete points during the flight of the spacecraft; using a numerical method to solve the initial value array of the Halo orbit to obtain the invariant manifold of the Halo orbit; the invariant manifold includes an unstable manifold and a stable manifold; using a point selection method to splice the invariant manifolds of the Halo orbit, and calculating and obtaining the optimal manifold splicing points of the comprehensive position deviation and speed deviation;

[0041] The optimal splicing orbit acquisition module is configured to construct a four-pulse dual-circle restricted four-body model; take the optimal manifold splicing points of the Halo orbits corresponding to the sun-earth system and the earth-moon system as initial values, use a particle swarm optimization algorithm and a sequential quadratic programming algorithm to solve the dual-circle restricted four-body model based on the four-pulse to obtain the optimal splicing orbit; specifically comprising: taking the optimal manifold splicing points of the Halo orbits corresponding to the sun-earth system and the earth-moon system as initial values, using a particle swarm optimization algorithm to solve the upper optimization model, and outputting the optimal perturbation factors of the Halo orbits corresponding to the sun-earth system and the earth-moon system; wherein, in the particle swarm optimization algorithm, During the iterative optimization process of the method, a perturbation factor threshold is set. At the end of each iteration, it is determined whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, the optimization result is punished, that is, its performance index is increased. The optimal manifold splicing point of the Halo orbit corresponding to the sun-earth system and the earth-moon system, as well as the optimal perturbation factor of the Halo orbit corresponding to the sun-earth system and the earth-moon system are taken as initial values, and the sequential quadratic programming algorithm is used to solve the lower optimization model, and the time taken for the spacecraft in the Halo orbit of the sun-earth system and the earth-moon system to run from the final corrected initial value to the optimal manifold splicing point is output.

[0042] The beneficial technical effects of the present invention are:

[0043] The low-energy transfer orbit design method under the high-precision ephemeris model proposed in the present invention uses two optimization models to design the low-energy transfer orbit, which can accelerate iteration to a certain extent and improve calculation efficiency. At the same time, the initial value has a low impact on the optimization result, which improves the convergence of the algorithm and realizes the rapid search for the low-energy transfer orbit in a specific launch window, meeting the needs of the space station for the transfer of the Earth-Moon low-energy orbit. The universal optimization method and model provided by the present invention can efficiently use the manifold under the three-body model to design the low-energy transfer orbit. This method can not only be used as a design method for the cargo transfer orbit of the future lunar base construction, but also can be used as a reference for the deployment of the future lunar space station, promoting the development and technological progress of these fields. BRIEF DESCRIPTION OF THE DRAWINGS

[0044] The present invention can be better understood by referring to the description given below in conjunction with the accompanying drawings, which together with the following detailed description are included in this specification and form a part of this specification, and are used to further illustrate the preferred embodiments of the present invention and explain the principles and advantages of the present invention.

[0045] Figure 1 It is a flow chart of the method for optimizing the Earth-Moon orbit under the four-pulse ephemeris model based on manifold splicing described in an embodiment of the present invention.

[0046] Figure 2 Schematic diagram of the transfer track design concept in an embodiment of the present invention.

[0047] Figure 3 Schematic diagram of the periodic Halo orbit generated in an embodiment of the present invention.

[0048] Figure 4 This is an example diagram of a special orbital invariant manifold in an embodiment of the present invention.

[0049] Figure 5 Schematic diagram of error deviation of optimal manifold splicing points in an embodiment of the present invention.

[0050] Figure 6 The present invention is a flowchart for solving a four-pulse based dual-circle restricted four-body model in an embodiment of the present invention.

[0051] Figure 7 Schematic diagram of the transfer orbit from the Sun-Earth Halo orbit to the Earth-Moon Halo orbit in an embodiment of the present invention.

[0052] Figure 8 Schematic diagram of the transfer orbit under the double-circle restricted four-body model based on four pulses in an embodiment of the present invention.

[0053] Fig. 9 Schematic diagram of low-energy transfer orbit under the high-precision ephemeris model in an embodiment of the present invention. DETAILED DESCRIPTION

[0054] In order to enable those skilled in the art to better understand the scheme of the present invention, exemplary implementations or embodiments of the present invention will be described below in conjunction with the accompanying drawings. Obviously, the described implementations or embodiments are only implementations or embodiments of a part of the present invention, not all of them. Based on the implementations or embodiments of the present invention, all other implementations or embodiments obtained by ordinary technicians in the field without creative work should fall within the scope of protection of the present invention.

[0055] To solve the problem of transferring from LEO to the Earth-Moon Halo orbit under the high-precision ephemeris model, multi-point shooting is usually used to correct the initial transfer orbit converted to the high-precision ephemeris model. [1] In order to save energy, the initial transfer orbit usually adopts the design method of invariant manifold. [2]However, the simple invariant manifold of the Earth-Moon system is difficult to approach large celestial bodies, and the local optimization design method and multi-point shooting have strict requirements on the initial values, low convergence, and it is difficult to select the number of points for multi-point shooting. Therefore, the present invention proposes an Earth-Moon orbit optimization design method and system based on manifold splicing under a four-pulse ephemeris model, which designs the initial value of the transfer orbit by splicing the Halo orbits of the Sun-Earth system and the Earth-Moon system, and then solves the proposed optimization model by combining the hybrid optimization method of the particle swarm algorithm and the sequential quadratic programming algorithm (SQP) to further optimize the initial value, thereby obtaining the transfer orbit under the double-circle restricted four-body model. The optimized orbit has high convergence and high computational efficiency under the multi-point shooting method, and the initial value conjecture of the method is simple, which can approach the transfer orbit with the minimum fuel, and has potential application value in the future design of the Earth-Moon transfer orbit and the future utilization and development of the moon.

[0056] The embodiment of the present invention provides a method for optimizing the Earth-Moon orbit under a four-pulse ephemeris model based on manifold splicing, such as Figure 1 As shown, the method comprises the following steps:

[0057] Step 1: Based on the invariant manifold, the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are processed to obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; the optimal manifold splicing points of the Halo orbits include the optimal position and optimal speed of the spacecraft; specifically, the following steps are included:

[0058] Step 11: Calculate and obtain the initial value array of the Halo orbit of the Sun-Earth system or the Earth-Moon system; the initial value array includes the positions and velocities of multiple discrete points during the flight of the spacecraft;

[0059] Step 12: using a numerical method to solve the initial value array of the Halo orbit to obtain an invariant manifold of the Halo orbit; the invariant manifold includes an unstable manifold and a stable manifold;

[0060] Step 13: Use the point selection method to splice the invariant manifold of the Halo orbit and calculate the optimal manifold splicing point of the comprehensive position deviation and velocity deviation;

[0061] Step 2: construct a four-pulse dual-circle restricted four-body model; take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, and use the particle swarm optimization algorithm and the sequential quadratic programming algorithm to solve the four-pulse dual-circle restricted four-body model to obtain the optimal splicing orbit; specifically including:

[0062] Step 21: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as the initial value, use the particle swarm optimization algorithm to solve the upper optimization model, and output the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; wherein, in the iterative optimization process of the particle swarm optimization algorithm, set a perturbation factor threshold, and at the end of each iteration, determine whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, punish the optimization result, that is, increase its performance index;

[0063] Step 22: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system, as well as the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, use the sequential quadratic programming algorithm to solve the lower-level optimization model, and output the time it takes for the spacecraft in the Halo orbits of the Sun-Earth system and the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

[0064] The method starts from step 1. In step 1, based on the invariant manifold, the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are processed to obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively.

[0065] According to an embodiment of the present invention, Figure 2 As shown, Figure 2 The transfer orbit design idea of ​​the present invention is shown. First, in step 11, the Halo orbit initial value array of the sun-earth system or the earth-moon system is calculated and obtained. The specific process includes: under the circular restricted three-body dynamics model, the Halo orbit starts from the xz plane, and its initial value is obtained by Richardson's third-order analytical solution; the initial value is subjected to a cycle of dynamic integration to obtain the final state; the error between the final state and the target state is calculated, and if the error is greater than the given error, the final state is substituted into the differential correction solution equation to calculate the initial value of the next iteration, that is, the corrected initial value; the above process is repeated until the error between the final state and the target state is less than the given error, and the iteration is terminated to obtain the final corrected initial value; the final corrected initial value is subjected to a cycle of dynamic integration to obtain the Halo orbit initial value array.

[0066] Specifically, in the circular restricted three-body dynamics model coordinate system, the center of mass motion of the spacecraft is:

[0067]

[0068] Where: x, y, z are the three-axis coordinates of the spacecraft in the normalized coordinate system; is the velocity of the spacecraft in the three-axis coordinate system in the normalized coordinate system; is the acceleration of the three-axis coordinates of the spacecraft in the normalized coordinate system; μ is the mass parameter of the system; r1 is the distance from the spacecraft to the earth; r2 is the distance from the spacecraft to the moon.

[0069] The Halo orbit is a kind of libration point orbit under the circular restricted three-body dynamics model. According to the third-order approximate analytical solution given by Richardson, the Halo orbit under the circular restricted three-body dynamics model is determined by three parameters: amplitude parameter, north-south family parameter and phase parameter. The amplitude parameter is generally represented by the normal amplitude, the x-axis amplitude A x and the y-direction amplitude A y Available through A z The phase parameter φ determines the specific phase of the initial point on the Halo orbit. Richardson's third-order analytical solution is:

[0070]

[0071] Due to the symmetry of the three-body problem, Halo orbits are divided into two families, the northern and southern families, corresponding to the two cases where n is 1 and 3 in the third-order analytical solution. For the L1 point, n = 1 corresponds to the northern family Halo orbit; for the L2 point, n = 3 corresponds to the northern family Halo orbit. The northern family Halo orbit moves above the motion plane of the main celestial body most of the time, so it has good coverage characteristics for the northern hemisphere of the main celestial body. Similarly, the southern family Halo orbit has good coverage characteristics for the southern hemisphere of the main celestial body. Therefore, in the design of the space mission orbit, the corresponding Halo orbit can be selected according to the actual situation, such as the landing point in the southern (northern) hemisphere.

[0072] In the conjunction coordinate system of the circular restricted three-body dynamics model, the Halo orbit is symmetric about the xz plane, so when the orbit passes through the xz plane, the velocity is perpendicular to the xz plane, that is, the velocities in the x and z directions are both 0. Based on this characteristic of the Halo orbit, the differential correction solution equation is given.

[0073] Assume that the Halo orbit starts from the xz plane, and the initial state quantity X0 has the following form:

[0074]

[0075] Take X0 as the initial value of the integration to perform orbital integration. When the orbit returns to the xz plane again, the integration time is considered to be half a period t f , assuming that the final state after half a period of integration is:

[0076]

[0077] Assumptions is the initial correction value of the initial state quantity X0. After half a cycle t f The change at the end is:

[0078]

[0079] Expand to:

[0080]

[0081] The integral termination condition is to integrate to the xz surface, which can be written as an expression:

[0082]

[0083] The differential correction solution equation is as follows:

[0084]

[0085] In the formula, Φ is the state transfer matrix, which is a 6×6 square matrix. ij Where i is the corresponding row number of the state transfer matrix, and j is the corresponding column number of the state transfer matrix. is the state transfer matrix, which represents the linear estimation relationship between the initial state and the final state of the system; x0 represents the position of the spacecraft on the X axis; represents the speed of the spacecraft on the Y axis; Indicates the velocity of the spacecraft on the X axis at the final state; represents the acceleration of the spacecraft on the X axis at the final state; f Indicates the position of the spacecraft on the Y axis at the final state; Indicates the velocity of the spacecraft on the Y axis at the final state; Indicates the velocity of the spacecraft on the Z axis at the final state; It represents the acceleration of the spacecraft on the Z axis at the final state; δ represents the correction amount, and the variable after it represents the correction amount of each iteration of the variable.

[0086] In summary, the initial value obtained by Richardson's third-order analytical solution is iterated using differential correction to solve the equation. The iteration process is the Newton iteration process until the iteration variable is less than the error. The specific process is as follows: the initial value obtained by Richardson's third-order analytical solution is the initial value of the Halo orbit X0. After dynamic integration to the integral termination condition, the final state is X f ; Calculate the final state X f If the error is greater than the given error, the iteration continues until the final state error is less than the given error, and the final corrected initial value X'0 is output. f |≤ε1, ε1 represents a small value, usually 1×10 -7 ; The constraint variables are and The control variables are x0, z0 and By adjusting x0, z0 and Make and If it is less than ε2, the orbit is considered to be periodic. ε2 is expressed as a small value, usually 1×10 -7 The initial value of the third-order analytical solution requires the z-axis amplitude parameter of Halo. The z-axis amplitude parameter is related to z0. Therefore, z0 is fixed, x0 and Figure 3 An example diagram of using differential correction to solve the equation for iteration is shown. The final corrected initial value X'0 is shown in Table 1.

[0087] Table 1 Examples of initial values ​​of Halo orbits corresponding to the Earth-Moon system

[0088]

[0089] Finally, the final corrected initial value is subjected to a period of dynamic integration to obtain the Halo orbit initial value array.

[0090] Then, in step 1 and 2, the Halo orbit initial value array is solved by numerical methods to obtain the invariant manifold of the Halo orbit; the invariant manifold includes unstable manifold and stable manifold. The specific process includes: the state transfer matrix corresponding to the final state obtained by the dynamic integration of the final corrected initial value after one cycle is expressed as a single-value matrix, and the single-value matrix is ​​solved to obtain the eigenvalue and eigenvector; the eigenvalue includes: λ1>1 and The initial value point of the unstable manifold is calculated by the following formula: In the formula, p i represents the column data in the Halo orbit initial value array; ε represents the perturbation factor; D u represents the eigenvector corresponding to the eigenvalue λ1; the initial value point of the integration of the stable manifold is obtained by the following calculation: Where D s Represents the eigenvector corresponding to the eigenvalue λ2; after obtaining the initial value points of the integration of the two manifolds, the circular restricted three-body dynamics model is used for recursion to obtain the invariant manifold array of the Halo orbit; the invariant manifold array includes an unstable manifold array and a stable manifold array.

[0091] When a spacecraft is performing orbit transfer, it needs to apply a velocity pulse to itself so that the spacecraft changes its velocity state and enters a transfer orbit. An embodiment of the present invention utilizes an invariant manifold and a splicing design of invariant manifolds to optimize a four-pulse transfer orbit, wherein the four pulses correspond to a pulse from an Earth parking orbit to a Sun-Earth stable manifold, a pulse from a Sun-Earth stable manifold to a Sun-Earth unstable manifold, a pulse from a Sun-Earth unstable manifold to an Earth-Moon stable manifold, and a pulse from an Earth-Moon stable manifold to a lunar orbit.

[0092] The set of Halo initial value arrays obtained in steps one by one after a period of dynamic integration is represented as P1. The dimension of P1 is 6×N, where N is the number of discrete points for discrete selection of the Halo orbit, that is, one period of Halo orbit selection is represented by N points. The calculation of the invariant manifold of the Halo orbit requires the integral initial value of the invariant manifold, and the integral initial value of the invariant manifold requires the use of a single-value matrix. The definition of the single-value matrix is: The state transfer matrix running along the Halo orbit for one cycle is called the single-value matrix M, that is:

[0093]

[0094] Where Φ(t0;T) is the state transfer matrix after running for a period T.

[0095] The next step is to calculate the eigenvalues ​​and eigenvectors of the single-value matrix. The eigenvalues ​​of the single-value matrix have 6 forms as follows. For the same Halo orbit, different points have the same single-value matrix: λ1>1; λ3=λ4=1; |λ5|=|λ6=1|.

[0096] The calculation of the unstable manifold requires the eigenvector corresponding to the eigenvalue λ1. Then any point p in the Halo initial value array obtained in step 1 is i , the initial value of the integral of the unstable manifold corresponding to i∈(0,1,2…,N) can be obtained by the following formula:

[0097]

[0098] The calculation of the stable manifold requires the eigenvector corresponding to the eigenvalue λ2. Then any point p in the Halo initial value array obtained in step 1 is i The corresponding initial value of the integral of the stable manifold can be obtained by the following formula:

[0099]

[0100] Among them, ε is the perturbation factor, which is usually designed through experience; the ± signs correspond to the two branches of the manifold, one flowing to the large celestial body and the other flowing to the small celestial body.

[0101] After obtaining the initial value point of the integration of the manifold, the invariant manifold array of the Halo orbit can be obtained by recursion according to the circular restricted three-body dynamics model (this process is the prior art); the dimension of the array is N×6×N1, where the array data of the manifold obtained at any point on the Halo orbit is 6×N1, and N1 is the number of points discretized on each manifold. Examples of invariant manifolds of the Halo orbit are as follows Figure 4 shown.

[0102] Then, in step 1-3, the invariant manifold of the Halo orbit is spliced ​​using the point selection method to calculate the optimal manifold splicing point with comprehensive position deviation and velocity deviation. The specific process includes: calculating the minimum position deviation Δr between the position in the unstable manifold array and the position in the stable manifold array i ; Calculate the minimum speed deviation Δv between the speed in the unstable manifold array and the speed in the stable manifold array i ; For the minimum position deviation Δr i Deviation from minimum speed Δv i Normalize and perform weighted summation to obtain the comprehensive deviation right Sort and select the point corresponding to the smallest comprehensive deviation, which is the optimal manifold splicing point of the comprehensive position deviation and velocity deviation.

[0103] According to steps one and two, the invariant manifolds of the Halo orbit of the Sun-Earth system and the Halo orbit of the Earth-Moon system are generated. The difference between the generation of the Earth-Moon Halo orbit and the Sun-Earth Halo orbit is that the normalized quantities such as distance, time and speed are converted into the normalized quantities of the Sun-Earth system, and the mass ratio in dynamics is converted into the mass ratio of the Sun-Earth system to obtain the Sun-Earth Halo orbit. The invariant manifold of the Halo orbit of the Sun-Earth system is spliced ​​with the Halo orbit of the Earth-Moon system. A very important point in the manifold splicing process is to select the splicing point of the two manifolds. Although the manifolds on both sides end at the Poincare section, the Poincare section is a two-dimensional plane under the circular restrictive three-body model of the Sun and the Earth, and the yz plane with x=0 is selected here, the manifolds on both sides of the Poincare section still cannot guarantee that the position deviation and velocity deviation can meet the requirements of low-energy orbit transfer. Therefore, a point selection method is given below to ensure that the obtained manifolds can be effectively spliced.

[0104] The Halo orbit has been discretized in the above, and the discrete points can be expressed as:

[0105] X=[X1,X2,...,X N ]

[0106] Perform perturbation integration on each starting discrete point, that is, solve the unstable manifold corresponding to each point, and use the Poincare section as the integral termination condition. The state matrix representation of each manifold on the Poincare section is as follows:

[0107] X u =[X u1 ,X u2 ,...,X uN ]

[0108] Where each column represents a point state of an unstable manifold:

[0109]

[0110] Similarly, for the target libration point orbit, the reverse integration after perturbation can be performed on each target discrete point to solve the stable manifold corresponding to each point. The state matrix representation of each manifold on the Poincare section is obtained by taking the Poincare section as the integral termination condition as follows:

[0111] X s =[X s1 ,X s2 ,...,X sN ]

[0112] Among them, each column represents a point state of the stable manifold:

[0113]

[0114] The minimum position and velocity deviations are calculated based on the manifold section state of the starting track, and the expressions are as follows:

[0115]

[0116] In the formula, r ui is the position of the end of the unstable manifold at any point on the Poincare section, r sj To stabilize the position of the end of the manifold at any point on the Poincare section, v ui is the velocity of any point on the Poincare section at the end of the unstable manifold, v sj is the velocity at any point on the Poincare section at the end of the stable manifold.

[0117] The deviation is normalized to the maximum value, so that the position and speed deviations are unified to the range of [0,1]:

[0118]

[0119] The position deviation and speed deviation are weighted and summed to obtain the comprehensive deviation to balance the influence of the two deviations on the state at this point. The mathematical expression is as follows:

[0120]

[0121] right Sort and select the point corresponding to the smallest comprehensive deviation, which is the optimal manifold splicing point obtained by combining the position deviation and velocity deviation. The error deviation example is as follows Figure 5 shown.

[0122] Then, step 2 is executed. In step 2, a four-pulse dual-circle restricted four-body model is constructed; the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are used as initial values, and the particle swarm optimization algorithm and sequential quadratic programming algorithm are used to solve the four-pulse dual-circle restricted four-body model to obtain the optimal splicing orbit. Specifically, it includes:

[0123] Step 21: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as the initial value, use the particle swarm optimization algorithm to solve the upper optimization model, and output the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; wherein, in the iterative optimization process of the particle swarm optimization algorithm, set a perturbation factor threshold, and at the end of each iteration, determine whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, punish the optimization result, that is, increase its performance index;

[0124] Step 22: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system, as well as the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, use the sequential quadratic programming algorithm to solve the lower-level optimization model, and output the time it takes for the spacecraft in the Halo orbits of the Sun-Earth system and the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

[0125] According to the embodiment of the present invention, the optimal manifold splicing point of the comprehensive position deviation and velocity deviation is obtained by step 1, and the time t0 of the Earth-Moon Halo orbit and the Sun-Earth Halo orbit from the starting point to the corresponding optimal splicing point can be obtained. 1 and t0 2 The result and the optimal manifold splicing point are used as the initial optimization input. The upper layer uses the particle swarm optimization algorithm for global optimization, and the optimization variables are the perturbation factors corresponding to the Earth-Moon manifold and the Sun-Earth manifold. The lower layer uses the sequential quadratic programming algorithm-SQP algorithm for local optimization, and the perturbation factor output by the upper layer is used as one of the inputs of the lower model.

[0126] The specific process is: select the initial value of the differentially corrected Halo orbit as the starting point in the Sun-Earth Halo orbit and the Earth-Moon Halo orbit. and The final state point reached by the Halo orbit after running from the starting point for t1 and t2 respectively is and and is the starting point of the Sun-Earth manifold and the Earth-Moon manifold. The position of the Sun-Earth manifold and the Earth-Moon manifold when they run to the Poincare section is and Speed ​​status is and The manifold of the Earth-Moon system determines the conversion matrix C(x) of the Earth-Moon three-body system to the Sun-Earth three-body system according to the given epoch start time. The optimization variables of this problem are the time t1 in the Sun-Earth Halo orbit and the time t2 in the Earth-Moon Halo orbit. The optimization goal is to minimize the velocity error Δv3 on the Poincare section. The constraints of the optimization model are bounded constraints on the absence of position error and deviation at the Poincare section. The optimal splicing orbit can be obtained by solving this optimization problem. The solution result of the lower model is affected by the output result of the upper model. When the solution value of the manifold deviation is unreasonable, the result of the lower model does not converge. A large index value is assigned to it as a penalty, and the next iterative solution is entered at the same time. The flowchart of the optimization model is as follows: Figure 6 shown.

[0127] The mathematical description of the optimization model (i.e., the four-pulse two-circle restricted four-body model) is as follows:

[0128] The upper optimization model is:

[0129]

[0130] The lower optimization model is:

[0131]

[0132] In the formula, min means minimization; and They represent the optimal speeds at the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; ε s represents the disturbance factor corresponding to the Sun-Earth system, ε m represents the disturbance factor corresponding to the Earth-Moon system; Indicates the preset disturbance factor threshold; and They respectively represent the optimal positions of the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; t1 represents the time it takes for the spacecraft in the Halo orbit of the Sun-Earth system to run from the final corrected initial value to the optimal manifold splicing point, and t2 represents the time it takes for the spacecraft in the Halo orbit of the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

[0133] Among them, the position and velocity update formula of the particle swarm optimization algorithm is:

[0134]

[0135] Where: is the velocity component of the particle at the kth iteration; is the position component of the particle at the kth iteration; ω is the inertia factor; c1 and c2 are acceleration factors; is the individual extreme value of particle i after the first k iterations; is the population extreme value of the entire particle population after the first k iterations; R is a random number. The update strategy of the acceleration factor is as follows: three different autonomous populations are designed, and the update strategies of the acceleration factors of the three populations are shown in Table 2.

[0136] Table 2 Update strategy for c1 and c2

[0137]

[0138] Select the number of particles, the maximum number of iterations, the maximum inertia weight coefficient, the minimum inertia weight coefficient, The double-layer optimization algorithm is used to optimize the transfer orbit from the Sun-Earth Halo orbit to the Earth-Moon Halo orbit, as Figure 7 As shown in the figure, the result of the transfer orbit is a 6×N array.

[0139] The point on the Sun-Earth Halo orbit corresponding to the initial point of the Sun-Earth system stable manifold in the obtained optimal spliced ​​orbit is taken as the starting point, and the unstable manifold of this point is calculated. The transfer orbit from the Sun-Earth Halo orbit to the LEO orbit is obtained by differential correction to satisfy the parking orbit height constraint. The obtained optimal spliced ​​orbit is spliced ​​with the corrected unstable manifold, as shown in Figure 8 The obtained transfer orbits are all 7×N arrays, the 7th dimension data is the time corresponding to the transfer, and the first 6 dimensions are the position and speed.

[0140] Furthermore, the method also includes step three: after obtaining the optimal stitched orbit, correcting the position and state of the data points in the optimal stitched orbit based on the four-pulse ephemeris model to obtain the optimal stitched orbit under the four-pulse ephemeris model.

[0141] According to the embodiment of the present invention, the optimized initial value is the transfer orbit of the double-circle restrictive four-body obtained in step 2, which is a 7×N array. Select an initial epoch time, transform the transfer orbit into the geocentric inertial coordinate system, convert the time into the corresponding ephemeris time, and discretize the array of the transfer orbit. For example, select the number of target points as k=20, divide N into 20 equal parts, and select the target points in turn, so as to obtain the position velocity information and the corresponding ephemeris time of each target point.

[0142] The specific process is as follows:

[0143] (1) Continuous position correction: The state quantity of the k points given by the current iteration value and the integration time of each segment are used for orbit recursion. At this time, the speed of the starting point of this segment is changed according to the position correction equation to make the end point of this segment coincide with the position of the starting point of the next segment. The control variable is the speed component of the starting point of this segment. The constraint variable is the position r of the end point if =[x if ,y if ,z if ] T , the termination condition is that the integration time is fixed. The position correction equation can be obtained:

[0144]

[0145] Repeat this correction operation for each segment of the trajectory until all segments of the trajectory are in continuous position.

[0146] (2) Speed ​​Continuity Correction: After the first level of differential correction, i.e. position continuous correction, the trajectory is continuous in position, but there will be a speed jump at each segment connection point (existing at every maneuvering point). In the second level of differential correction, the speed will be corrected to make the trajectory continuous in speed. The corresponding speed differential correction relationship is:

[0147]

[0148] After completing the continuous velocity correction, since the position of each segment end point is also changed, continuous position correction is required. This process is repeated iteratively until the position error and velocity error are within the allowable range. It can be considered that a transfer trajectory with continuous position and velocity is obtained.

[0149] Because the transfer orbit needs to meet engineering constraints and the starting and ending points need to be fixed during the correction process, the calculation efficiency of the transfer orbit is slow and the convergence speed decreases when the number of correction iterations is large. Therefore, a penalty function is set for the velocity termination condition of the corrected orbit to obtain an inaccurate convergence solution and then substitute it into the optimization model using the sequential quadratic programming algorithm for solution. The optimization variables of the optimization model are the velocity v4 required to leave the starting orbit and the velocity v5 required to enter the target orbit. The optimization target is the sum of the squares of the velocity increments required to leave the starting orbit and enter the target orbit. The constraint variables are the constraints of the orbital dynamics and the position constraints of the starting and target points. The mathematical description of the sequential quadratic programming algorithm is as follows:

[0150]

[0151] Where f(x) is the dynamics of the orbit, x c is the starting point of the starting track, x tis the starting point of the target orbit. The gravitational bodies added to the ephemeris model are: Earth, Sun, Moon, Venus, Mars, and Saturn. The three-segment transfer orbit, i.e., the stable manifold and unstable manifold of the Sun and the Earth, and the stable manifold of the Earth and the Moon, are corrected as a whole. The sequential quadratic programming algorithm is used to solve the optimization problem. The optimization variables are the three-dimensional components of the three-segment orbit entry pulses. The low-energy transfer orbit under the high-precision ephemeris model is obtained, i.e., the Earth-Moon orbit under the four-pulse ephemeris model based on manifold splicing. The orbital data is a 6×N3 array, such as Fig. 9 shown.

[0152] Another embodiment of the present invention provides a system for optimizing the Earth-Moon orbit based on a four-pulse ephemeris model with manifold splicing, the system comprising:

[0153] The optimal manifold splicing point calculation module is configured to process the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system based on the invariant manifold, and obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; the optimal manifold splicing points of the Halo orbits include the optimal position and the optimal speed of the spacecraft; the specific process includes: calculating and obtaining the initial value array of the Halo orbits of the Sun-Earth system or the Earth-Moon system; the initial value array includes the position and speed of multiple discrete points during the flight of the spacecraft; using a numerical method to solve the initial value array of the Halo orbit to obtain the invariant manifold of the Halo orbit; the invariant manifold includes an unstable manifold and a stable manifold; using a point selection method to splice the invariant manifolds of the Halo orbit, and calculating and obtaining the optimal manifold splicing points of the comprehensive position deviation and speed deviation;

[0154] The optimal splicing orbit acquisition module is configured to construct a four-pulse dual-circle restricted four-body model; take the optimal manifold splicing points of the Halo orbits corresponding to the sun-earth system and the earth-moon system as initial values, use a particle swarm optimization algorithm and a sequential quadratic programming algorithm to solve the dual-circle restricted four-body model based on the four-pulse to obtain the optimal splicing orbit; specifically comprising: taking the optimal manifold splicing points of the Halo orbits corresponding to the sun-earth system and the earth-moon system as initial values, using a particle swarm optimization algorithm to solve the upper optimization model, and outputting the optimal perturbation factors of the Halo orbits corresponding to the sun-earth system and the earth-moon system; wherein, in the particle swarm optimization algorithm, During the iterative optimization process of the method, a perturbation factor threshold is set. At the end of each iteration, it is determined whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, the optimization result is punished, that is, its performance index is increased. The optimal manifold splicing point of the Halo orbit corresponding to the sun-earth system and the earth-moon system, as well as the optimal perturbation factor of the Halo orbit corresponding to the sun-earth system and the earth-moon system are taken as initial values, and the sequential quadratic programming algorithm is used to solve the lower optimization model, and the time taken for the spacecraft in the Halo orbit of the sun-earth system and the earth-moon system to run from the final corrected initial value to the optimal manifold splicing point is output.

[0155] The functions of the Earth-Moon orbit optimization design system under the four-pulse ephemeris model based on manifold splicing in an embodiment of the present invention can be explained by the aforementioned Earth-Moon orbit optimization design method under the four-pulse ephemeris model based on manifold splicing. Therefore, for the parts not described in detail in the system embodiment, please refer to the above method embodiment and will not be repeated here.

[0156] Although the present invention has been described according to a limited number of embodiments, it will be apparent to those skilled in the art, with the benefit of the above description, that other embodiments are contemplated within the scope of the invention thus described. The disclosure of the present invention is intended to be illustrative rather than restrictive of the scope of the invention, which is defined by the appended claims.

[0157] The documents cited in the present invention are as follows:

[0158] [1] Lei Hanlun. Libration points, invariant currents and low-energy orbits[D]. Nanjing: Nanjing University, 2015.

[0159] [2] Li Mingtao. Design and optimization of energy-saving orbits for collinear libration missions[D]. Beijing: Chinese Academy of Sciences, 2010.

Claims

1. A method for optimizing the Earth-Moon orbit based on a four-pulse ephemeris model with manifold splicing, characterized in that: The following steps are involved: Step 1: Based on the invariant manifold, the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are processed to obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; the optimal manifold splicing points of the Halo orbits include the optimal position and optimal speed of the spacecraft; Step 2: Construct a two-circle restrictive four-body model based on four pulses; take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, and use the particle swarm optimization algorithm and sequential quadratic programming algorithm to solve the two-circle restrictive four-body model based on four pulses to obtain the optimal splicing orbit.

2. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 1 is characterized in that: In step 1, the process of obtaining the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system includes: Step 11: Calculate and obtain the initial value array of the Halo orbit of the Sun-Earth system or the Earth-Moon system; the initial value array includes the positions and velocities of multiple discrete points during the flight of the spacecraft; Step 12: using a numerical method to solve the initial value array of the Halo orbit to obtain an invariant manifold of the Halo orbit; the invariant manifold includes an unstable manifold and a stable manifold; Step 1: Use the point selection method to splice the invariant manifold of the Halo orbit and calculate the optimal manifold splicing point for the comprehensive position deviation and velocity deviation.

3. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 2 is characterized in that: The specific process of step one includes: Under the circular restricted three-body dynamics model, the Halo orbit starts from the xz plane, and its initial value is calculated by Richardson's third-order analytical solution; The initial value is subjected to a cycle of dynamic integration to obtain the final state; the error between the final state and the target state is calculated. If the error is greater than the given error, the final state is substituted into the differential correction solution equation to obtain the initial value of the next iteration, i.e. the corrected initial value; Repeat the above process until the error between the final state and the target state is less than the given error, and then terminate the iteration to obtain the final corrected initial value. The final corrected initial value is subjected to a period of dynamic integration to obtain the Halo orbit initial value array.

4. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 3 is characterized in that: The differential correction solution equation described in step 1 is expressed as: In the formula, is the state transfer matrix, which represents the linear estimation relationship between the initial state and the final state of the system; x0 represents the position of the spacecraft on the X axis; represents the speed of the spacecraft on the Y axis; Indicates the velocity of the spacecraft on the X axis at the final state; represents the acceleration of the spacecraft on the X axis at the final state; f Indicates the position of the spacecraft on the Y axis at the final state; Indicates the velocity of the spacecraft on the Y axis at the final state; Indicates the velocity of the spacecraft on the Z axis at the final state; Indicates the acceleration of the spacecraft on the Z axis at the final state; δ represents the correction amount, and the variable after it represents the correction amount of the variable at each iteration.

5. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 4 is characterized in that: The specific process of steps one and two includes: The state transfer matrix corresponding to the final state obtained by the dynamic integration of the final corrected initial value after one cycle is expressed as a single-value matrix, and the single-value matrix is ​​solved to obtain the eigenvalues ​​and eigenvectors; the eigenvalues ​​include: λ1>1 and The initial value point of the unstable manifold is calculated by the following formula: In the formula, p i represents the column data in the Halo orbit initial value array; ε represents the perturbation factor; D u represents the eigenvector corresponding to the eigenvalue λ1; The initial value point of the integration of the stable manifold is calculated by the following formula: Where D s represents the eigenvector corresponding to the eigenvalue λ2; After obtaining the initial value points of the integration of the two manifolds, the circular restricted three-body dynamics model is used for recursion to obtain the invariant manifold array of the Halo orbit; the invariant manifold array includes an unstable manifold array and a stable manifold array.

6. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 5 is characterized in that: The specific process of steps one and three includes: Calculate the minimum position deviation Δr between the position in the unstable manifold array and the position in the stable manifold array i ; Calculate the minimum speed deviation Δv between the speed in the unstable manifold array and the speed in the stable manifold array i ; For the minimum position deviation Δr i Deviation from minimum speed Δv i Normalize and perform weighted summation to obtain the comprehensive deviation right Sort and select the point corresponding to the smallest comprehensive deviation, which is the optimal manifold splicing point of the comprehensive position deviation and velocity deviation.

7. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 6 is characterized in that: Step 4: The four-pulse dual-circle restrictive four-body model includes an upper optimization model and a lower optimization model; wherein the upper optimization model is expressed as: The lower layer optimization model is expressed as: In the formula, min means minimization; and They represent the optimal speeds at the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; ε s represents the disturbance factor corresponding to the Sun-Earth system, ε m represents the disturbance factor corresponding to the Earth-Moon system; Indicates the preset disturbance factor threshold; and They respectively represent the optimal positions of the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; t1 represents the time it takes for the spacecraft in the Halo orbit of the Sun-Earth system to run from the final corrected initial value to the optimal manifold splicing point, and t2 represents the time it takes for the spacecraft in the Halo orbit of the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

8. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 7 is characterized in that: In step 2, the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system are used as initial values, and the particle swarm optimization algorithm and the sequential quadratic programming algorithm are used to solve the four-pulse based double-circle restricted four-body model to obtain the optimal splicing orbit, which includes: Step 21: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as the initial value, use the particle swarm optimization algorithm to solve the upper optimization model, and output the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system; wherein, in the iterative optimization process of the particle swarm optimization algorithm, set a perturbation factor threshold, and at the end of each iteration, determine whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, punish the optimization result, that is, increase its performance index; Step 22: Take the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system, as well as the optimal perturbation factors of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system as initial values, use the sequential quadratic programming algorithm to solve the lower-level optimization model, and output the time it takes for the spacecraft in the Halo orbits of the Sun-Earth system and the Earth-Moon system to run from the final corrected initial value to the optimal manifold splicing point.

9. The method for optimizing the Earth-Moon orbit based on the four-pulse ephemeris model based on manifold splicing according to claim 8 is characterized in that: The method further includes step three: after obtaining the optimal spliced ​​orbit, correcting the position and state of the data points in the optimal spliced ​​orbit based on the four-pulse ephemeris model to obtain the optimal spliced ​​orbit under the four-pulse ephemeris model.

10. The Earth-Moon orbit optimization design system based on the four-pulse ephemeris model based on manifold splicing is characterized by: include: The optimal manifold splicing point calculation module is configured to process the initial values ​​of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system based on the invariant manifold, and obtain the optimal manifold splicing points of the Halo orbits corresponding to the Sun-Earth system and the Earth-Moon system respectively; the optimal manifold splicing points of the Halo orbits include the optimal position and the optimal speed of the spacecraft; the specific process includes: calculating and obtaining the initial value array of the Halo orbits of the Sun-Earth system or the Earth-Moon system; the initial value array includes the position and speed of multiple discrete points during the flight of the spacecraft; using a numerical method to solve the initial value array of the Halo orbit to obtain the invariant manifold of the Halo orbit; the invariant manifold includes an unstable manifold and a stable manifold; using a point selection method to splice the invariant manifolds of the Halo orbit, and calculating and obtaining the optimal manifold splicing points of the comprehensive position deviation and speed deviation; The optimal splicing orbit acquisition module is configured to construct a four-pulse dual-circle restricted four-body model; take the optimal manifold splicing points of the Halo orbits corresponding to the sun-earth system and the earth-moon system as initial values, use a particle swarm optimization algorithm and a sequential quadratic programming algorithm to solve the dual-circle restricted four-body model based on the four-pulse to obtain the optimal splicing orbit; specifically comprising: taking the optimal manifold splicing points of the Halo orbits corresponding to the sun-earth system and the earth-moon system as initial values, using a particle swarm optimization algorithm to solve the upper optimization model, and outputting the optimal perturbation factors of the Halo orbits corresponding to the sun-earth system and the earth-moon system; wherein, in the particle swarm optimization algorithm, During the iterative optimization process of the method, a perturbation factor threshold is set. At the end of each iteration, it is determined whether the local optimal perturbation factor obtained each time is less than or equal to the perturbation factor threshold. If it is greater, the optimization result is punished, that is, its performance index is increased. The optimal manifold splicing point of the Halo orbit corresponding to the sun-earth system and the earth-moon system, as well as the optimal perturbation factor of the Halo orbit corresponding to the sun-earth system and the earth-moon system are taken as initial values, and the sequential quadratic programming algorithm is used to solve the lower optimization model, and the time taken for the spacecraft in the Halo orbit of the sun-earth system and the earth-moon system to run from the final corrected initial value to the optimal manifold splicing point is output.

Images