A spacecraft proximity maneuver data set construction method

Abstract

The application discloses a kind of spacecraft near approach maneuver data set construction method, belong to spacecraft orbit design and optimization technical field, including the following steps: step 1, determine the maneuver type of maneuvering satellite and set the input of machine generation algorithm;Step 2, calculate the maximum flight circle number of target satellite and mission time;Step 3, based on the current circle number, calculate mission time, and solve initial guess;Step 4, construct the nonlinear equation representing the successful execution of maneuver of maneuvering satellite, and solve convergent solution using the derivative-free local optimal algorithm Sbplx based on Nlopt;Step 5, according to the solution, judge whether the optimal solution needs to be updated;Step 6, according to the current circle number, judge whether the iteration reaches the termination condition;Step 7, according to the maneuver type of maneuvering satellite, judge whether the algorithm reaches the termination condition, and whether additional maneuver is needed;The method provided by the application can satisfy the rapid construction of multiple near approach maneuver data sets.

Description

Technical Field

[0001] This invention belongs to the field of spacecraft orbit design and optimization technology, and in particular relates to a method for constructing a spacecraft approach maneuver dataset. Background Technology

[0002] The current space situation is becoming increasingly complex, and the sources and severity of threats to the daily operation of spacecraft in orbit are constantly increasing. High-orbit satellites, such as geosynchronous satellites, are of great strategic value, and ensuring the safe operation of these high-value operational satellites in orbit has become an urgent problem to be solved. However, real-world orbital data of maneuvering satellites with clear intentions is very limited. Constructing a similar maneuvering dataset using maneuvering decision-making methods will be an effective measure for studying spacecraft protection issues, and will also provide guidance for research on maneuvering detection and intent recognition of maneuvering satellites.

[0003] While current orbit optimization design theory is mature, and some scholars combine orbit design with intent recognition in their research, existing studies generally only focus on motion-level intent, i.e., determining the relative motion relationship between maneuvering satellites and friendly satellites. Statistical and intelligent methods can achieve good recognition success rates at this level, and samples for this type of research can be quickly generated using backpropagation dynamics. However, research on mission-level intent recognition is currently limited. One of the challenges lies in correctly constructing orbital maneuver data for maneuvering satellites. This type of data cannot be easily obtained using backpropagation and random maneuvering methods; instead, it requires simulating real-world scenarios, considering the perspective of the maneuvering satellite, and completing the entire mission-level intent implementation process with optimal maneuvering within a suitable time window. Therefore, it is necessary to propose a method for constructing spacecraft approach maneuver datasets that can provide a large amount of orbital data with randomness and conforming to real-world scenarios with high success rate and efficiency. Summary of the Invention

[0004] The purpose of this invention is to provide a method for constructing a spacecraft approach maneuver dataset. Based on an orbit optimization design method, it can quickly generate maneuver satellite orbit maneuver data for four types of maneuvers: rapid approach maneuver, skimming maneuver, orbiting maneuver, and co-orbiting phase-lag maneuver. It supports large-scale random generation and solves the problem of difficulty in constructing mission-level maneuver datasets in existing technologies.

[0005] To achieve the above objectives, this invention provides a method for constructing a spacecraft approach maneuver dataset, comprising the following steps:

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

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

[0008] Step 3: Calculate the task time t based on the current lap number k. MT And solve the initial conjecture;

[0009] Step 4: Construct the nonlinear equations characterizing the successful maneuvering of the satellite, and use the Nlopt-based derivative-free local optimum algorithm Sbplx to solve for the convergent solution;

[0010] Step 5: Based on the solution results, determine whether the optimal solution needs to be updated;

[0011] Step 6: Based on the current lap count, determine if the iteration has reached the termination condition; compare the current lap count k with the maximum number of flight laps n. n The size; if the current lap number k is less than the maximum flight lap number n n If the current lap number k equals the maximum flight lap number n, then return to step 3; n Then proceed to step 7;

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

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

[0014] Scenario 1: If rapid approach maneuvers and skimming maneuvers are generated, there is no need to set up a virtual spacecraft; the initial six-digit number of the target satellite's orbit can be directly used as the algorithm input.

[0015] Scenario 2: If it is necessary to generate orbital maneuvers and phase-lag maneuvers in the same orbit, a virtual spacecraft needs to be generated. This is done by modifying the true anomaly angles of the initial orbital six elements of the target satellite to obtain the initial orbital six elements of the virtual spacecraft. The expression for calculating the true anomaly angles in the initial orbital six elements of the virtual spacecraft is as follows:

[0016]

[0017] Where a, e, and f are the semi-major axis, eccentricity, and true anomaly angle of the target satellite at the initial moment, respectively; f v It is the true perimeter angle of the virtual spacecraft, and Δr is the position offset. It needs to be set to a certain value when generating orbital maneuvers and phase-lag maneuvers in the same orbit.

[0018] Preferably, the input parameters of the motion generation algorithm include the maximum task time T. all The skimming altitude d and the initial maneuver value ΔV0.

