Spacecraft approaching maneuvering data set construction method

Abstract

The invention discloses a spacecraft approaching maneuvering data set construction method, and belongs to the technical field of spacecraft orbit design and optimization, and the method comprises the following steps: 1, determining the maneuvering type of a maneuvering satellite, and setting the input of a maneuvering generation algorithm; step 2, calculating the maximum number of flight turns and task time of the target satellite; step 3, calculating task time based on the current number of turns, and solving an initial guess; step 4, constructing a nonlinear equation representing successful maneuvering execution of the maneuvering satellite, and solving a convergence solution by using a non-derivative local optimal algorithm Sbplx based on Nlopt; 5, judging whether the optimal solution needs to be updated or not according to a solving result; 6, judging whether iteration reaches a termination condition or not according to the current number of turns; step 7, according to the maneuvering type of the maneuvering satellite, judging whether the algorithm achieves a termination condition and whether additional maneuvering is needed; the method provided by the invention can meet the rapid construction of various approaching maneuvering data sets.

Description

Technical Field

[0001] The present invention belongs to the technical field of spacecraft orbit design and optimization, and in particular relates to a method for constructing a spacecraft approach maneuvering data set. Background Art

[0002] Today's space situation is becoming increasingly complex, and the sources of threats faced by spacecraft in daily orbit are increasing and the degree of threats is increasing. High-orbit satellites such as synchronous belt satellites are of great strategic value. How to ensure the safety of our high-value business satellites in orbit has become a problem that needs to be solved urgently. However, the real orbit data of maneuvering satellites with obvious intentions are very limited. Constructing similar maneuvering data sets through maneuvering decision-making methods will be an effective measure to study spacecraft protection issues, and it will also be of guiding significance for the research on maneuvering detection and intention recognition of maneuvering satellites.

[0003] The current orbit optimization design theory system is mature, and some scholars have combined orbit design with intention recognition to conduct research, but existing research generally only stays at the motion level intention, that is, judging what kind of relative motion relationship the maneuvering satellite may form with our satellite. Statistical methods and intelligent methods can achieve a good recognition success rate in this level of research, and samples of this type of research can also be quickly generated through reverse recursive dynamics. However, there are currently few studies on mission-level intention recognition. One of the difficulties is how to correctly construct the orbital maneuvering data of the maneuvering satellite. This type of data cannot be simply realized by reverse recursion and random maneuvering methods, but requires simulating real scenes, starting from the perspective of the maneuvering satellite, and completing the entire mission-level intention implementation process with optimal maneuvers in the appropriate time window. Therefore, it is necessary to propose a method for constructing a spacecraft approach maneuvering data set that can provide a large amount of orbital data with randomness and in line with real scenes with high success rate and high efficiency. Summary of the invention

[0004] The purpose of the present invention is to provide a method for constructing a spacecraft approach maneuver data set. Based on the orbit optimization design method, the method can quickly generate maneuver satellite orbit maneuver data of four types of maneuvers, namely, rapid approach maneuver, plane-skimming maneuver, plane-circling maneuver and co-orbit lagging phase maneuver. The method supports a large amount of random generation and solves the problem of difficulty in constructing mission-level maneuver data sets in the prior art.

[0005] To achieve the above object, the present invention provides a method for constructing a spacecraft approach maneuver data set, comprising the following steps:

[0006] Step 1: determine the maneuver type of the maneuvering satellite and set the input of the maneuver generation algorithm;

[0007] Step 2: Calculate the maximum number of flight circles n of the target satellite n With task time t MT ;

[0008] Step 3: Calculate the task time t based on the current number of laps k MT , and solve for the initial guess;

[0009] Step 4: construct a nonlinear equation that characterizes the successful execution of the maneuver by the maneuvering satellite, and use the derivative-free local optimal algorithm Sbplx based on Nlopt to solve the converged solution;

[0010] Step 5: According to the solution result, determine whether the optimal solution needs to be updated;

[0011] Step 6: Based on the current number of laps, determine whether the iteration has reached the termination condition; compare the current number of laps k with the maximum number of flight laps n n If the current number of circles k is less than the maximum number of flight circles n n , then return to step 3; if the current number of circles k is equal to the maximum number of flight circles n n , then execute step 7;

[0012] Step 7: Based on the maneuver type of the maneuvering satellite, determine whether the algorithm has reached the termination condition and whether additional maneuvers are required.

[0013] Preferably, the maneuvering type of the maneuvering satellite in step 1 includes the following two cases:

[0014] Case 1: If a rapid approach maneuver and a swipe maneuver are generated, there is no need to set up a virtual spacecraft, and the six initial orbital numbers of the target satellite are directly used as algorithm inputs;

[0015] Case 2: If you need to generate orbiting maneuvers and co-orbital backward phase maneuvers, you need to generate a virtual spacecraft. By modifying the true anomaly of the six initial orbital elements of the target satellite, you can get the six initial orbital elements of the virtual spacecraft. The calculation expression of the true anomaly in the six initial orbital elements of the virtual spacecraft is as follows:

[0016]

[0017] Where a, e and f are the semi-major axis, eccentricity and true anomaly of the target satellite at the initial time, respectively; f v is the true anomaly angle of the virtual spacecraft, and Δr is the position offset, which needs to be set to a certain value when generating orbiting maneuvers and co-orbital lagging phase maneuvers.

[0018] Preferably, the input parameters of the maneuver generation algorithm include the maximum mission time T all , the flyby height d and the initial value of the maneuver ΔV 0 .

[0019] Preferably, in step 2, the maximum number of target satellite flight circles n is calculated. n With task time t MT The process is as follows:

[0020] S21. Calculate the position r and velocity v corresponding to the six numbers of the track. The calculation expression is as follows:

[0021]

[0022] Where a is the semi-major axis, e is the eccentricity, i is the orbital inclination, ω is the argument of perigee, Ω is the right ascension of the ascending node, and f is the true anomaly;

[0023] S22. Substitute the initial orbital elements of the maneuvering satellite and the target satellite into equations (2) and (3) to obtain their initial positions and velocities R. 1 , V 1 and R 2 , V 2 , and calculate the intersection vector L of the two satellite orbits 12 , the calculation expression is as follows:

[0024]

[0025] S23. Calculate the orbit pericenter direction of the target satellite as E 2 , the calculation expression is as follows:

[0026]

[0027] Among them, μ is the gravitational constant of the central celestial body;

[0028] S24, calculate the true anomaly angle θ corresponding to the intersection of the target satellite orbit and the initial orbit plane of the maneuvering satellite on the target satellite orbit mb , the calculation expression is as follows:

[0029]

[0030] S25, the eccentricity e of the six elements of the initial orbit of the target satellite 2 and true close angle f 2 Calculate the target satellite's anomaly angle ξ 2 , the calculation expression is as follows:

[0031]

[0032] S26, calculate the near-point angle ξ of the target satellite arriving at the flyby position mb , the calculation expression is as follows:

[0033]

[0034] S27, calculate the time ΔT taken by the target satellite to reach the flyby point from the initial position mb0 , the calculation expression is as follows:

[0035]

[0036] In the formula, a 2 represents the semi-major axis of the six elements of the initial orbit of the target satellite, and μ represents the gravitational constant of the central celestial body of the Earth;

[0037] S28. Calculate the maximum number of flight circles n of the target satellite n and task time t MT , the calculation expression is as follows:

[0038]

[0039] In the formula, floor means rounding down.

[0040] Preferably, in step 3, based on the current lap number k, if it is the first iteration, k = 1, otherwise k = k + 1, the task time t is calculated MT , and the process of solving the initial guess is as follows:

[0041] S31. According to equation (4), the intersection vector L of the two satellite orbits is obtained: 12 Then, the true anomaly angle in the six elements of the maneuvering satellite orbit is set to 0, and the position velocity R of the maneuvering satellite at the pericenter is obtained. 1e and V 1e , perform the first maneuver here to enter the phase adjustment orbit, and calculate the speed V after the maneuver 1e ′, the calculation expression is as follows:

[0042]

[0043] Where, ΔV 0 Represents the initial value of the maneuver, which is one of the input parameters of the maneuver generation algorithm;

[0044] S32. After the maneuver, the semi-major axis of the maneuvering satellite is a 12 , the eccentricity is e 12 ; Using the method of equations (4) to (5), the intersection vector L between the phase modulation orbit and the target satellite orbit is obtained 12 ′ and the direction of the pericenter of the phase modulation orbit E 2 ′, calculate the true anomaly angle θ of the flyby point on the phase modulation orbit mb ′, the calculation expression is as follows:

[0045]

[0046] S33, calculate the flying height r f , the calculation expression is as follows:

[0047]

[0048] Where d is the flyby altitude, θmb ′ represents the true anomaly angle of the flyby point on the phase modulation orbit, where θ mb,inv ′=θ mb '+π,θ mb,inv ′ represents the true anomaly angle when approaching in the opposite direction;

[0049] S34, the initial value of the computer momentum guess ΔV 0 ′, the calculation expression is as follows:

[0050]

[0051] in

[0052]

[0053] In the formula, r f Indicates the flyby altitude of the current solution.

[0054] Preferably, the process of constructing the nonlinear equation characterizing the successful execution of the maneuver by the maneuvering satellite in step 4 is as follows:

[0055] S41. According to the six initial orbital numbers of the maneuvering satellite and the target satellite, the intersection vector L of the two satellite orbits is solved by equation (4): 12 ; Calculate the anomaly angle ξ corresponding to the initial position of the target satellite 1,0 , the calculation expression is as follows:

[0056]

[0057] In the formula, f 1 It represents the true anomaly angle among the six elements of the initial orbit of the maneuvering satellite;

[0058] S42, calculate the time ΔT when the moving satellite reaches the pericenter from the initial position 1 , the calculation expression is as follows:

[0059]

[0060] S43. Calculate ΔT 4 ; The same process as (12) to (18) is used to obtain the true anomaly angle θ of the flyby point on the phase modulation orbit mb ′, the true anomaly angle θ when approaching in the opposite direction mb,inv ′ and the maneuvering amount ΔV of the maneuvering satellite at the pericenter 0 ′; let θ = θ mb,inv ′, ΔV=ΔV 0 ',make The tangential and radial velocities are obtained as follows:

[0061]

[0062] After the maneuvering satellite maneuvers at the pericenter, the speed changes to:

[0063]

[0064] The semi-major axis a of the maneuvering satellite at this time is calculated using the velocity, and the eccentricity e of the current maneuvering satellite is calculated from this. The calculation expression is as follows:

[0065]

[0066] After the maneuvering satellite moves to the pericenter on the phase modulation orbit, it enters the transition orbit; at the end of the transition orbit, it performs a second maneuver and enters the approach orbit, and approaches the target satellite at the end of the approach orbit; the pericenter angle ξ of the maneuvering point and the end point on the approach orbit is calculated. i1 and i2 , the calculation expression is as follows:

[0067]

[0068]

[0069] Among them, θ i1 =arccosθ,θ i2 =θ i1 +π, so the time ΔT for the maneuvering satellite from the transition orbit maneuvering point to the flyby point 4 for:

[0070]

[0071] S44. Calculate the running time ΔT on the transition orbit 3 , the expression is as follows:

[0072]

[0073] in

[0074]

