An emergency return orbit reachable domain calculation method based on interval analysis

Abstract

The present application relates to a kind of based on interval analysis emergency return orbit reachable domain calculation method, comprising: selecting the design variable of single emergency return orbit solution, interval variable and optimization variable for interval analysis are determined based on design variable, according to the search space of interval variable, interval variable is divided into corresponding multiple initial subinterval;Emergency return orbit is solved to the interval upper limit and interval lower limit of each initial subinterval;The initial subinterval of existence solution is obtained, and initial subinterval is contracted, and transition subinterval is obtained;Judge whether the interval boundary of transition subinterval exists the solution of emergency return orbit, if there is, then the foregoing step is repeated to further contract transition subinterval, and minimum subinterval is obtained;After all initial subinterval in search space is handled, all about the solution of emergency return orbit in minimum subinterval is obtained, and then the reachable domain of the emergency return orbit is obtained.

Description

Technical Field

[0001] This invention relates to the field of spacecraft space mission planning and design technology, and in particular to a method for calculating the reachability domain of an emergency return trajectory based on interval analysis. Background Technology

[0002] In manned lunar exploration missions, the lunar orbit phase is one of the most crucial stages of the entire mission. If the spacecraft malfunctions or encounters other emergencies during this phase, the mission must be aborted and the astronauts must return to Earth prematurely for their safety. Because the spacecraft is relatively far from Earth at this point, multiple maneuvers are not feasible in an emergency. Therefore, a single maneuver can be performed to enter an emergency return orbit and return to Earth. This single maneuver can employ either a coplanar orbit change or a non-coplanar orbit change. However, compared to a coplanar orbit change, a non-coplanar orbit change offers a wider return window and is more suitable for emergency needs requiring return at any time.

[0003] For near-Earth spacecraft, the reachability domain is generally defined as the set of all spatial locations that the spacecraft can reach under certain constraints. As an intuitive concept describing the reachable locations of a spacecraft, it can be directly used in space mission design and evaluation, including on-orbit servicing, formation flying, and trajectory safety assessment. The theoretical research on reachability domains mentioned above is usually referred to as the absolute reachability domain. Based on this, the relative reachability domain is derived and is often used in spacecraft collision avoidance planning research. For manned lunar emergency return orbits, the reachability domain concept is usually characterized as the Earth's reachable region of the reentry point, that is, the area of ​​the Earth's surface that the spacecraft can reach under certain orbital constraints, such as the latitude and longitude distribution of the landing point. However, the more general concept of Earth's reachability domain also refers to the distribution of characteristic state parameters of the emergency return orbit in high-dimensional space.

[0004] Research on the reachability region of emergency return orbits can serve not only as an engineering criterion for determining the reachability of spacecraft to Earth in emergency situations, but also potentially yield mathematical conclusions regarding the existence of solutions to the boundary value problem between two points in the Earth-Moon orbit. Currently, research on the reachability region problem of Earth-Moon orbits is limited both domestically and internationally, and a mature theory or methodology has not yet been established. Existing techniques mainly employ ergonomic calculation methods to obtain the set of all orbits satisfying the constraints, thereby acquiring the reachability region. However, this technique still suffers from drawbacks such as high computational complexity and long processing time. Summary of the Invention

[0005] The purpose of this invention is to provide a method for calculating the reachability domain of an emergency return trajectory based on interval analysis.

[0006] To achieve the above-mentioned objectives, this invention provides a method for calculating the reachability domain of an emergency return trajectory based on interval analysis, characterized by comprising:

[0007] S1. Select the design variables when solving a single emergency return trajectory, determine the interval variables and optimization variables for interval analysis based on the design variables, and divide the interval variables into multiple corresponding initial sub-intervals according to the search space of the interval variables;

[0008] S2. Based on the optimization variables, solve for the upper and lower bounds of each initial sub-interval to determine the emergency return trajectory;

[0009] S3. Obtain the initial sub-interval where a solution exists, and shrink the initial sub-interval to obtain the shrunken transition sub-interval;

[0010] S4. Determine whether there is a solution for the emergency return track at the boundary of the transition sub-interval. If there is, repeat step S3 to further shrink the transition sub-interval until the minimum interval length limit is met, and obtain the minimum sub-interval.