[0019] Preferably, in step 2, the maximum number of orbits n of the target satellite 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 elements of the track. The calculation expression is as follows:

[0021]

[0022] In the formula, 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 six initial orbital roots of the maneuvering satellite and the target satellite into equations (2) and (3) to obtain their initial positions and velocities R1, V1 and R2, V2, and calculate the intersection vector L of the two satellite orbits. 12 The calculation expression is as follows:

[0024]

[0025] S23. Calculate the direction of the target satellite's orbital pericenter as E2, using the following expression:

[0026]

[0027] Where μ is the gravitational constant of the central body;

[0028] S24. Calculate the true anomaly angle θ corresponding to the intersection point of the target satellite's orbit and the initial orbital plane of the maneuvering satellite on the target satellite's orbit. mb The calculation expression is as follows:

[0029]

[0030] S25. Calculate the target satellite's anomaly angle ξ2 from the eccentricity e2 and true anomaly angle f2 in the six roots of the target satellite's initial orbit. The calculation expression is as follows:

[0031]

[0032] S26. Calculate the angle of approach ξ of the target satellite as it approaches its flyby position. mb The calculation expression is as follows:

[0033]

[0034] S27. Calculate the time ΔT taken for the target satellite to reach the flyby point from its initial position. mb0 The calculation expression is as follows:

[0035]

[0036] In the formula, a2 represents the semi-major axis of the initial orbit of the target satellite, and μ represents the gravitational constant of the Earth's central body.

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

[0038]

[0039] In the formula, floor represents 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, and the task time t is calculated. MT The process of solving the initial conjecture is as follows:

[0041] S31. The intersection vector L of the two satellite orbits is obtained by solving equation (4). 12 Then, by setting the true anomaly angle in the six-element number of the maneuvering satellite's orbit to 0, the position and velocity R of the maneuvering satellite at its pericenter are obtained. 1e and V 1e The first maneuver is performed here to enter the phasing track, and the computer sets the speed V after the maneuver. 1e The calculation expression is as follows:

[0042]

[0043] In the formula, ΔV0 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 methods in equations (4) to (5), the intersection vector L between the phasing orbit and the target satellite orbit is obtained. 12 Calculate the true anterior angle θ of the flyby point on the phasing track, relative to the direction E2 of the pericenter point of the phasing track. mb The calculation expression is as follows:

[0045]

[0046] S33, Calculate the flyby height r f The calculation expression is as follows:

[0047]

[0048] In the formula, d represents the skimming altitude, and θ mb θ' represents the true anterior angle of the flyby point on the phasing trajectory, where θ mb,inv ′=θ mb +π, θ mb,inv ′ represents the true anterior angle when approaching from the opposite direction;

[0049] S34. The initial value of the computer's momentum is guessed as ΔV0′, and the calculation expression is as follows:

[0050]

[0051] in

[0052]

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

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

[0055] S41. Based on the initial orbital root 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 anomalous angle ξ corresponding to the initial position of the target satellite. 1,0 The calculation expression is as follows:

[0056]

[0057] In the formula, f1 represents the true anomaly angle in the six roots of the initial orbit of the maneuvering satellite;

[0058] S42. The time ΔT1 for the computer-controlled satellite to reach its pericenter from its initial position is calculated as follows:

[0059]

[0060] S43. Calculate ΔT4; following the same process as in equations (12) to (18), obtain the true near-point angle θ of the flyby point on the phase-adjusting track. mb ′, the true anterior angle θ when approaching from the opposite direction mb,inv And the maneuvering amount ΔV0′ of the maneuvering satellite at the pericenter; let θ = θ mb,inv Let ', ΔV=ΔV0', let The tangential and radial velocities are obtained using the following expressions:

[0061]

[0062] After maneuvering at the pericenter, the speed of the maneuvering satellite changes as follows:

[0063]

[0064] The semi-major axis 'a' of the maneuvering satellite is calculated using the velocity magnitude. From this, the eccentricity 'e' of the current maneuvering satellite is calculated using the following expression:

[0065]

[0066] After maneuvering to its pericenter in the phasing orbit, the satellite enters a transition orbit. At the end of the transition orbit, it performs a second maneuver to enter a close approach orbit, where it approaches the target satellite. The apogee angle ξ between the maneuver point and the terminal point in the close approach orbit is calculated. i1 and ξ i2 The calculation expression is as follows:

[0067]

[0068]

[0069] Where, θ i1 =arccosθ,θ i2 =θ i1 +π, therefore the time ΔT4 for the maneuvering satellite to move from the transition orbit maneuver point to the flyby point is:

[0070]

[0071] S44. Calculate the travel time ΔT3 on the transition track, as shown in the following expression:

[0072]

[0073] in

[0074]

[0075] S45. The time ΔT2 for the computer-controlled satellite to operate in phasing orbit is expressed as follows:

[0076] ΔT2=T 12 ·[(t MT -ΔT1-ΔT3-ΔT4) / T 12 ] round (29);

[0077] Here, "round" indicates rounding to the nearest integer. The period of the initial orbit;

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

[0079] ΔT1+ΔT2+ΔT3+ΔT4=t MT (30)

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