[0075] S45. Calculate the time ΔT of the dynamic satellite in the phase modulation orbit 2 , the expression is as follows:

[0076] ΔT 2 =T 12 ·[(t MT -ΔT 1 -ΔT 3 -ΔT 4 ) / T 12 ] round (29);

[0077] Among them, round means rounding to the nearest integer. is the period of the initial orbit;

[0078] S46. Obtain a nonlinear equation, the expression is as follows:

[0079] ΔT 1 +ΔT 2 +ΔT 3 +ΔT 4 =t MT (30).

[0080] Preferably, in step 5, the expression for determining whether the optimal solution needs to be updated according to the solution result is as follows:

[0081]

[0082] Where, ΔV * is the final solution, and ΔV represents the solution calculated in the current number of turns.

[0083] Preferably, in step 7, the process of determining whether the algorithm has reached the termination condition and whether additional maneuvers are required is as follows according to the maneuver type of the maneuvering satellite:

[0084] When the maneuver type performed by the maneuvering satellite is a rapid approach maneuver and a flyby maneuver, no virtual spacecraft is generated in step 1. After the maneuvering satellite enters the approach orbit, it implements a subsequent rapid approach maneuver or a flyby maneuver on the target satellite.

[0085] When the maneuver type performed by the maneuvering satellite is an orbiting maneuver and a co-orbital backward phase maneuver, according to step 1, a virtual spacecraft needs to be generated in the scene as the target satellite of the maneuvering satellite; steps 2 to 6 guide the maneuvering satellite to approach the virtual spacecraft, and at the end of the approach orbit, an additional maneuver is performed to transfer the target satellite of the maneuvering satellite from the virtual spacecraft to the real target satellite; the details are as follows:

[0086] For orbiting maneuvers, it is first necessary to make maneuver corrections so that the maneuvering satellite and the real target satellite are running on the same orbit; if at the end of the approaching orbit, the semi-major axis of the real target satellite is a 2f , the eccentricity is e 2f , the orbital inclination is i 2f , the argument of perigee is ω 2f , the right ascension of the ascending node is Ω 2f , the maneuvering satellite is now at position r 1f , the speed is v 1f , then the true anomaly angle f corresponding to the position of the maneuvering satellite at this time on the actual target satellite orbit is 1f for:

[0087]

[0088] To change the orbit of the maneuvering satellite to the actual orbit of the target satellite, the maneuver Δv 1 for:

[0089]

[0090] in

[0091]

[0092] λ=ω 2f +f 1f (35);

[0093] The maneuvering satellite needs to form a flying relationship with the real target satellite, and also needs to add a pulse maneuver Δv in the direction perpendicular to the orbital plane. 2 :

[0094]

[0095] Among them, ε is a small random quantity; the additional maneuver Δv that a maneuvering satellite performing an orbiting maneuver needs to do when approaching the end of the orbit is:

[0096] Δv=Δv 1 +Δv 2 (37);

[0097] For the co-orbital backward phase maneuver, it is only necessary to change the orbit of the maneuvering satellite to the actual orbit of the target satellite. The additional maneuver performed when approaching the end of the orbit is Δv 1 .

[0098] Therefore, the present invention adopts the above-mentioned method for constructing a spacecraft approach maneuvering data set, which has the following beneficial effects:

[0099] (1) The present invention discloses a method for constructing a spacecraft approach maneuvering dataset, which uses a nonlinear equation to characterize the terminal time error, and uses the derivative-free local optimal algorithm Sbplx based on Nlopt to quickly obtain the convergence solution of the equation group, thereby ensuring a high success rate of maneuvering generation. The proposed virtual spacecraft setting method and additional maneuvering processing of the terminal point can meet the needs of rapid construction of various approach maneuvering datasets.

[0100] (2) The present invention discloses a method for constructing a spacecraft approach maneuver dataset, which uses pulse maneuvers to complete orbital maneuvers, which conforms to the real scenario of spacecraft orbital maneuvers. The method uses a time window search strategy to find the most fuel-saving maneuvering points on the initial orbit and the transition orbit, which can achieve the purpose of saving fuel as much as possible in the complex construction process of multiple maneuvers, and is closer to the real maneuvering situation;

[0101] (3) The present invention discloses a method for constructing a spacecraft approach maneuver dataset, which uses the initial position of the spacecraft as input and calculates the maneuvering strategy of the maneuvering satellite from finding a time window to entering a transition orbit, approaching orbit, and the final mission implementation state after the mission begins. While conforming to the actual scenario of the maneuvering target from lurking to approaching, it has a strong randomness and is very suitable for the dataset construction work of spacecraft approach maneuver research.

[0102] The technical solution of the present invention is further described in detail below through the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0103] Figure 1 It is an overall flow chart of a method for constructing a spacecraft approach maneuvering data set according to the present invention;

[0104] Figure 2 A schematic diagram of the space trajectory of a rapid approach maneuver according to an embodiment of the present invention;

[0105] Figure 3 A schematic diagram of the space trajectory of a skimming maneuver according to an embodiment of the present invention;

[0106] Figure 4 A schematic diagram of a space trajectory of an aircraft orbiting an embodiment of the present invention;

[0107] Figure 5 A schematic diagram of the spatial trajectory of the co-orbital lagging phase maneuver according to an embodiment of the present invention;

[0108] Figure 6 It is a curve diagram of the relative distance between the maneuvering satellite and the target satellite during the execution of various maneuvering types in an embodiment of the present invention. DETAILED DESCRIPTION

[0109] The following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention claimed for protection, but merely represents selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0110] See also Figure 1-Figure 6 , a method for constructing a spacecraft approach maneuvering data set, comprising the following steps:

[0111] Step 1: Determine the maneuver type of the maneuver satellite and set the input of the maneuver generation algorithm; wherein the maneuver type of the maneuver satellite includes the following two cases:

[0112] Case 1: If a rapid approach maneuver and a swipe maneuver are generated, there is no need to set up a virtual spacecraft, and the six initial orbital numbers of the target satellite are directly used as algorithm inputs; in this example, the six initial orbital numbers of the target satellite are [29292.6362, 0.0100, 0.4700, 0.3491, 1.2352, 1.5776], and the six initial orbital numbers of the maneuvering satellite are [38407.1726, 0.0100, 0.4509, 0.3491, 5.1163, 1.5301], in order of semi-major axis, eccentricity, orbit inclination, argument of perigee, right ascension of ascending node, and true anomaly, with length units in km and angle units in rad.

[0113] Case 2: If you need to generate orbiting maneuvers and co-orbital backward phase maneuvers, you need to generate a virtual spacecraft. By modifying the true anomaly of the six initial orbital elements of the target satellite, you can get the six initial orbital elements of the virtual spacecraft. The calculation expression of the true anomaly in the six initial orbital elements of the virtual spacecraft is as follows:

[0114]

[0115] Where a, e and f are the semi-major axis, eccentricity and true anomaly of the target satellite at the initial time, respectively; f v is the true anomaly angle of the virtual spacecraft, Δr is the position offset, which needs to be set to a certain value when generating orbiting maneuvers and co-orbital phase-lag maneuvers. In this example, Δr is set to 100 for orbiting maneuvers and Δr is set to 1000 for co-orbital phase-lag maneuvers, in km.

[0116] The input parameters of the maneuver generation algorithm include the maximum mission time T all , the flyby height d and the initial value of the maneuver ΔV 0 In this example, the maximum task time T all The uniform setting is 72, in h. Except for the flyover height d of the flyover maneuver, which is set to 30, the flyover height d of other maneuver types is uniformly set to 0, in km. Maneuver initial value ΔV 0 The uniform setting is 0.015km / s.

[0117] Step 1 is mainly used to initialize the algorithm, clarify the input and constraints of the maneuverable construction, and complete the modeling of the track optimization design problem.

[0118] Step 2: Calculate the maximum number of flight circles n of the target satellite n With task time t MT ; The specific process is as follows:

[0119] S21. Calculate the position r and velocity v corresponding to the six numbers of the track. The calculation expression is as follows:

[0120]

[0121] Where a is the semi-major axis, e is the eccentricity, i is the orbital inclination, ω is the argument of perigee, Ω is the right ascension of the ascending node, and f is the true anomaly;

[0122] S22. Substitute the initial orbital elements of the maneuvering satellite and the target satellite into equations (2) and (3) to obtain their initial positions and velocities R. 1 , V 1 and R 2 , V 2 , and calculate the intersection vector L of the two satellite orbits 12 , the calculation expression is as follows:

[0123]

[0124] S23. Calculate the orbit pericenter direction of the target satellite as E 2 , the calculation expression is as follows:

[0125]

[0126] Among them, μ is the gravitational constant of the central celestial body;

[0127] S24, calculate the true anomaly angle θ corresponding to the intersection of the target satellite orbit and the initial orbit plane of the maneuvering satellite on the target satellite orbit mb , the calculation expression is as follows:

[0128]

[0129] S25, the eccentricity e of the six elements of the initial orbit of the target satellite 2 and true close angle f 2 Calculate the target satellite's anomaly angle ξ 2 , the calculation expression is as follows:

[0130]

[0131] S26, calculate the near-point angle ξ of the target satellite arriving at the flyby position mb , the calculation expression is as follows:

[0132]

[0133] S27, calculate the time ΔT taken by the target satellite to reach the flyby point from the initial position mb0 , the calculation expression is as follows:

[0134]

[0135] In the formula, a 2represents the semi-major axis of the six elements of the initial orbit of the target satellite, and μ represents the gravitational constant of the central celestial body of the Earth;

[0136] S28. Calculate the maximum number of flight circles n of the target satellite n and task time t MT , the calculation expression is as follows:

[0137]

[0138] In the formula, floor means rounding down.

[0139] Given the input in this example, the maximum number of flight circles n can be obtained n is 4, and the corresponding task time is t MT They are 87814.9978, 137709.0651, 187603.1325, 237497.1998 respectively, in seconds.

[0140] Step 2 mainly determines the initial positions of the maneuvering satellite and the target satellite based on the given maximum mission time T all , the time window for the maneuvering satellite to wait for maneuvers during the execution of the maneuver is obtained, and the number of cycles for subsequent iterative solutions is determined.

[0141] Step 3: Calculate the task time t based on the current number of laps k MT , and solve the initial guess; where, based on the current lap number k, if it is the first iteration, k = 1, otherwise k = k + 1, calculate the task time t MT , and the process of solving the initial guess is as follows:

[0142] S31. According to equation (4), the intersection vector L of the two satellite orbits is obtained: 12 Then, the true anomaly angle in the six elements of the maneuvering satellite orbit is set to 0, and the position velocity R of the maneuvering satellite at the pericenter is obtained. 1e and V 1e , perform the first maneuver here to enter the phase adjustment orbit, and calculate the speed V after the maneuver 1e ′, the calculation expression is as follows:

[0143]

[0144] Where, ΔV 0 Represents the initial value of the maneuver, which is one of the input parameters of the maneuver generation algorithm;

[0145] S32. After the maneuver, the semi-major axis of the maneuvering satellite is a 12 , the eccentricity is e 12 ; Using the method of equations (4) to (5), the intersection vector L between the phase modulation orbit and the target satellite orbit is obtained 12′ and the direction of the pericenter of the phase modulation orbit E 2 ′, calculate the true anomaly angle θ of the flyby point on the phase modulation orbit mb ′, the calculation expression is as follows:

[0146]

[0147] S33, calculate the flying height r f , the calculation expression is as follows:

[0148]

[0149] Where d is the flyby altitude, θ mb ′ represents the true anomaly angle of the flyby point on the phase modulation orbit, where θ mb,inv ′=θ mb ′+π,θ mb,inv ′ represents the true anomaly angle when approaching in the opposite direction;

[0150] S34, the initial value of the computer momentum guess ΔV 0 ′, the calculation expression is as follows:

[0151]

[0152] in

[0153]

[0154] In the formula, r f Indicates the flyby altitude of the current solution.

[0155] Step 3 assumes that the maneuvering satellite enters the phase-adjusting orbit at the pericenter of the initial orbit and solves the reasonable maneuver size of the first maneuver as the initial value guess for solving the subsequent nonlinear equations.

[0156] Step 4: Construct a nonlinear equation to characterize the successful execution of the maneuver by the maneuvering satellite, and use the derivative-free local optimal algorithm Sbplx based on Nlopt to find the converged solution; this equation characterizes the time from the initial moment to the time when the maneuvering satellite approaches the target satellite. When this time is equal to the mission time t MT When the difference is equal to zero, it means that the problem has a solution. To calculate the time of the maneuvering satellite, four time periods need to be calculated: the time ΔT from the initial position to the pericenter 1 , the time ΔT of the maneuverable satellite in the phase modulation orbit 2 , the time ΔT that the maneuverable satellite moves in the transition orbit 3 , the time ΔT of the maneuverable satellite in the approaching orbit 4 ; The process of constructing the nonlinear equation that characterizes the successful execution of the maneuver by the maneuvering satellite is as follows:

[0157] S41. According to the six initial orbital numbers of the maneuvering satellite and the target satellite, the intersection vector L of the two satellite orbits is solved by equation (4): 12 ; Calculate the anomaly angle ξ corresponding to the initial position of the target satellite 1,0 , the calculation expression is as follows:

[0158]

[0159] In the formula, f 1 It represents the true anomaly angle among the six elements of the initial orbit of the maneuvering satellite;

[0160] S42, calculate the time ΔT when the moving satellite reaches the pericenter from the initial position 1 , the calculation expression is as follows:

[0161]

[0162] S43. Calculate ΔT 4 ; The same process as (12) to (18) is used to obtain the true anomaly angle θ of the flyby point on the phase modulation orbit mb ′, the true anomaly angle θ when approaching in the opposite direction mb,inv ′ and the maneuvering amount ΔV of the maneuvering satellite at the pericenter 0 ′; let θ = θ mb,inv ′, ΔV=ΔV 0 ',make The tangential and radial velocities are obtained as follows:

[0163]

[0164] After the maneuvering satellite maneuvers at the pericenter, the speed changes to:

[0165]

[0166] The semi-major axis a of the maneuvering satellite at this time is calculated using the velocity, and the eccentricity e of the current maneuvering satellite is calculated from this. The calculation expression is as follows:

[0167]

[0168] After the maneuvering satellite moves to the pericenter on the phase modulation orbit, it enters the transition orbit; at the end of the transition orbit, it performs a second maneuver and enters the approach orbit, and approaches the target satellite at the end of the approach orbit; the pericenter angle ξ of the maneuvering point and the end point on the approach orbit is calculated. i1 and i2 , the calculation expression is as follows:

[0169]

[0170] Among them, θ i1=arccosθ,θ i2 =θ i1 +π, so the time ΔT for the maneuvering satellite from the transition orbit maneuvering point to the flyby point 4 for:

[0171]

[0172] S44. Calculate the running time ΔT on the transition orbit 3 , the expression is as follows:

[0173]

[0174] in

[0175]

[0176] S45. Calculate the time ΔT of the dynamic satellite in the phase modulation orbit 2 , the expression is as follows:

[0177] ΔT 2 =T 12 ·[(t MT -ΔT 1 -ΔT 3 -ΔT 4 ) / T 12 ] round (29);

[0178] Among them, round means rounding to the nearest integer. is the period of the initial orbit;

[0179] S46. Obtain a nonlinear equation, the expression is as follows:

[0180] ΔT 1 +ΔT 2 +ΔT 3 +ΔT 4 =t MT (30).

[0181] Based on the initial value given in step 3, the solution of the nonlinear equation can be found through the optimization algorithm. Since the gradient expression is not given, the local derivative-free optimization algorithm sbplx based on Nlopt is selected here to solve the nonlinear equation.

[0182] Step 4 calculates the maneuvering process of the maneuvering satellite, establishes the equation relationship between the flight time and the mission time, and finds the appropriate maneuvering amount ΔV by solving the nonlinear equation.

[0183] Step 5: Based on the solution, determine whether the optimal solution needs to be updated; the nonlinear equation obtained above may not have a solution; when the equation has no solution, go directly to step 6. When the equation has a solution, continue to execute the process of step 5. The specific expression for determining whether the optimal solution needs to be updated is as follows:

[0184]

[0185] Where, ΔV * is the final solution, and ΔV represents the solution calculated in the current number of turns.

[0186] When there is a solution in the first iteration, the solution obtained is the final solution. In subsequent judgments, if a solution is obtained that is better than the current ΔV * The smaller ΔV is, the final solution is updated to ensure that the maneuvering satellite performs the maneuver with the lowest fuel consumption.

[0187] Step 5 is mainly used to retain the solution with the lowest fuel consumption during multiple iterations, so that the maneuvering process of the maneuvering satellite is closer to the actual engineering application.