[0011] S5. Repeat steps S2 to S4. After processing all the initial sub-intervals in the search space, obtain all solutions related to the emergency return trajectory in the smallest sub-interval, and then obtain the reachable domain of the emergency return trajectory.

[0012] According to one aspect of the present invention, step S1, which involves selecting design variables for solving a single emergency return trajectory, determining interval variables and optimization variables for interval analysis based on the design variables, and dividing the interval variables into multiple corresponding initial sub-intervals according to the search space of the interval variables, includes:

[0013] S11. The design variable for solving the emergency return trajectory is determined to be (λ). ex φ ex ,Δt1), where λ ex φ represents the longitude at the point where the moon influences the sphere's exit. ex Δt1 represents the latitude of the exit point of the lunar influence sphere, and Δt1 represents the flight time from the emergency maneuver point to the lunar influence sphere.

[0014] S12. Select the longitude and latitude of the lunar influence sphere exit point from the design variables as the interval variables, and select the flight time from the emergency maneuver point to the lunar influence sphere from the design variables as the optimization variables;

[0015] S13. Based on the search space of the interval variables, each interval variable is divided into a plurality of corresponding initial sub-intervals; wherein the initial sub-intervals are a first initial sub-interval and a second initial sub-interval, the first initial sub-intervals are p in number and are obtained based on longitude division at the exit point of the lunar influence sphere, and the second initial sub-intervals are q in number and are obtained based on latitude division at the exit point of the lunar influence sphere.

[0016] According to one aspect of the present invention, step S2, which involves solving for the upper and lower bounds of each initial sub-interval based on the optimization variables, includes:

[0017] S21. Select one of the first initial sub-intervals and one of the second initial sub-intervals respectively, and obtain the upper or lower bound of the interval in the first initial sub-interval and the second initial sub-interval;

[0018] S22. Using the upper or lower bound of the first and second initial sub-intervals as the longitude and latitude of the emergency return trajectory at the lunar influence sphere exit point, construct the position vector of the lunar influence sphere exit point in the lunar-fixed system, which is expressed as:

[0019]

[0020] Wherein, the superscript MF represents the lunar fixed system, ρ M The radius of the sphere is affected by the moon;

[0021] S23. Based on the optimized variables and the position vector of the lunar influence sphere exit point in the lunar-fixed coordinate system, obtain the position vector of the lunar influence sphere exit point in the lunar-centered J2000 coordinate system, which is expressed as:

[0022]

[0023] Wherein, the superscript MJ represents the lunar center J2000 coordinate system, t ex At the exit point time, M(t) ex ) for t ex The coordinate transformation matrix from the fixed lunar coordinate system to the J2000 coordinate system at the lunar center.

[0024] S24. Based on the optimized variables, the position vector of the lunar influence sphere exit point in the lunar-fixed coordinate system, and the position vector of the lunar influence sphere exit point in the lunar center J2000 coordinate system, obtain the initial velocity vector V of the emergency return trajectory in the lunar center segment. cr0 Terminal velocity vector V ex and the emergency maneuver pulse vector ΔV; where the emergency maneuver pulse vector ΔV is expressed as:

[0025] ΔV=Vcr0 -V0

[0026] Where V0 represents the velocity vector before the emergency maneuver;

[0027] S25. Based on the exit point time, the position vector of the lunar influence sphere exit point in the lunar center J2000 coordinate system and the terminal velocity vector V. ex The position and velocity vector of the lunar influence sphere exit point in the geocentric J2000 coordinate system are obtained, and are expressed as follows:

[0028]

[0029]

[0030] in, and They are t ex The position and velocity vector of the Moon in the Earth-centered J2000 coordinate system at any given time can be obtained by solving the JPL DE430 ephemeris.

[0031] S26. By transforming the position and velocity vector of the lunar influence sphere exit point in the geocentric J2000 coordinate system, the orbital elements of the geocentric segment of the emergency return orbit are obtained. The transfer time from the lunar influence sphere exit point to the reentry point is then calculated, expressed as:

[0032]

[0033] Where, μ E Let a be the Earth's gravitational constant. gr e is the semi-major axis of the emergency return orbit for the geocentric segment. gr Let E be the eccentricity of the emergency return trajectory in the geocentric segment, and E be the angle of abduction at the reentry point, which can be obtained from the true angle of abduction at the reentry point.