[0081]

[0082] In the formula, ΔV *The final solution is represented by ΔV, which indicates the solution calculated for the current lap number.

[0083] Preferably, in step 7, the process of determining whether the algorithm has met the termination condition and whether additional maneuvers are needed based on the maneuver type of the maneuvering satellite is as follows:

[0084] When the maneuvering satellite performs a rapid approach maneuver or a skimming maneuver, no virtual spacecraft is generated in step 1. After the maneuvering satellite enters the approach orbit, it can then perform a rapid approach maneuver or a skimming maneuver on the target satellite.

[0085] When the maneuver performed by the maneuvering satellite is an orbital maneuver or a phase-lag maneuver within the same orbit, according to step 1, a virtual spacecraft needs to be generated in the scenario as the target satellite for the maneuvering satellite; steps 2-6 guide the maneuvering satellite to approach the virtual spacecraft, and at the end of the approach orbit, perform an additional maneuver to transfer the target satellite of the maneuvering satellite from the virtual spacecraft to the real target satellite; specifically as follows:

[0086] For orbiting satellites, the first step is to perform maneuvering corrections to ensure that the maneuvering satellite and the actual target satellite are in the same orbit. If, at the end of the approach orbit, the semi-major axis of the actual target satellite is a... 2f The eccentricity is e 2f The orbital inclination angle is i 2f The argument of perigee is ω 2f The right ascension of the ascending node is Ω. 2f The position of the mobile satellite at this time is r 1f The speed is v 1f Then the true anomaly angle f corresponding to the current position of the maneuvering satellite in 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 required maneuver Δv1 is:

[0089]

[0090] in

[0091]

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

[0093] The maneuvering satellite needs to establish an orbital relationship with the actual target satellite, and also needs to incorporate a pulse maneuver Δv2 in a direction perpendicular to the orbital plane:

[0094]

[0095] Where ε is a small random variable; the additional maneuver Δv required by the maneuvering satellite to perform orbital maneuvers at the end of its approach orbit is:

[0096] Δv=Δv1+Δv2(37);

[0097] For a phase-lag maneuver in the same orbit, it is only necessary to change the orbit of the maneuvering satellite to the current orbit of the actual target satellite, and the additional maneuver performed at the end of the approach orbit is Δv1.

[0098] Therefore, the above-mentioned method for constructing a spacecraft approach maneuver dataset has the following beneficial effects:

[0099] (1) The present invention discloses a method for constructing a spacecraft approach maneuver dataset. It uses a nonlinear equation to characterize the terminal time error and obtains the convergent solution of the equation system quickly through the derivativeless local optimal algorithm Sbplx based on Nlopt, which ensures a high success rate of maneuver generation. Through the proposed virtual spacecraft setting method and additional maneuver processing of the terminal point, it can meet the requirements for the rapid construction of various approach maneuver datasets.

[0100] (2) The present invention discloses a method for constructing a spacecraft approach maneuver dataset, which uses pulse maneuvering to complete orbital maneuvers, conforming to the real scenario of spacecraft orbital maneuvers. This method utilizes a time window search strategy to find the most fuel-efficient maneuvering points on the initial and transition orbits, achieving the goal of saving as much fuel as possible during the complex construction process of multiple maneuvers, and more closely resembling real maneuvering situations;

[0101] (3) The present invention discloses a method for constructing a spacecraft approach maneuver dataset. Taking the initial position of the spacecraft as input, it calculates the maneuver strategy of the maneuvering satellite from the start of the mission to the time window to enter the transition orbit, approach orbit and the final mission implementation state. It conforms to the real scenario of the maneuvering target from lurking to approach, while having strong randomness, which is very suitable for the dataset construction work of spacecraft approach maneuver research.

[0102] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description

[0103] Figure 1 This is an overall flowchart of a method for constructing a spacecraft approach maneuver dataset according to the present invention;

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

[0105] Figure 3This is a schematic diagram of the spatial trajectory of a grazing aircraft according to an embodiment of the present invention;

[0106] Figure 4 This is a schematic diagram of the spatial trajectory of an aircraft moving around an embodiment of the present invention;

[0107] Figure 5 This is a schematic diagram of the spatial trajectory of a co-track lagging phase maneuver according to an embodiment of the present invention;

[0108] Figure 6 This is a graph showing the relative distance between the maneuvering satellite and the target satellite during the execution of various maneuvering types in embodiments of the present invention. Detailed Implementation

[0109] The following detailed description of embodiments of the invention provided in the accompanying drawings is not intended to limit the scope of the claimed invention, but merely to illustrate selected embodiments of the invention. All other embodiments obtained by those skilled in the art based on the embodiments of the invention without inventive effort are within the scope of protection of the invention.

[0110] Please see Figures 1-6 A method for constructing a spacecraft approach maneuver dataset includes the following steps:

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

