Zero initial value guess earth-moon Halo orbit low-thrust intersection trajectory optimization method

Abstract

The invention discloses a zero-initial-value guess earth-moon Halo orbit low-thrust rendezvous trajectory optimization method, belongs to the field of low-thrust rendezvous trajectory optimization, and solves the problems that an optimal result cannot be ensured, the calculation time is long, the robustness is poor and the initial value sensitivity is poor. Comprising the steps of dispersing an original control trajectory and a state trajectory to obtain a target function meeting rendezvous task constraints; solving an intermediate variable of the target function, and calculating a discrete costate quantity through a costate mapping equation; in combination with the Pontryagin minimum principle, optimal control of a tracking spacecraft is obtained through Hamiltonian, and a two-point boundary value model is constructed; taking the discrete co-state quantity as an initial guess of the two-point boundary value model, and optimizing the discrete co-state quantity by using an SQP algorithm to obtain an output result; the kinetic equation and the co-state equation are positively integrated through the output result, and the optimal control rate and the small-thrust rendezvous trajectory of the tracking spacecraft are obtained; the method overcomes the problem of initial value sensitivity, shortens the calculation time, and ensures the optimal result and the robust calculation process.

Description

Technical Field

[0001] The present invention belongs to the technical field of small-thrust rendezvous trajectory optimization for the Earth-Moon Halo orbit, and particularly relates to a method for optimizing the small-thrust rendezvous trajectory of the Earth-Moon Halo orbit without initial value guessing. Background Art

[0002] With the exploration of space by humans, the Earth-Moon space has attracted extensive attention from scholars. Currently, the typical orbits existing in the Earth-Moon space are considered important transfer stations for human exploration of the Moon and entry into deep space. China has launched the Queqiao relay satellite, which operates on the Halo orbit at the Earth-Moon L2 point. This satellite has participated in the exploration of the far side of the Moon and serves as a relay communication satellite to maintain the communication link between the tracking station and the lunar rover. The Lunar GateWay planned to be built by the United States focuses on the 9:2-L2 southern family near-rectilinear halo orbit (NRHO) near the Moon. It aims to build a space station to support future deep-space exploration missions. In order to maintain the long-term operation of spacecraft, it is necessary to develop rendezvous and docking (RVD) technology on typical orbits in the Earth-Moon space. The rendezvous and docking mission in the Earth-Moon space was first proposed in the Apollo program, which considered a linear model of the circular restricted three-body problem. Subsequent studies have used relative dynamics models, such as the elliptical linear relative motion equation (ELERM) and circular linear relative motion equation (CLERM) constructed in the local vertical-local horizontal (LVLH) coordinate system centered on the Moon. After that, scholars have proposed strategies for rendezvous using impulses and continuous thrusts.

[0003] In recent years, compared with traditional chemical fuel propulsion, small-thrust propulsion technologies with higher specific impulse and smaller thrust have been developed to improve the utilization efficiency of propellants. As shown in Table 1 of typical propulsion technologies and performances, examples of small-thrust propulsion include electric propulsion and solar sail propulsion suitable for deep-space missions. Electric propulsion converts electrical energy into kinetic energy through electrothermal, electrostatic, or electromagnetic means. Solar sails use lightweight, highly reflective films to reflect photons and generate thrust from solar radiation pressure. For rendezvous problems, electric propulsion is more suitable for close-range missions, while solar sails have not been successfully applied due to their lower thrust and more complex sail control (including control of the sail and support mechanism). These new propulsion technologies correspond to small-thrust control strategies, such as linear quadratic regulator (LQR) control.

[0004] Table 1

[0005]

[0006]

[0007] Compared with the optimization problems formed by traditional impulsive maneuvers, the low-thrust propulsion mode increases the complexity of the optimization problem due to the higher dimension of the optimization domain. Generally speaking, the low-thrust rendezvous trajectory optimization problem is usually solved by direct methods or indirect methods;

[0008] For example, the publication number: CN112084581A, the invention name: A method and system for optimizing the low-thrust perturbed rendezvous trajectory of a spacecraft; this technical solution calculates the low-thrust rendezvous pulse velocity increment according to the given initial orbital position and target orbital position of the spacecraft, and then substitutes it as the initial value and the optimality criterion of the solution into the indirect method optimization model to obtain the optimal low-thrust transfer trajectory; although this method solves to a certain extent the problem of large dynamic recurrence calculations and low calculation efficiency when the transfer time is long in the traditional indirect method low-thrust rendezvous optimization algorithm, it still faces the problems of non-optimal solution set and sensitivity of the co-state quantity at the initial moment.

[0009] Although some techniques have been proposed in the prior art to reduce the sensitivity of the co-state quantity at the initial moment, such as the homotopy method and the continuation method, most researchers prefer a heuristic hybrid method. However, the most obvious problem is its robustness, that is, different solutions will be generated under the same initial condition settings.

[0010] In summary, the purpose of low-thrust rendezvous trajectory optimization is to quickly solve the optimal solution under the given scenario, and obtain the optimal control rate of the tracking spacecraft and the low-thrust rendezvous trajectory of the Earth-Moon Halo orbit. The direct method does not consider the necessary conditions and cannot guarantee the optimality of the result; the indirect method requires a guessed value of the co-state quantity and has the problem of initial value sensitivity; the heuristic hybrid method has a long calculation time and poor robustness. Summary of the Invention

[0011] To solve the above technical problems, the present invention proposes a method for optimizing the low-thrust rendezvous trajectory of the Earth-Moon Halo orbit without initial value guessing. This method is a method for calculating and optimizing the continuous low-thrust spacecraft / aircraft rendezvous trajectory design problem without initial value guessing. Comparing with the defects of traditional methods, the initial value sensitivity problem will make it difficult to obtain a convergent solution that meets the constraints in the calculation process, resulting in the spacecraft being unable to complete the rendezvous mission. Using traditional heuristic methods requires a large amount of computing resources (from a few minutes to dozens of minutes), and it is impossible to achieve fast calculation, thus missing a good mission window. The non-fuel optimal control result will cause the spacecraft to waste unnecessary fuel in actual control, thereby reducing the on-orbit life of the spacecraft. The non-robust calculation process will make the algorithm only applicable to specific rendezvous scenarios / tasks and is difficult to be competent for the calculation of the optimal orbital control rate of the spacecraft in multi-scenarios / tasks. The purpose of the present invention is to overcome the initial value sensitivity problem of the indirect method, shorten the calculation time, ensure the optimality of the result, and make the calculation process robust.

[0012] The object of the present invention is specifically achieved through the following technical solutions:

[0013] The present invention discloses a small-thrust rendezvous trajectory optimization method for a lunar Halo orbit with zero initial value guessing, and the method includes:

[0014] Step 1: After discretizing the received original small-thrust rendezvous control trajectory and state trajectory of the lunar Halo orbit, determine the control quantity and state quantity at the discrete Legendre-Gauss-Radau collocation points, and obtain an objective function that satisfies the rendezvous mission constraints;

[0015] Step 2: Use the open-source nonlinear programming solver IPOPT to solve the intermediate variables of the objective function;