[0034]

[0035] The true anterior angle is calculated using the following formula:

[0036]

[0037] Where, r vcp and r RP These are the geocentric distances from the vacuum perigee and the reentry point, respectively, from which the orbital elements at the reentry point can be obtained;

[0038] S27. Further constraints are set on the emergency return orbit terminal reentry point, which are expressed as follows:

[0039]

[0040] Among them, h RP and γ RP h represents the altitude and reentry angle of the reentry point aimed at by the emergency return orbit terminal. f and γ f These are the altitude and reentry angle of the emergency return orbit terminal reentry point, respectively.

[0041] S28. The objective function of the emergency return trajectory is set and optimized using a sequential quadratic programming algorithm to obtain the emergency return trajectory that satisfies the constraints; wherein, the objective function is to minimize the emergency pulse, which is expressed as:

[0042] J=(||ΔV||) min

[0043] Where J represents the objective function and ΔV represents the emergency maneuver pulse vector.

[0044] According to one aspect of the present invention, in step S3, the initial sub-interval with a solution is obtained, and the initial sub-interval is shrunk to obtain the shrunk transition sub-interval. The initial sub-interval with a solution is shrunk using the midpoint of the initial sub-interval as the boundary to obtain two interval variables to obtain two transition sub-intervals.

[0045] According to one aspect of the invention, in step S3, if the initial sub-interval does not contain the solution, then the initial sub-interval is eliminated.

[0046] According to one aspect of the present invention, in step S4, it is determined whether there is a solution for the emergency return track at the boundary of the transition sub-interval. If there is, proceed to step S4; if there is no solution, eliminate it.

[0047] According to one aspect of the present invention, the method of interval analysis is introduced into the solution of the reachable domain of the emergency return trajectory, which effectively improves the computational efficiency and reduces the computational load.

[0048] According to one aspect of the present invention, the calculation method is based on interval analysis to solve for the reachability domain of an emergency return orbit. It improves upon the large-scale calculation method of traversing and searching design variables by using interval variables for calculation. Compared to the existing method of traversing and solving for the reachability domain, this significantly improves the calculation speed. The method of the present invention is simple to implement in engineering and has good prospects for widespread application. It can not only solve for the reachability domain of emergency return orbits but also be applied to the reachability domain of other Earth-Moon space orbits. Attached Figure Description

[0049] Figure 1This is a schematic diagram illustrating the steps of an emergency return trajectory reachability domain calculation method according to an embodiment of the present invention;

[0050] Figure 2 This is a flowchart schematically illustrating an emergency return trajectory reachability domain calculation method according to an embodiment of the present invention;

[0051] Figure 3 This is a schematic representation of the emergency return trajectory reachability region diagram obtained by the emergency return trajectory reachability region calculation method according to an embodiment of the present invention;

[0052] Figure 4 This is a schematic representation of the reachable region of the emergency return trajectory obtained using a traversal solution method. Detailed Implementation

[0053] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the accompanying drawings used in the embodiments will be briefly described below.

[0054] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments. The embodiments cannot be described in detail here, but the embodiments of the present invention are not limited to the following embodiments.

[0055] Combination Figure 1 and Figure 2 As shown, according to one embodiment of the present invention, under the premise of known relevant orbital constraints, the set of orbital parameters reachable from Earth by an emergency return orbit is solved. Furthermore, the present invention provides a method for calculating the reachability domain of an emergency return orbit based on interval analysis, comprising:

[0056] S1. Select the design variables when solving a single emergency return trajectory, determine the interval variables and optimization variables for interval analysis based on the design variables, and divide the interval variables into multiple corresponding initial sub-intervals according to the search space of the interval variables;

[0057] S2. Based on the optimization variables, solve for the upper and lower bounds of each initial sub-interval to determine the emergency return trajectory;

[0058] S3. Obtain the initial subinterval where a solution exists, and shrink the initial subinterval to obtain the shrunken transition subinterval;

[0059] S4. Determine whether there is a solution for the emergency return track at the boundary of the transition sub-interval. If there is, repeat step S3 to further shrink the transition sub-interval until the minimum interval length limit is met, and obtain the minimum sub-interval.