[0112] Scenario 1: If generating rapid approach maneuvers and skimming maneuvers, there is no need to set up virtual spacecraft; directly use the initial six roots of the target satellite's orbit as the algorithm input. In this example, the initial six roots of the target satellite's orbit are [29292.6362,0.0100,0.4700,0.3491,1.2352,1.5776], and the initial six roots of the maneuvering satellite's orbit are [38407.1726,0.0100,0.4509,0.3491,5.1163,1.5301], in the following order: semi-major axis, eccentricity, orbital inclination, argument of perigee, right ascension of the ascending node, and true anomaly. The length unit is km, and the angle unit is rad.

[0113] Scenario 2: If it is necessary to generate orbital maneuvers and phase-lag maneuvers in the same orbit, a virtual spacecraft needs to be generated. This is done by modifying the true anomaly angles of the initial orbital six elements of the target satellite to obtain the initial orbital six elements of the virtual spacecraft. The expression for calculating the true anomaly angles in the initial orbital six elements of the virtual spacecraft is as follows:

[0114]

[0115] Where a, e, and f are the semi-major axis, eccentricity, and true anomaly angle of the target satellite at the initial moment, respectively; fv This is the true perimeter angle of the virtual spacecraft, and Δr is the position offset. These values ​​need to be set when generating orbital maneuvers and phase-lag maneuvers in the same orbit. In this example, Δr for orbital maneuvers is set to 100, and Δr for phase-lag maneuvers in the same orbit is set to 1000, in km.

[0116] The input parameters for the motion generation algorithm include the maximum task time T. all The flyby altitude d and the initial maneuver value ΔV0. In this example, the maximum mission time T. all The uniform value is set to 72, in hours (h). Except for the ski-flying altitude d of the ski-flying aircraft maneuver, which is set to 30, the ski-flying altitude d of all other maneuvers is uniformly set to 0, in kilometers (km). The initial maneuver value ΔV0 is uniformly set to 0.015 km / s.

[0117] Step 1 is mainly used for algorithm initialization, which clarifies the input and constraints of the maneuver construction and completes the modeling of the trajectory optimization design problem.

[0118] Step 2: Calculate the maximum number of orbits n for 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 elements of the track. The calculation expression is as follows:

[0120]

[0121] In the formula, 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 six initial orbital roots of the maneuvering satellite and the target satellite into equations (2) and (3) to obtain their initial positions and velocities R1, V1 and R2, V2, and calculate the intersection vector L of the two satellite orbits. 12 The calculation expression is as follows:

[0123]

[0124] S23. Calculate the direction of the target satellite's orbital pericenter as E2, using the following expression:

[0125]

[0126] Where μ is the gravitational constant of the central body;

[0127] S24. Calculate the true anomaly angle θ corresponding to the intersection point of the target satellite's orbit and the initial orbital plane of the maneuvering satellite on the target satellite's orbit. mb The calculation expression is as follows:

[0128]

[0129] S25. Calculate the target satellite's anomaly angle ξ2 from the eccentricity e2 and true anomaly angle f2 in the six roots of the target satellite's initial orbit. The calculation expression is as follows:

[0130]

[0131] S26. Calculate the angle of approach ξ of the target satellite as it approaches its flyby position. mb The calculation expression is as follows:

[0132]

[0133] S27. Calculate the time ΔT taken for the target satellite to reach the flyby point from its initial position. mb0 The calculation expression is as follows:

[0134]

[0135] In the formula, a2 represents the semi-major axis of the initial orbit of the target satellite, and μ represents the gravitational constant of the Earth's central body.

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

[0137]

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

[0139] Given the input in this example, the maximum number of flight laps n can be calculated. n The value is 4, corresponding to the task time t. MT The values ​​are 87814.9978, 137709.0651, 187603.1325, and 237497.1998, respectively, in seconds.

[0140] Step 2 mainly involves determining the initial positions of the maneuvering satellite and the target satellite, based on the given maximum mission time T. all The time window during which the maneuvering satellite waits to maneuver during the maneuver is obtained by solving the problem, which determines the number of iterations for subsequent solutions.

[0141] Step 3: Calculate the task time t based on the current lap number 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 The process of solving the initial conjecture is as follows:

[0142] S31. The intersection vector L of the two satellite orbits is obtained by solving equation (4). 12 Then, by setting the true anomaly angle in the six-element number of the maneuvering satellite's orbit to 0, the position and velocity R of the maneuvering satellite at its pericenter are obtained. 1e and V 1e The first maneuver is performed here to enter the phasing track, and the computer sets the speed V after the maneuver. 1e The calculation expression is as follows:

[0143]

[0144] In the formula, ΔV0 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 methods in equations (4) to (5), the intersection vector L between the phasing orbit and the target satellite orbit is obtained. 12 Calculate the true anterior angle θ of the flyby point on the phasing track, relative to the direction E2 of the pericenter point of the phasing track. mb The calculation expression is as follows:

[0146]

[0147] S33, Calculate the flyby height r f The calculation expression is as follows:

[0148]

[0149] In the formula, d represents the skimming altitude, and θ mb θ' represents the true anterior angle of the flyby point on the phasing trajectory, where θ mb,inv ′=θ mb +π, θ mb,inv ′ represents the true anterior angle when approaching from the opposite direction;