[0016] Step 3: Take the intermediate variables as inputs, and calculate the discrete co-state quantities through the co-state mapping equation;

[0017] Step 4: Combining with the Pontryagin minimum principle, obtain the optimal control of the chasing spacecraft through the Hamiltonian, and form a two-point boundary value model from the optimal control of the chasing spacecraft, the dynamic equation of the chasing spacecraft, the co-state equation, and the terminal constraints;

[0018] Step 5: Take the discrete co-state quantities as the initial guess of the co-state quantities at the initial moment of the two-point boundary value model, use the SQP algorithm to optimize the co-state quantities at the initial moment, obtain the output result of the two-point boundary value model, and obtain the minimum value of the objective function through the output result;

[0019] Step 6: Use the output result to integrate the dynamic equation and co-state equation of the chasing spacecraft forward, and obtain the optimal control rate of the chasing spacecraft and the small-thrust rendezvous trajectory of the lunar Halo orbit.

[0020] In Step 1, the Legendre-Gauss-Radau collocation points refer to the roots of the P N (τ)+P N-1 polynomial;

[0021] In the formula, P N (τ) is the Nth-order Legendre polynomial, P N-1 (τ) is the (N - 1)th-order Legendre polynomial, where P N (τ) is calculated by the following formula:

[0022]

[0023] Replace N in P N (τ) with N - 1 to obtain the expression of P N-1 (τ), where τ is the independent variable, d is the differential symbol, and N is the specified order.

[0024] In Step 1, the objective functions satisfying the rendezvous mission constraints include: the objective function of time-optimal control and the objective function of fuel-optimal control; among them,

[0025] (1) The objective function of time-optimal control is: Among them, time-optimal control means satisfying the control rate of the chaser spacecraft in the rendezvous mission to minimize the time required for rendezvous;

[0026] (2) The objective function of fuel-optimal control is: Among them, fuel-optimal control means satisfying the control rate of the chaser spacecraft in the rendezvous mission to minimize the fuel required for rendezvous;

[0027] In the formula, J 1 is the time required for rendezvous, t 0 is the initial time, t f is the terminal time, J 2 is the fuel required for rendezvous, and U is the control vector of the chaser spacecraft.

[0028] In Step 2, the intermediate variables include: Gaussian quadrature weights, Lagrange multiplier matrix, and differential matrix.

[0029] In Step 3, taking the intermediate variables as inputs, the discrete co-states calculated through the co-state mapping equation include: the discrete co-states of time-optimal control and the discrete co-states of fuel-optimal control; among them, the co-state mapping equation is:

[0030]

[0031] In the formula, the discrete co-states of time-optimal control and the discrete co-states of fuel-optimal control calculated through the co-state mapping equation are composed of the discrete co-states at the Legendre-Gauss-Radau collocation points and the discrete co-states at the endpoints. is the (M + 1)-th column of the transposed matrix D T of the differential matrix D, Λ k represents the k-th row of the Lagrange multiplier matrix Λ, ω k is the value of the Gaussian quadrature weight ω k at the k-th Legendre-Gauss-Radau collocation point, k is the number of the collocation point, and M is the total number of collocation points.

[0032] In Step 4, combining the Pontryagin's minimum principle, the optimal control of the chaser spacecraft obtained through the Hamiltonian includes: the time-optimal control U 1 of the chaser spacecraft obtained through the Hamiltonian H 1 of time-optimal control and the Hamiltonian H 2Obtain the optimal fuel control U of the tracking spacecraft 2 ; where

[0033] (1) The Hamiltonian H of the time-optimal control 1 is as follows:

[0034]