[0060] S5. Repeat steps S2 to S4. After processing all the initial subintervals in the search space, obtain all solutions for the emergency return trajectory in the smallest subinterval, and then obtain the reachable domain of the emergency return trajectory.

[0061] Combination Figure 1 and Figure 2 As shown, according to one embodiment of the present invention, step S1, which involves selecting design variables for solving a single emergency return trajectory, determining interval variables and optimization variables for interval analysis based on the design variables, and dividing the interval variables into multiple corresponding initial sub-intervals according to the search space of the interval variables, includes:

[0062] S11. When solving for a single free return trajectory, the design variable for solving the emergency return trajectory is determined to be (λ). ex φ ex ,Δt1), where λ ex φ represents the longitude at the point where the moon influences the sphere's exit. ex Δt1 represents the latitude of the exit point of the lunar influence sphere, and Δt1 represents the flight time from the emergency maneuver point to the lunar influence sphere. In this embodiment, given the initial state of the emergency return, an emergency return trajectory can be determined by a set of design variables.

[0063] S12. Select the longitude and latitude of the lunar influence sphere exit point as interval variables, which are respectively represented as [λ]. ex ] and [φ ex The flight time Δt1 from the emergency maneuver point to the lunar influence sphere is selected as the optimization variable.

[0064] S13. Based on the search space of the interval variables, each interval variable (i.e., [λ]) ex ] and [φ ex The region is divided into multiple initial sub-intervals; the initial sub-intervals are the first initial sub-interval and the second initial sub-interval, the first initial sub-interval consists of p regions and is obtained based on the longitude division at the exit point of the lunar influence sphere, and the second initial sub-interval consists of q regions and is obtained based on the latitude division at the exit point of the lunar influence sphere.

[0065] Combination Figure 1 and Figure 2 As shown, according to one embodiment of the present invention, in step S2, in the step of solving for the upper and lower bounds of each initial sub-interval based on the optimization variables, the interval variable [λ] ex The upper and lower bounds of a certain initial subinterval (i.e., the first initial subinterval) of are respectively and Interval variable [φex The upper and lower bounds of a certain initial subinterval (i.e., the second initial subinterval) of ] are respectively and include:

[0066] S21. Select a first initial sub-interval and a second initial sub-interval respectively, and obtain the upper or lower bound of the interval in the first and second initial sub-intervals; in this embodiment, when solving for the emergency return trajectory of the interval boundary, two processes are used: solving for the upper bound of the interval and solving for the lower bound of the interval. Taking the solution for the lower bound of the interval as an example, the lower bound of the first initial sub-interval and the second initial sub-interval is selected as the input for the solution, that is, the input is... and Correspondingly, when solving using the upper bound of an interval, the upper bound of a first initial subinterval and a second initial subinterval is selected as the input for the solution. That is, the input is... and

[0067] S22. Using the upper or lower bound of either the first or second initial sub-interval as the longitude and latitude of the emergency return trajectory at the lunar influence sphere exit point, construct the position vector of the lunar influence sphere exit point in the lunar-solid system, which is expressed as:

[0068]

[0069] Wherein, the superscript MF represents the lunar fixed system, ρ M The radius of the sphere is affected by the moon;

[0070] S23. Based on the optimized variables and the position vector of the lunar influence sphere exit point in the lunar-fixed coordinate system, the position vector of the lunar influence sphere exit point in the lunar-centered J2000 coordinate system is obtained, which is expressed as:

[0071]

[0072] Wherein, the superscript MJ represents the lunar center J2000 coordinate system, t ex At the exit point time, M(t) ex ) for t ex The coordinate transformation matrix from the fixed lunar coordinate system to the J2000 coordinate system at the lunar center.

[0073] S24. Based on the optimized variables, the position vector of the lunar influence sphere exit point in the lunar-fixed coordinate system, and the position vector of the lunar influence sphere exit point in the lunar center J2000 coordinate system, obtain the initial velocity vector V of the lunar center segment of the emergency return orbit. cr0 Terminal velocity vector V ex And the emergency maneuver pulse vector ΔV; in this embodiment, the Lambert algorithm is used to solve for the initial velocity vector V.cr0 Terminal velocity vector V ex The emergency maneuver pulse vector ΔV is represented as:

[0074] ΔV=V cr0 -V0

[0075] Wherein, V0 represents the velocity vector before the emergency maneuver, which is a known quantity;

[0076] S25. Based on the exit point time, the position vector and terminal velocity vector V of the lunar influence sphere's exit point in the lunar center J2000 coordinate system. ex The position and velocity vector of the lunar influence sphere exit point in the geocentric J2000 coordinate system are obtained, and are expressed as follows:

[0077]

[0078]

[0079] in, and They are t ex The position and velocity vector of the Moon in the Earth-centered J2000 coordinate system at any given time can be obtained by solving the JPL DE430 ephemeris.

[0080] S26. By transforming the position and velocity vector of the lunar influence sphere exit point in the geocentric J2000 coordinate system, the orbital elements of the geocentric segment of the emergency return orbit are obtained. The transfer time from the lunar influence sphere exit point to the reentry point is then calculated, expressed as:

[0081]

[0082] Where, μ E Let a be the Earth's gravitational constant. gr e is the semi-major axis of the emergency return orbit for the geocentric segment. gr Let E be the eccentricity of the emergency return trajectory in the geocentric segment, and E be the angle of abduction at the reentry point, which can be obtained from the true angle of abduction at the reentry point.

[0083]

[0084] The true anterior angle is calculated using the following formula:

[0085]

[0086] Where, r vcp and r RP These are the geocentric distances from the vacuum perigee and the reentry point, respectively, from which the orbital elements at the reentry point can be obtained;

[0087] S27. Further constraints are set for the emergency return orbit terminal reentry point, which are expressed as follows:

[0088]

[0089] Among them, h RP and γ RP These are the altitude and reentry angle of the reentry point aimed at by the emergency return orbit terminal, h. f and γ f These are the altitude and reentry angle of the emergency return orbit terminal reentry point, respectively. In this step, after obtaining the orbital elements at the reentry point in the previous steps, these two quantities can be obtained.

[0090] S28. The objective function for the emergency return trajectory is set and optimized using a sequential quadratic programming algorithm to obtain an emergency return trajectory that satisfies the constraints. In this embodiment, fuel may be scarce in emergency situations, so it is often required that emergency maneuvers consume as little energy as possible. Therefore, the objective function for calculating the emergency return trajectory is set to minimize the emergency impulse, which is expressed as:

[0091] J=(||ΔV||) min

[0092] Where J represents the objective function and ΔV represents the emergency maneuver pulse vector.

[0093] In this embodiment, in step S2, the corresponding computational model is obtained in steps S21 to S27, and the corresponding solution algorithm is executed in step S28. The sequential quadratic programming algorithm is executed based on the explicit computational model, constraints and objective function.

[0094] Combination Figure 1 and Figure 2 As shown, according to one embodiment of the present invention, in step S3, in the step of obtaining the initial sub-interval with a solution and shrinking the initial sub-interval to obtain the shrunken transition sub-interval, the midpoint segmentation method is used to shrink the initial sub-interval with a solution; wherein, the initial sub-interval is further divided into two interval variables with the midpoint of the initial sub-interval as the boundary to obtain two transition sub-intervals.

[0095] like Figure 1 As shown, according to one embodiment of the present invention, in step S3, if there is no solution in the initial sub-interval, the initial sub-interval is eliminated.

[0096] like Figure 1 As shown, according to one embodiment of the present invention, in step S4, determining whether there is a solution for an emergency return track at the boundary of the transition sub-interval is a step in which, if no solution exists, it is eliminated.

[0097] Combination Figure 1 and Figure 2 As shown, according to one embodiment of the present invention, in step S5, steps S2 to S4 are repeated. After processing all initial sub-intervals in the search space, all solutions for the emergency return trajectory within the smallest sub-interval are obtained, thereby obtaining the reachable domain of the emergency return trajectory. Since the reachable domain has been obtained based on the aforementioned process, in this step, some important parameters are summarized based on the obtained reachable domain to obtain a corresponding set. The reachable domain of the emergency return trajectory refers to a high-dimensional parameter set consisting of the trajectory parameters of the emergency return terminal re-entry point that satisfy the constraints. In practical tasks, some parameters are often of particular interest, such as the trajectory inclination and right ascension of the ascending node of the emergency return trajectory re-entry point.