[0150] S34. The initial value of the computer's momentum is guessed as ΔV0′, and the calculation expression is as follows:

[0151]

[0152] in

[0153]

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

[0155] Step 3 assumes that the maneuvering satellite maneuvers into a phasing orbit at the pericenter of the initial orbit. Solve for the reasonable maneuver size of the first maneuver, which will serve as an initial value guess for solving the subsequent nonlinear equations.

[0156] Step 4: Construct the nonlinear equation characterizing the successful maneuver of the maneuvering satellite, and solve for the convergent solution using the Nlopt-based derivative-free local optimum algorithm Sbplx. This equation characterizes the time from the initial moment to the successful approach of the target satellite. When this time is equal to the mission time t... MT When the difference equals zero, it indicates that the problem has a solution. To calculate the time it takes for the maneuvering satellite to travel, four time segments need to be calculated: ΔT1 (time for the maneuvering satellite to reach its pericenter from its initial position), ΔT2 (time for the maneuvering satellite to travel in the phasing orbit), ΔT3 (time for the maneuvering satellite to travel in the transition orbit), and ΔT4 (time for the maneuvering satellite to travel in the approach orbit). The process of constructing the nonlinear equations characterizing the successful execution of the maneuvering satellite is as follows:

[0157] S41. Based on the initial orbital root 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 anomalous angle ξ corresponding to the initial position of the target satellite. 1,0 The calculation expression is as follows:

[0158]

[0159] In the formula, f1 represents the true anomaly angle in the six roots of the initial orbit of the maneuvering satellite;

[0160] S42. The time ΔT1 for the computer-controlled satellite to reach its pericenter from its initial position is calculated as follows:

[0161]

[0162] S43. Calculate ΔT4; following the same process as in equations (12) to (18), obtain the true near-point angle θ of the flyby point on the phase-adjusting track. mb ′, the true anterior angle θ when approaching from the opposite direction mb,inv And the maneuvering amount ΔV0′ of the maneuvering satellite at the pericenter; let θ = θ mb,inv Let ', ΔV=ΔV0', let The tangential and radial velocities are obtained using the following expressions:

[0163]

[0164] After maneuvering at the pericenter, the speed of the maneuvering satellite changes as follows:

[0165]

[0166] The semi-major axis 'a' of the maneuvering satellite is calculated using the velocity magnitude. From this, the eccentricity 'e' of the current maneuvering satellite is calculated using the following expression:

[0167]

[0168] After maneuvering to its pericenter in the phasing orbit, the satellite enters a transition orbit. At the end of the transition orbit, it performs a second maneuver to enter a close approach orbit, where it approaches the target satellite. The apogee angle ξ between the maneuver point and the terminal point in the close approach orbit is calculated. i1 and ξ i2 The calculation expression is as follows:

[0169]

[0170] Where, θ i1 =arccosθ,θ i2 =θ i1 +π, therefore the time ΔT4 for the maneuvering satellite to move from the transition orbit maneuver point to the flyby point is:

[0171]

[0172] S44. Calculate the travel time ΔT3 on the transition track, as shown in the following expression:

[0173]

[0174] in

[0175]

[0176] S45. The time ΔT2 for the computer-controlled satellite to operate in phasing orbit is expressed as follows:

[0177] ΔT2=T 12 ·[(t MT -ΔT1-ΔT3-ΔT4) / T 12 ] round (29);

[0178] Here, "round" indicates rounding to the nearest integer. The period of the initial orbit;

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

[0180] ΔT1+ΔT2+ΔT3+ΔT4=t MT (30)

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

[0182] Step 4 establishes an equation between flight time and mission time by tracking the satellite's maneuvering process using a computer, and finds a suitable maneuvering amount ΔV by solving nonlinear equations.

[0183] Step 5: Based on the solution results, determine whether the optimal solution needs to be updated. The nonlinear equation obtained above may not have a solution. If the equation has no solution, proceed directly to Step 6. If the equation has a solution, continue executing the process in Step 5. The specific expression for determining whether the optimal solution needs to be updated is as follows:

[0184]

[0185] In the formula, ΔV * The final solution is represented by ΔV, which indicates the solution calculated for the current lap number.

[0186] When a solution is found in the first iteration, that solution is the final solution. In subsequent judgments, if a solution is obtained that is greater than the current ΔV... * A smaller ΔV updates the final solution, ensuring that the maneuvering satellite performs the maneuver with the lowest possible 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 satellite is closer to the actual engineering application.

[0188] Step 6: Based on the current lap count, determine if the iteration has reached the termination condition; compare the current lap count k with the maximum number of flight laps n. n The size; if the current number of laps k is less than the maximum number of flight laps n n If the current lap number k equals the maximum flight lap number n, then return to step 3; n Then proceed to step 7;