[0188] Step 6: Based on the current number of laps, determine whether the iteration has reached the termination condition; compare the current number of laps k with the maximum number of flight laps n n If the current number of circles k is less than the maximum number of flight circles n n , then return to step 3; if the current number of circles k is equal to the maximum number of flight circles n n , then execute step 7;

[0189] Step 6 ensures the smooth connection of the loop and ensures that the loop ends at the appropriate time when the loop enters the maximum number of flight circles. In this example, in the process of solving the rapid approach maneuver, the plane-skimming maneuver, the plane-circling maneuver and the co-orbital lagging phase maneuver, the solution of the nonlinear equation at different flight circles is shown in Table 1. Since the input is consistent and the difference between the various maneuver types is not reflected in the nonlinear equation solution process, the solution results of each maneuver type are consistent. The first and second circles did not converge, and the third and fourth circles had solutions. Finally, the solution with the minimum fuel consumption was selected to obtain ΔV * =0.0150km / s.

[0190] Table 1 Nonlinear equation solutions for various maneuvers with different numbers of turns

[0191]

[0192] Step 7: According to the maneuver type of the maneuvering satellite, determine whether the algorithm has reached the termination condition and whether additional maneuvers are required. The specific process is as follows:

[0193] When the maneuver type performed by the maneuvering satellite is a rapid approach maneuver or a flyby maneuver, no virtual spacecraft is generated in step 1. After the maneuvering satellite enters the approach orbit, it can subsequently perform a rapid approach maneuver or a flyby maneuver on the target satellite. When performing a rapid approach maneuver or a flyby maneuver, the space trajectory of the maneuvering satellite is as follows: Figure 2 and Figure 3 As shown, after the maneuvering satellite moves to the pericenter of the initial orbit, it performs the first maneuver to enter the phase adjustment orbit and the transition orbit, and after waiting for a suitable number of circles, it performs the second maneuver to enter the approach orbit;

[0194] When the maneuver type performed by the maneuvering satellite is an orbiting maneuver and a co-orbital backward phase maneuver, according to step 1, a virtual spacecraft needs to be generated in the scene as the target satellite of the maneuvering satellite; steps 2 to 6 guide the maneuvering satellite to approach the virtual spacecraft, and at the end of the approach orbit, an additional maneuver is performed to transfer the target satellite of the maneuvering satellite from the virtual spacecraft to the real target satellite; the details are as follows:

[0195] For orbiting maneuvers, it is first necessary to make maneuver corrections so that the maneuvering satellite and the real target satellite are running on the same orbit; if at the end of the approaching orbit, the semi-major axis of the real target satellite is a 2f , the eccentricity is e 2f , the orbital inclination is i 2f , the argument of perigee is ω 2f , the right ascension of the ascending node is Ω 2f , the position of the maneuvering satellite is now r 1f , the speed is v 1f , then the true anomaly angle f corresponding to the position of the maneuvering satellite at this time on the actual target satellite orbit is 1f for:

[0196]

[0197] To change the orbit of the maneuvering satellite to the actual orbit of the target satellite, the maneuver Δv 1 for:

[0198]

[0199] in

[0200]

[0201] λ=ω 2f +f 1f (35);

[0202] The maneuvering satellite needs to form a flying relationship with the real target satellite, and also needs to add a pulse maneuver Δv in the direction perpendicular to the orbital plane. 2 :

[0203]

[0204] Among them, ε is a small random quantity; the additional maneuver Δv that a maneuvering satellite performing an orbiting maneuver needs to do when approaching the end of the orbit is:

[0205] Δv=Δv 1 +Δv 2 (37);

[0206] For the co-orbital backward phase maneuver, it is only necessary to change the orbit of the maneuvering satellite to the actual orbit of the target satellite. The additional maneuver performed when approaching the end of the orbit is Δv 1 .

[0207] When performing orbiting maneuvers and co-orbital phase-lag maneuvers, the space trajectory of the maneuvering satellite is as follows: Figure 4 and Figure 5 As shown in the figure, after the maneuvering satellite reaches the pericenter of the initial orbit, it performs the first maneuver to enter the phase adjustment orbit and transition orbit, and after waiting for the appropriate number of circles, it performs the second maneuver to enter the approach orbit. At the end of the approach orbit, the maneuvering satellite performs the third maneuver, that is, the additional maneuver, and the orbit after the maneuver is the flyby orbit or the co-orbital lagging phase orbit.

[0208] Step 7 provides additional maneuvering methods for two types of maneuvers: orbiting maneuvers and co-orbital backward phase maneuvers, ensuring the integrity of maneuver execution.

[0209] During the four maneuvering processes, the relative distance curves between the maneuvering satellite and the target satellite are as follows: Figure 6 As shown. By observing the end of the relative distance curve, the difference between the four types of maneuvers can be known. The relative distance of the rapid approach maneuver drops rapidly at the end of the approach orbit and eventually reaches near 0 km, indicating that the maneuvering satellite has achieved a close approach to the target satellite. The relative distance of the grazing maneuver reaches a minimum value of around 30 km at the end of the approach orbit and then rises, indicating that the maneuvering satellite has grazed the target satellite. The relative distances of the two types of maneuvers, the orbiting maneuver and the co-orbital lagging phase maneuver, finally remain at around 100 km and 1000 km, respectively, indicating that the maneuvering satellite and the target satellite have achieved a relatively stable relative motion relationship.