[0098]

[0099] Among them, i RP Ω represents the orbital inclination at the reentry point. RP The right ascension of the ascending node represents the re-entry point. This represents a partial set of parameters of the reachable domain.

[0100] To better understand the above-described solutions of the present invention, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0101] Example 1

[0102] The emergency return time is set at 00:00 on February 17, 2028. In the lunar center J2000 coordinate system, the lunar orbital altitude is 100km, the eccentricity is 0, the orbital inclination is 156°, the right ascension of the ascending node is 194°, and the lunar perigee's latitude argument and true anomaly angle are 0°. The maximum permissible emergency return duration is 3 days. The reentry point altitude is 122km, and the reentry angle is -6°. The interval variable [λ]... ex The search space is -180° to 180°, and the interval variable [φ] is... ex The search space is -90° to 90°. The interval variable [λ] is... ex ] and [φ ex The intervals are divided into 40 initial subintervals (i.e., p = 40, q = 40). To verify the accuracy and superiority of the calculation method of this invention, the reachable region of the emergency return trajectory is calculated using a traditional traversal solution method. During the traversal solution, λ is searched in steps of 1°. ex and φ ex .

[0103] Figure 3 The reachable region of the emergency return trajectory obtained by the calculation method of this invention is presented. Figure 4The reachable region of the emergency return trajectory obtained by the traversal solution method is presented. A comparison of the two figures shows that the geometric characteristics of the reachable region obtained by the calculation method of this invention are basically similar to those obtained by the traversal solution method, thus verifying the accuracy of the calculation method of this invention. Regarding solution time, the solution time of the calculation method of this invention is 2424.261s, while the solution time of the traversal solution method is 20324.349s. Therefore, compared with the traversal solution method, the calculation method of this invention can shorten the solution time by 88.1%, greatly improving computational efficiency.

[0104] The above description is merely an example of a specific solution of the present invention. For any devices and structures not described in detail herein, it should be understood that they are implemented using common devices and methods already available in the art.

[0105] The above description is merely one embodiment of the present invention and is not intended to limit the invention. Various modifications and variations can be made to the invention by those skilled in the art. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for calculating the reachability domain of an emergency return trajectory based on interval analysis, characterized in that, include: S1. Select the design variables when solving a single emergency return trajectory, determine the interval variables and optimization variables for interval analysis based on the design variables, and divide the interval variables into multiple corresponding initial sub-intervals according to the search space of the interval variables; These include: S11. Determine the design variables for solving the emergency return trajectory as ( ). λ ex , φ ex Δ t 1), among which, λ ex This indicates the longitude at the point where the moon influences the Earth's exit. φ ex Δ represents the latitude of the point where the moon influences the sphere's exit. t 1 indicates the flight time from the emergency mobile location to the lunar influence sphere; S12. Select the longitude and latitude of the lunar influence sphere exit point from the design variables as the interval variables, and select the flight time from the emergency maneuver point to the lunar influence sphere from the design variables as the optimization variables; S13. Based on the search space of the interval variables, each interval variable is divided into multiple corresponding initial sub-intervals; wherein the initial sub-intervals are a first initial sub-interval and a second initial sub-interval, and the first initial sub-interval is... p The number of intervals is determined based on the longitude division at the exit point of the lunar influence sphere, and the second initial sub-interval is... q The number of samples is based on the latitudinal segmentation of the exit point of the lunar influence sphere. S2. Based on the optimization variables, solve for the upper and lower bounds of each initial sub-interval to determine the emergency return trajectory; S3. Obtain the initial sub-interval where a solution exists, and shrink the initial sub-interval to obtain the shrunken transition sub-interval; S4. Determine whether there is a solution for the emergency return track at the boundary of the transition sub-interval. If there is, repeat step S3 to further shrink the transition sub-interval until the minimum interval length limit is met, and obtain the minimum sub-interval. S5. Repeat steps S2 to S4. After processing all the initial sub-intervals in the search space, obtain all solutions related to the emergency return trajectory in the smallest sub-interval, and then obtain the reachable domain of the emergency return trajectory.