[0189] Step 6 ensures smooth loop transitions and guarantees that the loop ends at the appropriate time when it reaches the maximum number of flight laps. Table 1 shows the solution results of the nonlinear equations for rapid approach maneuvers, skimming maneuvers, orbiting maneuvers, and co-track lagging phase maneuvers at different numbers of flight laps in this example. Since the input is consistent and the differences between the various maneuver types are not reflected in the nonlinear equation solution process, the solution results for each maneuver type are consistent. The first and second laps did not converge, while the third and fourth laps had solutions. Finally, the solution with the minimum fuel consumption was selected, yielding ΔV. * = 0.0150km / s.

[0190] Table 1. Solution results of nonlinear equations for various maneuvers with different revolution counts.

[0191]

[0192] Step 7: Based on the maneuver type of the maneuvering satellite, determine whether the algorithm has met the termination condition and whether additional maneuvers are needed. The specific process is as follows:

[0193] When the maneuvering satellite performs 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 then perform a rapid approach maneuver or a flyby maneuver on the target satellite. During rapid approach maneuvers and flyby maneuvers, the spatial trajectory of the maneuvering satellite is as follows: Figure 2 and Figure 3 As shown, after the maneuvering satellite reaches the pericenter of its initial orbit, it performs its first maneuver to enter the phasing orbit and transition orbit. After waiting for the appropriate number of orbits, it performs its second maneuver to enter the approach orbit.

[0194] When the maneuver performed by the maneuvering satellite is an orbital maneuver or a phase-lag maneuver within the same orbit, according to step 1, a virtual spacecraft needs to be generated in the scenario as the target satellite for the maneuvering satellite; steps 2-6 guide the maneuvering satellite to approach the virtual spacecraft, and at the end of the approach orbit, perform an additional maneuver to transfer the target satellite of the maneuvering satellite from the virtual spacecraft to the real target satellite; specifically as follows:

[0195] For orbiting satellites, the first step is to perform maneuvering corrections to ensure that the maneuvering satellite and the actual target satellite are in the same orbit. If, at the end of the approach orbit, the semi-major axis of the actual target satellite is a... 2f The eccentricity is e 2f The orbital inclination angle is i 2f The argument of perigee is ω 2f The right ascension of the ascending node is Ω. 2f The position of the mobile satellite at this time is r 1f The speed is v 1f Then the true anomaly angle f corresponding to the current position of the maneuvering satellite in 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 Δv1 required is:

[0198]

[0199] in

[0200]

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

[0202] The maneuvering satellite needs to establish an orbital relationship with the actual target satellite, and also needs to incorporate a pulse maneuver Δv2 in a direction perpendicular to the orbital plane:

[0203]

[0204] Where ε is a small random quantity; the additional maneuver Δv required by the maneuvering satellite to perform orbital maneuvers at the end of its approach orbit is:

[0205] Δv=Δv1+Δv2(37);

[0206] For a phase-lag maneuver in the same orbit, it is only necessary to change the orbit of the maneuvering satellite to the current orbit of the actual target satellite, and the additional maneuver performed at the end of the approach orbit is Δv1.

[0207] When performing orbital maneuvers and phase-lag maneuvers in the same orbit, the spatial trajectory of the maneuvering satellite is as follows: Figure 4 and Figure 5 As shown, after the maneuvering satellite reaches the pericenter of its initial orbit, it performs its first maneuver to enter the phasing orbit and transition orbit. After waiting for the appropriate number of orbits, it performs a second maneuver to enter the approach orbit. At the end of the approach orbit, the maneuvering satellite performs a third maneuver, namely the additional maneuver. The orbit after the maneuver is either a flyby orbit or a phase-lag orbit.

[0208] Step 7 provides additional maneuvering methods for two maneuvering types: orbital maneuvering and phase-lag maneuvering on the same track, 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. The differences between the four maneuver types can be discerned by observing the ends of the relative distance curves. The relative distance of the rapid approach maneuver decreases rapidly at the end of the approach orbit, eventually reaching around 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 of around 30 km at the end of the approach orbit, and then increases, indicating that the maneuvering satellite has passed the target satellite. The relative distances of the orbital maneuver and the phase-lag maneuver ultimately remain 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, this invention employs the aforementioned method for constructing a spacecraft approach maneuver dataset. First, it determines whether to construct a virtual spacecraft based on the generated maneuver type. This process allows the orbit optimization design algorithm for approach maneuvers to be flexibly applied to solving other maneuver types. Then, based on a given mission time, it solves for the number of orbits the maneuvering satellite needs and the actual total mission time. Subsequently, it solves for different numbers of orbits to ensure that the maneuvering satellite's maneuver is the optimal solution with the lowest fuel consumption. Next, it solves for initial value guesses based on data from different orbit numbers, and uses the NLopt local derivative-free optimization algorithm sbplx to solve nonlinear equations, achieving rapid convergence of the end-time error. Finally, it determines whether additional maneuvers are needed based on the generated maneuver type, achieving the generation of multiple maneuver types while maintaining a simple algorithm structure and ensuring the correctness of maneuver generation.

[0211] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit them. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the technical solutions of the present invention, and these modifications or equivalent substitutions cannot cause the modified technical solutions to deviate from the spirit and scope of the technical solutions of the present invention.

Claims