[0210] Therefore, the present invention adopts the above-mentioned method for constructing a spacecraft approach maneuver data set. First, it is determined whether to construct a virtual spacecraft according to the generated maneuver type. This process can flexibly apply the orbit optimization design algorithm of the approach maneuver to the solution of other maneuver types; then, based on the given mission time, the number of orbits of the maneuvering satellite and the actual total mission time are solved; then, the solution is performed within different numbers of orbits to ensure that the maneuver of the maneuvering satellite is the optimal solution with the lowest fuel consumption. Then, the initial value guess is solved according to the data of different numbers of circles, and the nonlinear equation is solved based on the local derivative-free optimization algorithm sbplx of NLopt to achieve rapid convergence of the terminal time error; finally, it is determined whether additional maneuvers need to be performed according to the generated maneuver type, and the generation of various types of maneuvers is achieved while keeping the algorithm structure simple, and the correctness of the maneuver generation is guaranteed.

[0211] Finally, it should be noted that the above embodiments are only used to illustrate the technical solution of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that they can still modify or replace the technical solution of the present invention with equivalents, and these modifications or equivalent replacements cannot cause the modified technical solution to deviate from the spirit and scope of the technical solution of the present invention.

Claims

1. A method for constructing a spacecraft approach maneuvering data set, characterized in that: The following steps are involved: Step 1: determine the maneuver type of the maneuvering satellite and set the input of the maneuver generation algorithm; Step 2: Calculate the maximum number of flight circles n of the target satellite n With task time t MT ; Step 3: Calculate the task time t based on the current number of laps k MT , and solve for the initial guess; Step 4: construct a nonlinear equation that characterizes the successful execution of the maneuver by the maneuvering satellite, and use the derivative-free local optimal algorithm Sbplx based on Nlopt to solve the converged solution; Step 5: According to the solution result, determine whether the optimal solution needs to be updated; Step 6: According to the current number of circles, determine whether the iteration has reached the termination condition; Compare the current number of laps k with the maximum number of flight laps n n If the current number of circles k is less than the maximum number of flight circles n n , then return to step 3; if the current number of circles k is equal to the maximum number of flight circles n n , then execute step 7; Step 7: Based on the maneuver type of the maneuvering satellite, determine whether the algorithm has reached the termination condition and whether additional maneuvers are required.

2. The method for constructing a spacecraft approach maneuvering data set according to claim 1, characterized in that: The maneuvering types of the maneuvering satellite in step 1 include the following two cases: Case 1: If a rapid approach maneuver and a swipe maneuver are generated, there is no need to set up a virtual spacecraft, and the six initial orbital numbers of the target satellite are directly used as algorithm inputs; Case 2: If you need to generate orbiting maneuvers and co-orbital backward phase maneuvers, you need to generate a virtual spacecraft. By modifying the true anomaly of the six initial orbital elements of the target satellite, you can get the six initial orbital elements of the virtual spacecraft. The calculation expression of the true anomaly in the six initial orbital elements of the virtual spacecraft is as follows: Where a, e and f are the semi-major axis, eccentricity and true anomaly of the target satellite at the initial time, respectively; f v is the true anomaly angle of the virtual spacecraft, and Δr is the position offset, which needs to be set to a certain value when generating orbiting maneuvers and co-orbital lagging phase maneuvers.

3. The method for constructing a spacecraft approach maneuvering data set according to claim 2, characterized in that: The input parameters of the maneuver generation algorithm include the maximum mission time T all , flyby altitude d and maneuver initial value ΔV0.

4. The method for constructing a spacecraft approach maneuvering data set according to claim 3, characterized in that: In step 2, calculate the maximum number of target satellite flight circles n n With task time t MT The process is as follows: S21. Calculate the position r and velocity v corresponding to the six numbers of the track. The calculation expression is as follows: Where a is the semi-major axis, e is the eccentricity, i is the orbital inclination, ω is the argument of perigee, Ω is the right ascension of the ascending node, and f is the true anomaly; S22, Substitute the six initial orbital elements of the maneuvering satellite and the target satellite into equations (2) and (3), obtain the initial positions and velocities R1, V1 and R2, V2 of the two, and calculate the intersection vector L of the two satellite orbits 12 , the calculation expression is as follows: S23. Calculate the orbit pericenter direction of the target satellite as E2. The calculation expression is as follows: Among them, μ is the gravitational constant of the central celestial body; S24, calculate the true anomaly angle θ corresponding to the intersection of the target satellite orbit and the initial orbit plane of the maneuvering satellite on the target satellite orbit mb , the calculation expression is as follows: S25, calculate the eccentric anomaly ξ2 of the target satellite from the eccentricity e2 and the true anomaly f2 in the six elements of the initial orbit of the target satellite, and the calculation expression is as follows: S26, calculate the near-point angle ξ of the target satellite arriving at the flyby position mb , the calculation expression is as follows: S27, calculate the time ΔT taken by the target satellite to reach the flyby point from the initial position mb0 , the calculation expression is as follows: In the formula, a2 represents the semi-major axis of the six elements of the initial orbit of the target satellite, and μ represents the gravitational constant of the central celestial body of the Earth; S28. Calculate the maximum number of flight circles n of the target satellite n and task time t MT , the calculation expression is as follows: In the formula, floor means rounding down.

5. The method for constructing a spacecraft approach maneuvering data set according to claim 4, characterized in that: In step 3, based on the current lap number k, if it is the first iteration, k = 1, otherwise k = k + 1, calculate the task time t MT , and the process of solving the initial guess is as follows: S31. According to equation (4), the intersection vector L of the two satellite orbits is obtained: 12 Then, the true anomaly angle in the six elements of the maneuvering satellite orbit is set to 0, and the position velocity R of the maneuvering satellite at the pericenter is obtained. 1e and V 1e , perform the first maneuver here to enter the phase adjustment orbit, and calculate the speed V after the maneuver 1e ′, the calculation expression is as follows: Where ΔV0 represents the initial value of the maneuver, which is one of the input parameters of the maneuver generation algorithm; S32. After the maneuver, the semi-major axis of the maneuvering satellite is a 12 , the eccentricity is e 12 ; Using the method of equations (4) to (5), the intersection vector L between the phase modulation orbit and the target satellite orbit is obtained 12 ′ and the pericenter direction E2′ of the phase modulation orbit, calculate the true anomaly angle θ of the flyby point on the phase modulation orbit mb ′, the calculation expression is as follows: S33, calculate the flying height r f , the calculation expression is as follows: Where d is the flyby altitude, θ mb ′ represents the true anomaly angle of the flyby point on the phase modulation orbit, where θ mb,inv ′=θ mb ′+π,θ mb,inv ′ represents the true anomaly angle when approaching in the opposite direction; S34. The initial value of the computer momentum is guessed ΔV0′, and the calculation expression is as follows: in In the formula, r f Indicates the flyby altitude of the current solution.