2. The emergency return trajectory reachability domain calculation method according to claim 1, characterized in that, Step S2, which involves solving for the emergency return trajectory based on the optimization variables for the upper and lower bounds of each initial sub-interval, includes: S21. Select one of the first initial sub-intervals and one of the second initial sub-intervals respectively, and obtain the upper or lower bound of the interval in the first initial sub-interval and the second initial sub-interval; S22. Using the upper or lower bound of the first and second initial sub-intervals as the longitude and latitude of the emergency return trajectory at the lunar influence sphere exit point, construct the position vector of the lunar influence sphere exit point in the lunar-fixed system, which is expressed as: The superscript MF indicates the lunar fixed system. ρ M The radius of the sphere is affected by the moon; S23. Based on the optimized variables and the position vector of the lunar influence sphere exit point in the lunar-fixed coordinate system, obtain the position vector of the lunar influence sphere exit point in the lunar-centered J2000 coordinate system, which is expressed as: Wherein, the superscript MJ indicates the lunar center J2000 coordinate system, t ex For the exit point time, M ( t ex )for t ex The coordinate transformation matrix from the lunar fixed coordinate system to the lunar center J2000 coordinate system; S24. Based on the optimized variables, the position vector of the lunar influence sphere exit point in the lunar-fixed coordinate system, and the position vector of the lunar influence sphere exit point in the lunar center J2000 coordinate system, obtain the initial velocity vector of the emergency return trajectory in the lunar center segment. V cr0 Terminal velocity vector V ex and emergency mobility pulse vector ΔV Among them, emergency mobility pulse vector ΔV Represented as: ΔV = V cr0 - V 0 in, V 0 represents the velocity vector before the emergency maneuver; S25. Based on the exit point time, the position vector of the lunar influence sphere exit point in the lunar center J2000 coordinate system and the terminal velocity vector. V ex The position and velocity vector of the lunar influence sphere exit point in the geocentric J2000 coordinate system are obtained, and are expressed as follows: in, and They are t ex The position and velocity vector of the Moon in the Earth-centered J2000 coordinate system at any given time are obtained by solving the JPL DE430 ephemeris. S26. By transforming the position and velocity vector of the lunar influence sphere exit point in the geocentric J2000 coordinate system, the orbital elements of the geocentric segment of the emergency return orbit are obtained. The transfer time from the lunar influence sphere exit point to the reentry point is then calculated, expressed as: in, μ E The gravitational constant of Earth, a gr e is the semi-major axis of the emergency return orbit for the geocentric segment. gr The eccentricity of the emergency return trajectory in the geocentric segment. E The deviated angle of reentry is obtained from the true angle of reentry: The true anterior angle is calculated using the following formula: in, r vcp and r RP These are the distances from the Earth's center to the vacuum perigee and the reentry point, respectively, from which the orbital elements at the reentry point can be obtained; S27. Further constraints are set on the emergency return orbit terminal reentry point, which are expressed as follows: in, h RP and γ RP These are the altitude and reentry angle of the reentry point aimed at by the emergency return orbit terminal, respectively. h f and γ f These are the altitude and reentry angle of the emergency return orbit terminal reentry point, respectively. S28. The objective function of the emergency return trajectory is set and optimized using a sequential quadratic programming algorithm to obtain the emergency return trajectory that satisfies the constraints; wherein, the objective function is to minimize the emergency pulse, which is expressed as: in, J Describe the objective function. ΔV This represents the emergency maneuver pulse vector.

3. The emergency return trajectory reachability domain calculation method according to claim 2, characterized in that, In step S3, the initial sub-interval with a solution is obtained, and the initial sub-interval is shrunk to obtain the shrunk transition sub-interval. The midpoint segmentation method is used to shrunk the initial sub-interval with a solution. Specifically, the initial sub-interval is further divided into two interval variables with the midpoint of the initial sub-interval as the boundary to obtain two transition sub-intervals.

4. The emergency return trajectory reachability domain calculation method according to claim 3, characterized in that, In step S3, if the initial sub-interval does not contain the solution, then the initial sub-interval is eliminated.

5. The emergency return trajectory reachability domain calculation method according to claim 4, characterized in that, In step S4, it is determined whether there is a solution for the emergency return track at the boundary of the transition sub-interval. If there is, proceed to step S4; otherwise, eliminate the solution.

Images