1. A method for constructing a proximity maneuver data set for a spacecraft, the method comprising: Comprising the following steps: Step 1, determine the type of maneuver of the maneuverable satellite and set the input of the algorithm of the maneuver generation; Step 2, calculate the maximum flight circle number n of the target satellite n With the mission time t MT ; Step 3, calculate the task time t based on the current round k MT and solve the initial guess; Step 4, construct a nonlinear equation representing the successful execution of the maneuver of the maneuverable satellite, and solve the convergent solution by using the derivative-free local optimal algorithm Sbplx based on Nlopt; wherein the process of constructing the nonlinear equation representing the successful execution of the maneuver of the maneuverable satellite is as follows: S41、According to the initial orbit elements of the maneuverable satellite and the target satellite, the intersection vector L of the two satellite orbits is solved by formula (4) 12 ; the eccentric anomaly angle ξ corresponding to the initial position of the target satellite is calculated 1,0 , and the calculation expression is as follows: In the formula, f1 represents the true anomaly in the initial orbit six elements of the maneuverable satellite, R1 and V1 represent the initial position and velocity of the maneuverable satellite, and R2 and V2 represent the initial position and velocity of the target satellite; S42, calculate the time ΔT1 of the maneuverable satellite from the initial position to the pericenter, and the calculation expression is as follows: In the formula, μ represents the central celestial body gravitational constant of the earth; S43, calculate ΔT4; obtain true anomaly θ at the perigee of the flyby on the osculating orbit mb , true anomaly θ at the perigee when approaching in the opposite direction mb,inv , and the maneuvering amount ΔV0' of the maneuvering satellite at the perigee; let θ = θ mb,inv , ΔV = ΔV0', let obtain tangential and radial velocities v r1 , v t1 , the expressions are as follows: θ mb,inv ′ = θ mb ′ + π; where L 12 E2' represents the perigee direction of the phasing orbit, r f represents the flyby height of the current solution, a 12 and e 12 respectively represent the semi-major axis and eccentricity of the maneuvered satellite, After the maneuver of the maneuverable satellite at the pericenter, the velocity is changed to: The semi-major axis a of the maneuverable satellite at this time is calculated by using the velocity, and the eccentricity e of the current maneuverable satellite is calculated, and the calculation expression is as follows: The maneuverable satellite runs on the phase modulation orbit to the pericenter, enters the transition orbit, carries out the second maneuver at the end of the transition orbit, enters the approaching orbit, realizes the approaching to the target satellite at the end of the approaching orbit, and calculates the eccentric angle ξ of the maneuver point and the end point on the approaching orbit i1 and ξ i2 , and the calculation expression is as follows: where θ i1 = arccos θ, θ i2 = θ i1 + π, so the time ΔT4 for the maneuvering satellite to move from the transition orbit maneuver point to the flyby point is: S44, calculate the running time ΔT3 on the transition orbit, and the expression is as follows: Wherein S45, calculate the time ΔT2 of the maneuverable satellite running on the phase modulation orbit, and the expression is as follows: ΔT2 = T 12 • [(t MT - ΔT1- ΔT3- ΔT4) / T 12 ] round (29); wherein round denotes rounding to the nearest integer, is the period of the initial orbit; S46, obtain the nonlinear equation, and the expression is as follows: ΔT1+ ΔT2+ ΔT3+ ΔT4= t MT (30); Step 5, according to the solving result, judge whether the optimal solution needs to be updated; Step 6, judging whether the iteration reaches the termination condition according to the current number of circles; comparing the current number of circles k with the maximum number of flight circles n n ; if the current number of circles k is less than the maximum number of flight circles n n , returning to step 3; if the current number of circles k is equal to the maximum number of flight circles n n , executing step 7; Step 7, according to the type of maneuver of the maneuverable satellite, judge whether the algorithm reaches the termination condition and whether additional maneuver is needed.

2. The method of claim 1, wherein: The type of maneuver of the maneuverable satellite in step 1 includes the following two cases: Case one: if the fast approach maneuver and the skimming maneuver are generated, the virtual spacecraft does not need to be set, and the initial orbit six elements of the target satellite are directly taken as the input of the algorithm; Case two: if the fly-around maneuver and the same orbit trailing phase maneuver need to be generated, the virtual spacecraft needs to be generated, the initial orbit six elements of the virtual spacecraft are obtained by modifying the true anomaly of the initial orbit six elements of the target satellite, and the calculation expression of the true anomaly in the initial orbit six 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; f v is the true anomaly of the virtual spacecraft, and Δr is the position offset, which is set to a certain value when generating the orbiting maneuver and the same-orbit trailing phase maneuver.

3. The method of claim 2, wherein: The input parameters of the machine generation algorithm include the maximum mission time T all , the flight path angle d and the maneuver initial value AV0.