[0035] In the formula, Λ = [Λ r , Λ v , λ m represents the co-state quantity, Λ r is the position co-state quantity, Λ v is the velocity co-state quantity, λ m is the mass co-state quantity; v is the velocity of the tracking spacecraft, g(r) is the term related to the position r of the tracking spacecraft in the dynamic equation of the tracking spacecraft, h(v) is the term related to the velocity v of the tracking spacecraft in the dynamic equation of the tracking spacecraft, U is the control vector of the tracking spacecraft, T max is the maximum available thrust, m is the mass of the tracking spacecraft, u is the modulus of the control vector of the tracking spacecraft, I sp is the specific impulse, g 0 is the gravitational acceleration at sea level;

[0036] According to the Pontryagin minimum principle, the time-optimal control of the tracking spacecraft is obtained through the Hamiltonian H of the time-optimal control 1 as follows: In the formula, λ v is the modulus of the velocity co-state quantity;

[0037] (2) The Hamiltonian H of the fuel-optimal control 2 is as follows:

[0038]

[0039] According to the Pontryagin minimum principle, the fuel-optimal control of the tracking spacecraft is obtained through the Hamiltonian H of the fuel-optimal control 2 as follows: where, u 2 ∈[0, T max , representing the modulus of the fuel-optimal control vector of the tracking spacecraft.

[0040] In step four, the two-point boundary value model includes: the two-point boundary value model of the time-optimal control and the two-point boundary value model of the fuel-optimal control; where

[0041] The two-point boundary value model of the time-optimal control consists of the dynamic equation of the tracking spacecraft, the time-optimal control of the tracking spacecraft, the co-state equation of the time-optimal control, and the terminal constraint of the time-optimal control;

[0042] The two-point boundary value model of fuel optimal control consists of the dynamic equation of the chaser spacecraft, the fuel optimal control of the chaser spacecraft, the co-state equation of the fuel optimal control, and the terminal constraints of the fuel optimal control. Among them, the dynamic equation of the chaser spacecraft is:

[0043]

[0044] In the formula, is the derivative of the state quantity of the chaser spacecraft, f(X,U) is the dynamic equation obeyed by the state quantity of the chaser spacecraft, X = [r, v] is the state quantity of the chaser spacecraft, r is the position of the chaser spacecraft, v is the velocity of the chaser spacecraft, is the derivative of the position of the chaser spacecraft, is the derivative of the velocity of the chaser spacecraft, is the derivative of the mass of the chaser spacecraft, u = |U| represents the modulus value of the control vector U of the chaser spacecraft;

[0045] (1) The co-state equation of time optimal control is:

[0046]

[0047] In the formula, is the derivative of the co-state quantity of time optimal control, H 1 is the Hamiltonian of time optimal control;

[0048] The terminal constraints of time optimal control are:

[0049]

[0050] In the formula, φ 1 ([X, Λ 1 T , t 0 , t f ) is the terminal constraint value of time optimal control, Λ 1 is the co-state quantity of time optimal control of the chaser spacecraft, t 0 is the initial time, t f is the terminal time, r(t f ) is the position of the chaser spacecraft at the terminal time, r f is the position of the target spacecraft at the terminal time, v(t f ) is the velocity of the chaser spacecraft at the terminal time, v f is the velocity of the target spacecraft at the terminal time, λ m (t f ) is the mass co-state quantity of the chaser spacecraft at the terminal time, H 1 (t f ) is the Hamiltonian at the terminal time, 0 8 is an eight-dimensional zero vector;​

[0051] (2) The co - state equation of the fuel - optimal control is as follows:

[0052]

[0053] In the formula, is the derivative of the co - state quantity of the fuel - optimal control, and H 2 is the Hamiltonian of the fuel - optimal control;

[0054] The terminal constraint of the fuel - optimal control is:

[0055]

[0056] In the formula, φ 2 ([X, Λ 2 ) T , t 0 , t f ) is the terminal constraint value of the fuel - optimal control, Λ 2 is the co - state quantity of the fuel - optimal control of the tracking spacecraft, and 0 7 is the seven - dimensional zero vector.

[0057] In step 5, use the SQP algorithm to optimize the co - state quantity at the initial time, obtain the output result of the two - point boundary - value model, and obtain the minimum value of the objective function from the output result, including: the minimum value of the objective function of the time - optimal control and the minimum value of the objective function of the fuel - optimal control; among them,

[0058] (1) The calculation method of the minimum value of the objective function of the time - optimal control includes:

[0059] Use the SQP algorithm to optimize the discrete co - state quantity of the time - optimal control; stop the calculation when the control rate calculated by the two - point boundary - value model of the time - optimal control satisfies the terminal constraint of the time - optimal control, output the co - state quantity at the initial time of the time - optimal control and the terminal time of the time - optimal control, and obtain the minimum value of the objective function of the time - optimal control from the co - state quantity at the initial time of the time - optimal control and the terminal time of the time - optimal control;

[0060] (2) The calculation method of the minimum value of the objective function of the fuel - optimal control includes:

[0061] Construct the objective function J of the energy - optimal control that satisfies the rendezvous mission constraints;

[0062] Use the open - source nonlinear programming solver IPOPT to solve the intermediate variables of the objective function of the energy - optimal control;

[0063] Take the intermediate variables as the input and calculate the discrete co - state quantity of the energy - optimal control through the co - state mapping equation;

[0064] Combined with the Pontryagin minimum principle, the Hamiltonian H of the energy optimal control 3 obtains the energy optimal control U of the tracking spacecraft 3 ;

[0065] Taking the discrete co-state quantity of the energy optimal control as the initial guess of the co-state quantity at the initial moment of the two-point boundary value model of the fuel optimal control, using the SQP algorithm to optimize the co-state quantity at the initial moment, based on the homotopy iteration process, the Hamiltonian of the fuel optimal control and the fuel optimal control of the tracking spacecraft are obtained from the Hamiltonian of the energy optimal control and the energy optimal control of the tracking spacecraft; when the control rate calculated by the two-point boundary value model of the fuel optimal control satisfies the terminal constraint of the fuel optimal control, the calculation is stopped, and the co-state quantity at the initial moment of the fuel optimal control is output, and the minimum value of the objective function of the fuel optimal control is obtained from the co-state quantity at the initial moment of the fuel optimal control; among them,

[0066] The objective function of the energy optimal control that satisfies the rendezvous mission constraint is:

[0067] The Hamiltonian H of the energy optimal control 3 is:

[0068] The energy optimal control of the tracking spacecraft where, u 3 ∈[0, T max , representing the modulus value of the energy optimal control vector of the tracking spacecraft.

[0069] In step five, the iterative calculation method of the homotopy iteration process is:

[0070]

[0071] In the formula, is the homotopy parameter for each step; i = 10, 9,..., 1, corresponding to the homotopy parameter ∈ = 1 to ∈ = 0; among them, ∈ = 1 corresponds to the Hamiltonian of the energy optimal control, ∈ = 0 corresponds to the Hamiltonian of the fuel optimal control, and H 4 represents the Hamiltonian in the homotopy iteration process;

[0072] In each iteration step, the SQP algorithm is used to optimize the co-state quantity at the initial moment. After satisfying the terminal constraint of the fuel optimal control, the calculation stops and enters the next iteration step; when the homotopy parameter ∈ = 0, the Hamiltonian of the energy optimal control calculated is the same as the Hamiltonian of the fuel optimal control, and the Hamiltonian of the fuel optimal control is obtained.

[0073] In step six, the integration process of the forward integration uses a Runge-Kutta 7 / 8 order variable step size integrator for forward integration.

[0074] The beneficial effects of the present invention are:

[0075] 1. After discretizing the received original Earth-Moon Halo orbit low-thrust rendezvous control trajectory and state trajectory, the control quantities and state quantities at the discrete Legendre-Gauss-Radau collocation points are determined, and a technical solution for the objective function that meets the rendezvous mission constraints is obtained. The conversion from the original continuous low-thrust rendezvous problem to the nonlinear programming process expands the convergence domain of the original problem and reduces the difficulty of solving the continuous low-thrust rendezvous problem.

[0076] 2. The open-source nonlinear programming solver IPOPT is used to solve the intermediate variables of the objective function. Taking the intermediate variables as inputs, the discrete co-state quantities are calculated through the co-state mapping equation, and the discrete co-state quantities are used as the initial guess of the co-state quantities at the initial time of the two-point boundary value model; a reasonable initial value guess is provided for the co-state quantities at the initial time of the two-point boundary value model. Compared with the traditional indirect method of selecting random co-state quantities, the problem of initial value sensitivity is overcome, and the calculation speed is greatly improved compared with the heuristic method.

[0077] 3. The Pontryagin minimum principle is applied to ensure that the calculation result is the optimal solution.

[0078] 4. In the calculation of the homotopy iteration process, homotopy parameters with gradually decreasing modulus and gradually decreasing step size are used, which improves the robustness of the calculation process.

[0079] 5. The optimal control rate of the tracking spacecraft can be obtained by integrating the minimum value of the objective function forward.

[0080] 6. The technical solution disclosed in the present invention uses different convergence tolerances respectively, and greatly improves the calculation speed compared with the traditional indirect method and heuristic algorithm with a single convergence tolerance. BRIEF DESCRIPTION OF THE DRAWINGS

[0081] The present invention will be further described in detail below with reference to the drawings and embodiments.

[0082] Figure 1 It is a schematic diagram of the Earth-Moon rotating coordinate system provided by an embodiment of the present invention.

[0083] Figure 2 It is a schematic diagram of the rendezvous trajectory on the Halo orbit near L1 provided by Example 1 of the present invention.

[0084] Figure 3 It is a schematic diagram of the rendezvous trajectory on the Halo orbit near L2 provided by Example 1 of the present invention.

[0085] Figure 4 It is a schematic diagram of the rendezvous trajectory on the near-rectilinear halo orbit in the south of L2 provided by Example 1 of the present invention.

[0086] Figure 5It is a schematic diagram of the Halo orbit near L2 and the near-rectilinear halo orbit of the L2 southern family provided in Example 3 of the present invention.

[0087] Figure 6 It is a schematic diagram of the running time of three algorithms in the MT problem provided in Example 4 of the present invention.

[0088] Figure 7 It is a schematic diagram of the running time of three algorithms in the ME problem provided in Example 4 of the present invention.

[0089] Figure 8 It is a schematic diagram of the running time of three algorithms in the MF problem provided in Example 4 of the present invention.

[0090] Figure 9 It is a schematic diagram of the fuel consumption of three algorithms for calculating the MT problem provided in Example 4 of the present invention.

[0091] Figure 10 It is a schematic diagram of the fuel consumption of three algorithms for calculating the ME problem provided in Example 4 of the present invention.

[0092] Figure 11 It is a schematic diagram of the fuel consumption of three algorithms for calculating the MF problem provided in Example 4 of the present invention. Detailed implementation manners

[0093] The embodiment of the present invention provides a small-thrust rendezvous trajectory optimization method for the Earth-Moon Halo orbit without initial value guessing, and the method includes:

[0094] Step 1: After discretizing the received original small-thrust rendezvous control trajectory and state trajectory of the Earth-Moon Halo orbit, determine the control quantity and state quantity at the discrete Legendre-Gauss-Radau collocation points to obtain an objective function that satisfies the rendezvous mission constraints;

[0095] Step 2: Use the open-source nonlinear programming solver IPOPT to solve the intermediate variables of the objective function;

[0096] Among them, the open-source nonlinear programming solver IPOPT can be downloaded from the corresponding website, and the preferred download address is: https: / / coin-or.github.io / Ipopt / .

[0097] Step 3: Take the intermediate variables as inputs and calculate the discrete adjoint variables through the adjoint mapping equation;

[0098] Step 4: Combine the Pontryagin minimum principle, obtain the optimal control of the tracking spacecraft through the Hamiltonian, and form a two-point boundary value model from the optimal control of the tracking spacecraft, the dynamic equation of the tracking spacecraft, the adjoint equation, and the terminal constraints;

[0099] The Pontryagin minimum principle is a principle used to determine the optimal control rate, which is specifically stated as follows: the value of the optimal control rate should be such that the system Hamiltonian is minimized. Among them, the system Hamiltonian is determined by the dynamic equation and the objective function.

[0100] Step 5: Use the discrete co-state quantity as the initial guess of the co-state quantity at the initial time of the two-point boundary value model, optimize the co-state quantity at the initial time using the SQP algorithm, obtain the output result of the two-point boundary value model, and obtain the minimum value of the objective function through the output result;

[0101] Step 6: Use the output result to integrate forward to track the spacecraft dynamic equation and the co-state equation, and obtain the optimal control rate of the tracking spacecraft and the low-thrust rendezvous trajectory of the Earth-Moon Halo orbit.

[0102] In Step 1, the Legendre-Gauss-Radau collocation points refer to the roots of the P N (τ) + P N-1 (τ) polynomials;

[0103] In the formula, P N (τ) is the Nth-order Legendre polynomial, and P N-1 (τ) is the (N - 1)th-order Legendre polynomial. Among them, P N (τ) is calculated by the following formula:

[0104]

[0105] Replace N in P N (τ) with N - 1 to obtain the expression of P N-1 (τ). τ is the independent variable, d is the differential symbol, and N is the specified order.

[0106] In Step 1, the objective functions that satisfy the rendezvous mission constraints include: the objective function of time-optimal control and the objective function of fuel-optimal control; among them,

[0107] (1) The objective function of time-optimal control is: Among them, time-optimal control refers to the control rate of the tracking spacecraft that satisfies the rendezvous mission, making the time required for rendezvous the shortest;

[0108] (2) The objective function of fuel-optimal control is: Among them, fuel-optimal control refers to the control rate of the tracking spacecraft that satisfies the rendezvous mission, making the fuel required for rendezvous the least;

[0109] In the formula, J 1 is the time required for rendezvous, t 0 is the initial time, t f is the terminal time, and J 2For the fuel required for rendezvous, U is the control vector of the chaser spacecraft.

[0110] In step two, the intermediate variables include: Gauss quadrature weights, Lagrange multiplier matrix, and differential matrix.

[0111] In step three, taking the intermediate variables as inputs, the discrete co-states calculated by the co-state mapping equation include: the discrete co-states of time-optimal control and the discrete co-states of fuel-optimal control; where the co-state mapping equation is:

[0112]

[0113] In the formula, the discrete co-states of time-optimal control and the discrete co-states of fuel-optimal control calculated by the co-state mapping equation are composed of the discrete co-states at the Legendre-Gauss-Radau collocation points and the discrete co-states at the endpoints, is the (M + 1)-th column of the transpose matrix D T of the differential matrix D, Λ k represents the k-th row of the Lagrange multiplier matrix Λ, ω k is the value of the Gauss quadrature weight ω k at the k-th Legendre-Gauss-Radau collocation point, k is the number of the collocation point, and M is the total number of collocation points.

[0114] In step four, combining the Pontryagin minimum principle, the optimal control of the chaser spacecraft obtained through the Hamiltonian includes: the time-optimal control U 1 of the chaser spacecraft obtained through the Hamiltonian H 1 of time-optimal control and the fuel-optimal control U 2 of the chaser spacecraft obtained through the Hamiltonian H 2 of fuel-optimal control; where,

[0115] (1) The Hamiltonian H 1 of time-optimal control is:

[0116]

[0117] In the formula, Λ = [Λ r , Λ v , λ m represents the co-state, Λ r is the position co-state, Λ v is the velocity co-state, λ mis the mass co-state quantity; v is the velocity of the chasing spacecraft, g(r) is the term related to the position r of the chasing spacecraft in the dynamic equation of the chasing spacecraft, h(v) is the term related to the velocity v of the chasing spacecraft in the dynamic equation of the chasing spacecraft, U is the control vector of the chasing spacecraft, T max is the maximum available thrust, m is the mass of the chasing spacecraft, u is the magnitude of the control vector of the chasing spacecraft, I sp is the specific impulse, g 0 is the gravitational acceleration at sea level. r = [r x , r y , r z , v = [v x , v y , v z ;

[0118]

[0119] Among them, and respectively represent the distances between the chasing spacecraft and the Earth and the Moon. μ = 0.01215 is the normalized value of the lunar mass. Among them, the definition of the Earth-Moon rotating coordinate system is: the x-axis of this coordinate system points from the Earth to the Moon, the z-axis points along the angular momentum direction of the lunar orbit, the y-axis is determined according to the orthogonal rule, and the center of the coordinate system is the position of the Earth-Moon barycenter. The schematic diagram of the Earth-Moon rotating coordinate system is as Figure 1 shown.

[0120] As Figure 1 shown, the motion of the chasing spacecraft in the Earth-Moon space can be approximated as a circular restricted three-body problem (CR3BP). The Earth and the Moon move in a circular plane around their common barycenter. The motion of the chasing spacecraft is affected by the combined gravitational force of the two, and the mass of the chasing spacecraft is negligible.

[0121] According to the Pontryagin minimum principle, the time-optimal control Hamiltonian H 1 of the chasing spacecraft is obtained as follows: In the formula, λ b is the magnitude of the velocity co-state quantity;

[0122] (2) The fuel-optimal control Hamiltonian H 2 is:

[0123]

[0124] According to the Pontryagin minimum principle, the fuel-optimal control Hamiltonian H 2 of the chasing spacecraft is obtained as follows: Among them, u 2 ∈[0, T max, which represents the magnitude of the optimal control vector of the tracking spacecraft fuel.

[0125] In step 4, the two-point boundary value model includes: the two-point boundary value model of time-optimal control and the two-point boundary value model of fuel-optimal control. The known parameters of the model are the state quantities at the initial time and the state quantities at the terminal time, and the parameters to be solved are the co-state quantities at the initial time; the state quantities at the initial time are composed of the state quantities of the tracking spacecraft and the target spacecraft at the start time of the rendezvous, and the state quantities at the terminal time are composed of the state quantities of the tracking spacecraft and the target spacecraft at the end time of the rendezvous. Due to the rendezvous constraints, the state quantities of the tracking spacecraft and the target spacecraft are equal at the terminal time. The present invention uses r f = [x, y, z] and represents the state quantity of the target spacecraft at the terminal time; where,

[0126] The two-point boundary value model of time-optimal control is composed of the dynamic equation of the tracking spacecraft, the time-optimal control of the tracking spacecraft, the co-state equation of time-optimal control, and the terminal constraint of time-optimal control;

[0127] The two-point boundary value model of fuel-optimal control is composed of the dynamic equation of the tracking spacecraft, the fuel-optimal control of the tracking spacecraft, the co-state equation of fuel-optimal control, and the terminal constraint of fuel-optimal control; where, the dynamic equation of the tracking spacecraft is:

[0128]

[0129] In the formula, is the derivative of the state quantity of the tracking spacecraft, f(X,U) is the dynamic equation obeyed by the state quantity of the tracking spacecraft, X = [r, v] is the state quantity of the tracking spacecraft, r is the position of the tracking spacecraft, v is the velocity of the tracking spacecraft, is the derivative of the position of the tracking spacecraft, is the derivative of the velocity of the tracking spacecraft, is the derivative of the mass of the tracking spacecraft, u = |U| represents the magnitude of the control vector U of the tracking spacecraft;

[0130] (1) The co-state equation of time-optimal control is:

[0131]

[0132] In the formula, is the derivative of the co-state quantity of time-optimal control, H 1 is the Hamiltonian of time-optimal control;

[0133] The terminal constraint of time-optimal control is:

[0134]

[0135] where φ 1 ([X, Λ 1 T , t 0 , t f ) is the terminal constraint value of the time-optimal control, Λ 1 is the co-state quantity of the time-optimal control of the chaser spacecraft, t 0 is the initial time, t f is the terminal time, r(t f ) is the position of the chaser spacecraft at the terminal time, r f is the position of the target spacecraft at the terminal time, v(t f ) is the velocity of the chaser spacecraft at the terminal time, v f is the velocity of the target spacecraft at the terminal time, λ m (t f ) is the co-state quantity of the mass of the chaser spacecraft at the terminal time, H 1 (t f ) is the Hamiltonian at the terminal time, 0 8 is the eight-dimensional zero vector;

[0136] (2) The co-state equation of the fuel-optimal control is:

[0137]

[0138] where is the derivative of the co-state quantity of the fuel-optimal control, H 2 is the Hamiltonian of the fuel-optimal control; X = [r, v] is the state quantity of the chaser spacecraft;

[0139] (2) The co-state equation of the fuel-optimal control is:

[0140]

[0141] where is the derivative of the co-state quantity of the fuel-optimal control, H 2 is the Hamiltonian of the fuel-optimal control;

[0142] The terminal constraint of the fuel-optimal control is:

[0143]

[0144] where φ 2 ([X, Λ 2 T , t 0 , t f ) is the terminal constraint value of the fuel-optimal control, Λ 2 is the co-state quantity of the fuel-optimal control of the chaser spacecraft, 0 7 is the seven-dimensional zero vector. ​​

[0145] In Step 5, the SQP algorithm is used to optimize the co - state variables at the initial moment, and the output result of the two - point boundary - value model is obtained. The minimum value of the objective function obtained from the output result includes: the minimum value of the objective function of time - optimal control and the minimum value of the objective function of fuel - optimal control. Among them,

[0146] (1) In time - optimal control, the terminal time is free. Therefore, the variables to be solved in the two - point boundary - value model of time - optimal control are [Λ, t f , a total of 8 dimensions, where t f is the terminal time; the calculation method of the minimum value of the objective function of time - optimal control includes:

[0147] Use the SQP algorithm to optimize the discrete co - state variables of time - optimal control; stop the calculation when the control rate calculated by the two - point boundary - value model of time - optimal control satisfies the terminal constraints of time - optimal control, and output the co - state variables at the initial moment and the terminal time of time - optimal control. The minimum value of the objective function of time - optimal control is obtained from the co - state variables at the initial moment and the terminal time of time - optimal control;

[0148] (2) In fuel - optimal control, the terminal time is fixed. Therefore, the variable to be solved in the two - point boundary - value model of fuel - optimal control is Λ with 7 dimensions. Since it is difficult to directly solve fuel - optimal control, first construct and solve the minimum - energy (ME) control, and then solve fuel - optimal control through a homotopy iteration process. The calculation method of the minimum value of the objective function of fuel - optimal control includes:

[0149] Construct the objective function of the minimum - energy (ME) control that satisfies the rendezvous mission constraints:

[0150] Use the open - source nonlinear programming solver IPOPT to solve the intermediate variables of the objective function J of the minimum - energy control;

[0151] Take the intermediate variables: Gaussian quadrature weights, differential matrices, and Lagrange multipliers as inputs, and calculate the discrete co - state variables of the minimum - energy control through the co - state mapping equation;

[0152] Combined with the Pontryagin minimum principle, through the Hamiltonian H 3 : Obtain the minimum - energy control of the chaser spacecraft where u 3 ∈[0, T max , representing the modulus of the minimum - energy control vector of the chaser spacecraft;

[0153] Take the discrete co - state quantity of the energy - optimal control as the initial guess of the co - state quantity at the initial moment of the two - point boundary - value model of the fuel - optimal control. Use the SQP algorithm to optimize the co - state quantity at the initial moment. Based on the homotopy iteration process, stop the calculation when the control rate calculated by the two - point boundary - value model of the fuel - optimal control satisfies the terminal constraint of the fuel - optimal control, output the co - state quantity at the initial moment of the fuel - optimal control, and obtain the minimum value of the objective function of the fuel - optimal control from the co - state quantity at the initial moment of the fuel - optimal control.

[0154] In step five, the iterative calculation method of the homotopy iteration process is as follows:

[0155]

[0156] In the formula, is the homotopy parameter for each step; i = 10, 9, …, 1, corresponding to the homotopy parameter ∈ from ∈ = 1 to ∈ = 0; where ∈ = 1 corresponds to the Hamiltonian of the energy - optimal control, ∈ = 0 corresponds to the Hamiltonian of the fuel - optimal control, and H 4 represents the Hamiltonian in the homotopy iteration process;

[0157] In each iteration step, the SQP algorithm is used to optimize the initial - moment co - state quantity Λ. Stop after satisfying the terminal constraint of the fuel - optimal control and enter the next iteration step; when the homotopy parameter ∈ = 0, the calculated Hamiltonian of the energy - optimal control is the same as the Hamiltonian of the fuel - optimal control, and the Hamiltonian of the fuel - optimal control is obtained. The corresponding result at this time is the fuel - optimal control, that is, the output result of the two - point boundary - value model of the fuel - optimal control. This result makes the corresponding objective function of the fuel - optimal control take the minimum value. This iterative process corresponds to the transition from the energy - optimal control to the fuel - optimal control and is called the homotopy iteration process.

[0158] In step six, the output results of the two - point boundary - value model of the time - optimal control are the co - state quantity at the initial moment of the time - optimal control and the terminal moment of the time - optimal control, and the output result of the two - point boundary - value model of the fuel - optimal control is the co - state quantity at the initial moment of the fuel - optimal control. The forward - integration process uses a Runge - Kutta 7 / 8 - order variable - step - size integrator for forward integration.

[0159] To enable those skilled in the art to understand the technical solution of the present invention, five specific examples are provided below to elaborate in detail on the technical solution ZIGA disclosed in the present invention:

[0160] ZIGA is applied to the continuous low - thrust rendezvous problem of the Earth - Moon Halo orbit

[0161] The zero-initial-value guessing algorithm is used to solve three typical low-thrust rendezvous problems and compare the results with those of traditional two-impulse maneuver rendezvous (referred to as the three-body Lambert problem (TBLP)). The initial orbits of the target spacecraft and the chaser spacecraft are selected from three typical Earth-Moon orbits.

[0162] Example 1: Halo orbit near L1

[0163] The Halo orbit around L1 is usually selected for scientific missions. In the rendezvous mission, the target and chaser spacecraft are located on the L1 southern family Halo orbit, with an amplitude of 5000 km along the z-axis and a period of 11.99 days. It is assumed that the chaser spacecraft is at the 0-degree (r y = 0) phase at the start, while the target is at the 12.5-degree phase. In the ME and MF problems, the rendezvous time is assumed to be 0.25 of a Halo orbit period. In step two of ZIGA, the convergence tolerance is set to 10 -5 for quickly solving the MT and ME problems, and in step five, the tolerance is set to 10 -15 for calculating high-precision results. In addition, T max = 30 mN, I sp = 3000 s and m = 1500 kg.

[0164] In Figure 2 , Figure 3 and Figure 4 , I represents the initial position, F represents the final position, O represents the rendezvous orbit, the subscript C represents the chaser spacecraft, and T represents the target spacecraft. The words in parentheses represent the corresponding rendezvous problems. Since the final time is fixed, the final positions of the ME problem, MF problem, and TBLP coincide. However, the rendezvous position of the MT problem is ahead of their phases on the orbit, which means less time is consumed for the rendezvous.

[0165] Example 2: Halo orbit near L2

[0166] Queqiao satellite is the first communication relay satellite placed on the Earth-Moon L2-Halo orbit, connecting the Earth and the far side of the Moon. This example uses the L2 northern family Halo orbit, with an amplitude of 5000 km along the z-axis and a period of 14.91 days. The rendezvous time is 0.25 of an orbit period. It is assumed that the chaser spacecraft is at the apogee (r y = 0), and the phase of the target spacecraft is 12.5 degrees ahead of the chaser spacecraft.

[0167] To solve the rendezvous problem on the Halo orbit, it is usually necessary to solve for λ(t 0 ) and t fSelect a suitable value. In the MT problem, it is well known that the engine will always use the maximum available thrust. The introduction of this prior information greatly reduces the dimension of the NLP problem, thus facilitating its rapid convergence with arbitrary precision. However, this comes at the cost of an increase in CPU running time. In the ME problem, low-order numerical integration methods can be used to calculate the dynamic equations. However, during the RPM calculation process, the thrust curve must be approximated by interpolation, which results in a slightly increased solution time for the ME problem compared to the MT problem. The MF problem is the most difficult to solve and requires determining the appropriate rendezvous time and reasonable initial values. To accelerate the solution of the MF problem, a larger tolerance is used in Step 2 of ZIGA, which enables the rapid acquisition of the co-state mapping results and initiates the subsequent steps of the algorithm. Considering the rapid change of the right side of the dynamic equation near the switching point of the optimal bang-bang control, a high-order variable-step numerical integration method is usually required in Step 5. Considering the terminal relative distance requirement of the rendezvous and docking, the tolerance in Step 5 is usually smaller than that in Step 2, resulting in a longer running time for the MF problem than for the MT and ME problems. If the rendezvous time is less than the solution of the MT problem, or if random initial-time co-state quantities are adopted, the result is non-convergent.

[0168] Example 3: L2 Southern Family Near-Rectilinear Halo Orbit

[0169] The near-rectilinear halo orbit (NRHO) is a special subset of the north-south branch halo orbits around the L1 and L2 points. LunarGateWay utilizes the L2 Southern Family Near-Rectilinear Halo Orbit because of its excellent performance, including low maintenance costs, easy access to the lunar surface and return to Earth, and convenience for round trips to deep space. To asymptotically generate the L2 Southern Family Near-Rectilinear Halo Orbit from the L2 Lyapunov orbit, another state correction function is constructed in this example, which assumes δx = 0.

[0170]

[0171] Figure 5 Shows the process of constructing the L2 Lyapunov orbit to the L2 Southern Family Near-Rectilinear Halo Orbit by alternately using these two correction functions. The ratio in parentheses represents the ratio of the orbital period to the lunar rotation period. Table 2 shows the comparison of the small-thrust rendezvous results and the impulsive rendezvous results.

[0172] Table 2

[0173]

[0174]

[0175] The present invention adopts the same 9:2 NRHO as Lunar GateWay, with a perilune of 3250 km and an orbital period of 6.562 days. Assuming that the chasing spacecraft is at the apolune at the beginning, the target spacecraft has a phase lead of 12.5 degrees, and the rendezvous is set after 0.25 orbital periods. Table 2 shows that compared with the three-body Lambert problem (TBLP) (also known as double-impulse rendezvous), using low-thrust propulsion can reduce fuel consumption by 42.36% to 84.62%. Different control strategies of low-thrust propulsion also affect the final mass, and the final mass is usually negatively correlated with the rendezvous time.

[0176] Example 4: Comparison with other algorithms

[0177] The low-thrust trajectory optimization problem usually encounters two main challenges: the sensitivity of the costate at the initial moment and the optimality of the result. To solve the sensitivity problem, the costate at the initial moment is usually estimated first, and then the indirect shooting method (ISM) is applied to obtain the optimal solution. Two common algorithms for estimating the costate at the initial moment to solve the low-thrust trajectory optimization problem are: particle swarm optimization (PSO) and ISM based on random initial costates.

[0178] The present invention uses an AMD Ryzen5-5600H CPU @ 3.3 GHz and 16 GB of RAM for simulation and conducts it independently 10 times. The CPU running time is calculated using the tic-toc function in MATLAB software, as Figure 6 、 Figure 7 、 Figure 8 shown. The fuel consumption results are as Figure 9 、 Figure 10 、 Figure 11 shown.

[0179] The PSO algorithm consumes a large amount of computing resources and has a low convergence speed, especially in the fuel-optimal two-point boundary value problem. This is due to the continuous generation of new and different information by particles and frequent information exchange. In addition, due to the new information generated near the last generation, the PSO algorithm may converge to a suboptimal solution. In contrast, the ISM algorithm produces optimal results by using Newton iteration to optimize along the gradient direction. However, the costate at the initial moment is highly sensitive, which may make it difficult for the algorithm to converge and require more trial and error for different random initial guesses. In the most sensitive fuel-optimal two-point boundary value problem, compared with the PSO algorithm, the ISM algorithm may require more computing time.

[0180] The comparison of the three algorithms shows that although the PSO algorithm overcomes the sensitivity problem, it has the problems of long running time and non-optimal solutions. The ISM algorithm performs best in the MT problem, but has convergence problems in the ME and MF problems. The ZIGA method produces optimal results in the MT and MF problems, with a slightly longer running time than ISM in the MT problem and much shorter in the MF problem. In addition, the result of the ME problem is an intermediate step in solving the MF problem and does not affect the final optimality. For the MT problem, ZIGA converges in about 2 seconds, the ME problem converges in 0.5 seconds, and the MF problem converges in 15 seconds. Therefore, ZIGA is suitable for calculating the optimal solution of the small-thrust rendezvous problem on the typical Earth-Moon Halo orbit, overcoming the sensitivity problem of the adjoint variables at the initial moment and having no convergence problem.

[0181] The beneficial effects of the embodiments of the present invention are as follows:

[0182] 1. After discretizing the received original small-thrust rendezvous control trajectory and state trajectory of the Earth-Moon Halo orbit, the control quantities and state quantities at the discrete Legendre-Gauss-Radau collocation points are determined, and the technical solution of the objective function that satisfies the rendezvous mission constraints is obtained. The conversion from the original continuous small-thrust rendezvous problem to the nonlinear programming process expands the convergence domain of the original problem and reduces the difficulty of solving the continuous small-thrust rendezvous problem.

[0183] 2. The open-source nonlinear programming solver IPOPT is used to solve the intermediate variables of the objective function, and the intermediate variables are used as inputs. The discrete adjoint variables are calculated through the adjoint mapping equation, and the technical solution of using the discrete adjoint variables as the initial guess of the adjoint variables at the initial moment of the two-point boundary value model; it provides a reasonable initial value guess for the adjoint variables at the initial moment of the two-point boundary value model, overcomes the initial value sensitivity problem compared with the traditional indirect method of selecting random adjoint variables, and greatly improves the calculation speed compared with the heuristic method.

[0184] 3. The Pontryagin minimum principle is applied to ensure that the calculation result is the optimal solution.

[0185] 4. In the calculation of the homotopy iteration process, homotopy parameters with gradually decreasing modulus values and gradually decreasing step sizes are used, which improves the robustness of the calculation process.

[0186] 5. The optimal control rate of the tracking spacecraft can be obtained by integrating the minimum value of the objective function forward.

[0187] 6. The technical solutions disclosed in the present invention use different convergence tolerances respectively, and greatly improve the calculation speed compared with the traditional indirect method and heuristic algorithm with a single convergence tolerance.

[0188] As described above, it is only the specific implementation manner of the present invention, but the protection scope of the present invention is not limited thereto. Any person skilled in the art within the technical scope disclosed by the present invention can easily think of changes or substitutions, which should all be covered within the protection scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the protection scope of the claims described above.

Claims

1. A zero-initial value guessing method for low-thrust rendezvous trajectory optimization of Earth-Moon Halo orbit, characterized in that: The method includes: Step 1: After discretizing the received original Earth-Moon Halo orbit low-thrust rendezvous control trajectory and state trajectory, determine the control quantity and state quantity at the discrete Legendre-Gauss-Radau collocation point, and obtain the objective function that satisfies the rendezvous mission constraints; Step 2: Use the open source nonlinear programming solver IPOPT to solve the intermediate variables of the objective function; Step 3: Take the intermediate variable as input and calculate the discrete co-state quantity through the co-state mapping equation; Step 4: Combined with the Pontryagin minimum principle, the optimal control of the tracking spacecraft is obtained through the Hamiltonian. The two-point boundary value model is composed of the optimal control of the tracking spacecraft, the dynamic equation of the tracking spacecraft, the co-state equation and the terminal constraint. Step 5: Use the discrete costate as the initial guess of the costate at the initial moment of the two-point boundary value model, use the SQP algorithm to optimize the costate at the initial moment, obtain the output result of the two-point boundary value model, and obtain the minimum value of the objective function through the output result; Step 6: Use the output results to forward integrate the tracking spacecraft dynamic equations and co-state equations to obtain the optimal control rate of the tracking spacecraft and the low-thrust rendezvous trajectory of the Earth-Moon Halo orbit.

2. The method according to claim 1, characterized in that In step 1, Legendre-Gauss-Radau collocation, P N (τ)+P N-1 (τ) roots of the polynomial; Where P N (τ) is the Nth order Legendre polynomial, P N-1 (τ) is the N-1 order Legendre polynomial, where P N (τ) is calculated by the following formula: P N Replace N in (τ) with N-1, and we get P N-1 (τ), where τ is the independent variable, d is the differential symbol, and N is the specified order.

3. The method according to claim 1, characterized in that In step 1, the objective function that satisfies the rendezvous task constraints includes: the objective function of time optimal control and the objective function of fuel optimal control; wherein, (1) The objective function of time optimal control is: Among them, time optimal control refers to satisfying the tracking spacecraft control rate in the rendezvous mission so that the time required for rendezvous is the shortest; (2) The objective function of fuel optimal control is: Among them, optimal fuel control refers to satisfying the tracking spacecraft control rate in the rendezvous mission so that the fuel required for rendezvous is minimized; Where J1 is the time required for the intersection, t0 is the initial time, t f is the terminal time, J2 is the fuel required for rendezvous, and U is the tracking spacecraft control vector.

4. The method according to claim 3, characterized in that In step 2, the intermediate variables include: Gaussian quadrature weights, Lagrange multiplier matrix and differential matrix.

5. The method according to claim 4, characterized in that In step 3, the intermediate variables are used as inputs, and the discrete co-state quantities calculated by the co-state mapping equation include: the discrete co-state quantity of time optimal control and the discrete co-state quantity of fuel optimal control; wherein the co-state mapping equation is: In the formula, the discrete co-state quantity of time optimal control and the discrete co-state quantity of fuel optimal control calculated by the co-state mapping equation are given by the discrete co-state quantity at the Legendre-Gauss-Radau collocation point and discrete co-state quantities at the endpoints composition, is the transposed matrix D of the differential matrix D T The M+1th column, Λ k represents the kth row of the Lagrange multiplier matrix Λ, ω k is the Gaussian quadrature weight ω k The value at the kth Legendre-Gauss-Radau collocation point, where k is the number of the collocation point and M is the total number of collocation points.

6. The method according to claim 1, characterized in that In step 4, in combination with the Pontryagin minimum principle, the optimal control of the tracking spacecraft is obtained through the Hamiltonian, including: obtaining the time optimal control U1 of the tracking spacecraft through the Hamiltonian H1 of time optimal control and obtaining the fuel optimal control U2 of the tracking spacecraft through the Hamiltonian H2 of fuel optimal control; wherein, (1) The Hamiltonian H1 of time optimal control is: In the formula, Λ=[Λ r , Λ v ,λ m ] represents the co-state quantity, Λ r is the position covariance, Λ v is the velocity covariate, λ m is the mass costate; v is the velocity of the tracking spacecraft, g(r) is the term related to the position r of the tracking spacecraft in the tracking spacecraft dynamics equation, h(v) is the term related to the velocity v of the tracking spacecraft in the tracking spacecraft dynamics equation, U is the tracking spacecraft control vector, T max is the maximum available thrust, m is the mass of the tracking spacecraft, u is the modulus of the tracking spacecraft control vector, I sp is the specific impulse, g0 is the gravitational acceleration at sea level; According to the Pontryagin minimum principle, the optimal time control of the tracking spacecraft is obtained through the Hamiltonian H1 of the time optimal control: In the formula, λ v is the modulus of the velocity co-state; (2) The Hamiltonian H2 for optimal fuel control is: According to the Pontryagin minimum principle, the optimal fuel control of the tracking spacecraft is obtained through the Hamiltonian H2 of the optimal fuel control: Among them, u2∈[0,T max ], which represents the modulus of the optimal control vector for tracking spacecraft fuel.

7. The method according to claim 6, characterized in that In step 4, the two-point boundary value model includes: a two-point boundary value model for time optimal control and a two-point boundary value model for fuel optimal control; wherein, The two-point boundary value model of time optimal control consists of the tracking spacecraft dynamic equation, the tracking spacecraft time optimal control, the co-state equation of time optimal control and the terminal constraint of time optimal control. The two-point boundary value model of fuel optimal control consists of the tracking spacecraft dynamic equation, the tracking spacecraft fuel optimal control, the co-state equation of fuel optimal control and the terminal constraint of fuel optimal control; among them, the tracking spacecraft dynamic equation is: In the formula, is the derivative of the tracking spacecraft state quantity, f(X,U) is the dynamic equation obeyed by the tracking spacecraft state quantity, X=[r,v] is the tracking spacecraft state quantity, r is the tracking spacecraft position, v is the tracking spacecraft speed, To track the derivative of the spacecraft position, To track the derivative of the spacecraft velocity, is the derivative of the tracking spacecraft mass, u = |U| represents the modulus of the tracking spacecraft control vector U; (1) The co-state equation for time optimal control is: In the formula, is the derivative of the co-state quantity of time optimal control, H1 is the Hamiltonian of time optimal control; The terminal constraints for time-optimal control are: In the formula, φ1([X,Λ1] T ,t0,t f ) is the terminal constraint value of time optimal control, Λ1 is the co-state quantity of time optimal control of tracking spacecraft, t0 is the initial time, t f is the terminal moment, r(t f ) is to track the terminal position of the spacecraft, r f is the terminal position of the target spacecraft, v(t f ) is the terminal velocity of the spacecraft, v f is the terminal velocity of the target spacecraft, λ m (t f ) is the mass co-state quantity tracking the spacecraft at the terminal time, H1(t f ) is the Hamiltonian at the terminal moment, 08 is the eight-dimensional zero vector; (2) The co-state equation for optimal fuel control is: In the formula, is the derivative of the co-state quantity of the optimal fuel control, H2 is the Hamiltonian of the optimal fuel control; The terminal constraint for optimal fuel control is: In the formula, φ2([X,Λ2] T ,t0,t f ) is the terminal constraint value of the optimal fuel control, Λ2 is the co-state quantity of the optimal fuel control of the tracking spacecraft, and 07 is the seven-dimensional zero vector.

8. The method according to claim 7, characterized in that In step 5, the SQP algorithm is used to optimize the initial co-state quantity to obtain the output result of the two-point boundary value model. The minimum value of the objective function obtained through the output result includes: the minimum value of the objective function of time optimal control and the minimum value of the objective function of fuel optimal control; wherein, (1) The calculation method of the minimum value of the objective function of time optimal control includes: The SQP algorithm is used to optimize the discrete co-state quantity of the time optimal control; the calculation is stopped when the control rate calculated by the two-point boundary value model of the time optimal control satisfies the terminal constraint of the time optimal control, and the initial co-state quantity of the time optimal control and the terminal time of the time optimal control are output, and the minimum value of the objective function of the time optimal control is obtained from the initial co-state quantity of the time optimal control and the terminal time of the time optimal control; (2) The method for calculating the minimum value of the objective function of the optimal fuel control includes: Construct the objective function J of energy optimal control that meets the rendezvous mission constraints; Use the open source nonlinear programming solver IPOPT to solve the intermediate variables of the objective function of energy optimal control; Taking the intermediate variables as input, the discrete co-state quantity of energy optimal control is calculated through the co-state mapping equation; Combined with the Pontryagin minimum principle, the energy optimal control U3 of the tracking spacecraft is obtained through the Hamiltonian H3 of energy optimal control; The discrete co-state quantity of energy optimal control is used as the initial guess of the co-state quantity of the two-point boundary value model of fuel optimal control at the initial moment, and the SQP algorithm is used to optimize the co-state quantity at the initial moment. Based on the homotopy iteration process, the Hamiltonian of fuel optimal control and the energy optimal control of the tracking spacecraft are used to obtain the Hamiltonian of fuel optimal control and the fuel optimal control of the tracking spacecraft; when the control rate calculated by the two-point boundary value model of fuel optimal control meets the terminal constraint of fuel optimal control, the calculation is stopped, and the co-state quantity of the initial moment of fuel optimal control is output, and the minimum value of the objective function of fuel optimal control is obtained from the co-state quantity of the initial moment of fuel optimal control; wherein, The objective function of energy optimal control that meets the rendezvous mission constraints is: The Hamiltonian H3 for energy optimal control is: Energy-optimal control of tracking spacecraft Among them, u3∈[0,T max ], which represents the modulus of the energy-optimal control vector of the tracking spacecraft.

9. The method according to claim 8, characterized in that In step 5, the iterative calculation method of the homology iteration process is: In the formula, is the homotopy parameter of each step; i = 10, 9, …, 1, corresponding to the homotopy parameter ∈ = 1 to ∈ = 0; where ∈ = 1 corresponds to the Hamiltonian of energy optimal control, ∈ = 0 corresponds to the Hamiltonian of fuel optimal control, and H4 represents the Hamiltonian in the homotopy iteration process; In each iteration step, the SQP algorithm optimizes the initial co-state quantity, stops after satisfying the terminal constraint of the fuel optimal control, and enters the next iteration step; when the homology parameter ∈=0, the calculated Hamiltonian of the energy optimal control is the same as the Hamiltonian of the fuel optimal control, and the Hamiltonian of the fuel optimal control is obtained.

10. The method according to claim 1, characterized in that In step 6, the integration process of the forward integration uses a Runge-Kutta 7 / 8th order variable step-size integrator for forward integration.

Images