6. The method for constructing a spacecraft approach maneuvering data set according to claim 5, characterized in that: The process of constructing the nonlinear equation that characterizes the successful execution of the maneuvering satellite in step 4 is as follows: S41. According to the six initial orbital numbers of the maneuvering satellite and the target satellite, the intersection vector L of the two satellite orbits is solved by equation (4): 12 ; Calculate the anomaly angle ξ corresponding to the initial position of the target satellite 1,0 , the calculation expression is as follows: Where, f1 represents the true anomaly angle among the six elements of the initial orbit of the maneuvering satellite; S42. Calculate the time ΔT1 when the moving satellite reaches the pericenter from the initial position. The calculation expression is as follows: S43, calculate ΔT4; the same process as equations (12) to (18) to obtain the true anomaly angle θ of the flyby point on the phase modulation orbit mb ′, the true anomaly angle θ when approaching in the opposite direction mb,inv ′ and the maneuvering amount ΔV0′ of the maneuvering satellite at the pericenter; let θ = θ mb,inv ′, ΔV=ΔV0′, let The tangential and radial velocities are obtained as follows: After the maneuvering satellite maneuvers at the pericenter, the speed changes to: The semi-major axis a of the maneuvering satellite at this time is calculated using the velocity, and the eccentricity e of the current maneuvering satellite is calculated from this. The calculation expression is as follows: After the maneuverable satellite moves to the pericenter in the phase modulation orbit, it enters the transition orbit; At the end of the transition orbit, a second maneuver is performed to enter the approach orbit, and the approach to the target satellite is achieved at the end of the approach orbit; the approach angle ξ between the maneuver point and the end point on the approach orbit is calculated. i1 and i2 , the calculation expression is as follows: Among them, θ i1 =arccosθ,θ i2 =θ i1 +π, so the time ΔT4 from the maneuvering point of the transition orbit to the flyby point of the maneuvering satellite is: S44, calculate the running time ΔT3 on the transition orbit, the expression is as follows: in S45. Calculate the time ΔT2 that the dynamic satellite runs in the phase modulation orbit. The expression is as follows: ΔT2=T 12 ·[(t MT -ΔT1-ΔT3-ΔT4) / T 12 ] round (29); Among them, round means rounding to the nearest integer. is the period of the initial orbit; S46. Obtain a nonlinear equation, the expression is as follows: ΔT1+ΔT2+ΔT3+ΔT4=t MT (30).

7. The method for constructing a spacecraft approach maneuvering data set according to claim 6, characterized in that: In step 5, based on the solution result, the expression for determining whether the optimal solution needs to be updated is as follows: Where, ΔV * is the final solution, and ΔV represents the solution calculated in the current number of turns.

8. The method for constructing a spacecraft approach maneuvering data set according to claim 7, characterized in that: In step 7, the process of judging whether the algorithm has reached the termination condition and whether additional maneuvers are required is as follows based on the maneuver type of the maneuvering satellite: When the maneuver type performed by the maneuvering satellite is a rapid approach maneuver and a flyby maneuver, no virtual spacecraft is generated in step 1. After the maneuvering satellite enters the approach orbit, it implements a subsequent rapid approach maneuver or a flyby maneuver on the target satellite. When the maneuver type performed by the maneuvering satellite is an orbiting maneuver and a co-orbital backward phase maneuver, according to step 1, a virtual spacecraft needs to be generated in the scene as the target satellite of the maneuvering satellite; steps 2 to 6 guide the maneuvering satellite to approach the virtual spacecraft, and at the end of the approach orbit, an additional maneuver is performed to transfer the target satellite of the maneuvering satellite from the virtual spacecraft to the real target satellite; the details are as follows: For orbiting maneuvers, it is first necessary to make maneuver corrections so that the maneuvering satellite and the real target satellite are running on the same orbit; if at the end of the approaching orbit, the semi-major axis of the real target satellite is a 2f , the eccentricity is e 2f , the orbital inclination is i 2f , the argument of perigee is ω 2f , the right ascension of the ascending node is Ω 2f , the position of the maneuvering satellite is now r 1f , the speed is v 1f , then the true anomaly angle f corresponding to the position of the maneuvering satellite on the actual target satellite orbit at this time is 1f for: To change the orbit of the maneuvering satellite to the actual orbit of the target satellite, the maneuver Δv1 to be performed is: in λ=ω 2f +f 1f (35); The maneuvering satellite needs to form a flying relationship with the real target satellite, and also needs to add a pulse maneuver Δv2 in the direction perpendicular to the orbital plane: Among them, ε is a random quantity; the additional maneuver Δv that a maneuvering satellite performing an orbiting maneuver needs to do when approaching the end of the orbit is: Δv=Δv1+Δv2(37); For the co-orbital backward phase maneuver, it is only necessary to change the orbit of the maneuvering satellite to the actual orbit of the target satellite, and the additional maneuver performed when approaching the end of the orbit is Δv1.

Images