4. The method of claim 3, wherein: The maximum flight circle number n of the target satellite is calculated in step 2 n The mission time t MT The process is as follows: S21, calculate the position r and velocity v corresponding to the orbit six elements, and the calculation expression is as follows: In the formula, a is the semi-major axis, e is the eccentricity, i is the orbit inclination, ω is the perigee amplitude, Ω is the ascending node right ascension, and f is the true anomaly; S22, the initial orbit elements of the maneuverable satellite and the target satellite are substituted into the equations (2) and (3) to obtain the initial positions and velocities R1, V1 and R2, V2 of the two satellites, and the intersection vector L of the orbits of the two satellites is calculated 12 The calculation expression is as follows: S23, calculate the orbit pericenter direction of the target satellite as E2, and the calculation expression is as follows: Wherein, μ is the central celestial body gravitational constant; S24, calculate the true anomaly θ corresponding to the intersection point of the target satellite orbit and the maneuvering satellite initial orbit plane on the target satellite orbit mb The calculation expression is as follows: S25, calculate the pericenter angle ξ2 of the target satellite from the eccentricity e2 and the true anomaly f2 in the initial orbit six elements of the target satellite, and the calculation expression is as follows: S26, calculate the argument of perigee ξ of the target satellite at the flyby location mb The calculation expression is as follows: S27, calculating 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 in the initial orbit six elements of the target satellite, and μ represents the central celestial body gravitational constant of the earth; S28、Calculate the maximum flight number n of the target satellite n and the mission time t MT , the calculation expression is as follows: In the formula, floor represents the down rounding.

5. The method of claim 4, wherein: Step 3: Calculate the task time t based on the current cycle k, if it is the first iteration, k = 1, otherwise k = k + 1 MT The process of solving the initial guess is as follows: S31, the intersection vector L of two satellite orbits is solved according to formula (4) 12 Then, the true anomaly of the six elements of the maneuvering satellite orbit is set to 0, and the position and velocity R of the maneuvering satellite at the pericenter are obtained 1e and V 1e At this point, the first maneuver into the phase modulation orbit is performed, and the post-maneuver velocity V 1e ′ is calculated, and the calculation expression is as follows: In the formula, ΔV0 represents the initial value of the maneuver, which is one of the input parameters of the algorithm of the maneuver generation; S32, after the maneuver, the semi-major axis of the maneuvered satellite is a 12 , and the eccentricity is e 12 ; using the method of equations (4) to (5), the intersection vector L 12 ′ of the phase-adjusted orbit and the target satellite orbit is obtained, and the true anomaly angle θ mb ′ of the flyby point on the phase-adjusted orbit is calculated, and the calculation expression is as follows: S33, calculate the fly height r f The calculation expression is as follows: where d represents the grazing altitude, θ mb represents the true anomaly of the flyby point on the phase curve, where θ mb,inv ′ = θ mb ′ + π, θ mb,inv ′ represents the true anomaly at the counter-approach. S34, calculate the initial value guess ΔV0' of the maneuver, and the calculation expression is as follows: Wherein where r f represents the current solution's flight altitude.

6. The method of claim 5, wherein: In step 5, according to the solving result, whether the optimal solution needs to be updated is judged, and the expression is as follows: where ΔV * is the final solution, and ΔV represents the solution calculated for the current round.

7. The method of claim 6, wherein: In step 7, according to the type of maneuver of the maneuverable satellite, whether the algorithm reaches the termination condition and whether additional maneuver is needed is judged, and the process is as follows: When the type of the maneuver performed by the maneuvering satellite is the fast-approaching maneuver and the grazing maneuver, no virtual space vehicle is generated in step 1, and the maneuvering satellite enters the approaching orbit to implement the fast-approaching maneuver or the grazing maneuver on the target satellite subsequently; When the type of the maneuver performed by the maneuvering satellite is the grazing maneuver and the same-orbit trailing-phase maneuver, according to step 1, a virtual space vehicle needs to be generated in the scene as the target satellite of the maneuvering satellite; steps 2-6 guide the maneuvering satellite to implement the approaching on the virtual space vehicle, and at the end of the approaching orbit, an additional maneuver is performed to transfer the target satellite of the maneuvering satellite from the virtual space vehicle to the real target satellite; specifically as follows: For the plane maneuver, first, the maneuvering satellite needs to be corrected to the same orbit as the real target satellite. Assuming that at the end of the approach orbit, the real target satellite has semi-major axis a 2f , eccentricity e 2f , orbital inclination i 2f , argument of perigee ω 2f , longitude of ascending node Ω 2f , the maneuvering satellite is at position r 1f , and the velocity is v 1f , then the true anomaly f 1f corresponding to the position of the maneuvering satellite on the orbit of the real target satellite is: To change the orbit of the maneuvering satellite into the orbit of the real target satellite, the maneuvering Δv1 needs to be performed: wherein λ = ω 2f +f 1f (35); The maneuvering satellite needs to form a grazing relationship with the real target satellite, and also needs to add a pulse maneuver Δv2 in the direction perpendicular to the orbital plane: wherein ε is a random quantity; the additional maneuver Δv performed by the maneuvering satellite at the end of the approaching orbit is: Δv = Δv1 + Δv2 (37); For the same-orbit trailing-phase maneuver, only the orbit of the maneuvering satellite needs to be changed into the current orbit of the real target satellite, and the additional maneuver performed at the end of the approaching orbit is Δv1.

Images