A fast calculation method for the reachable area of ​​space-based reentry vehicles considering no-fly zone constraints

Abstract

The present invention discloses a method for quickly calculating the reachable domain of a space-based re-entry vehicle taking into account the constraints of a no-fly zone, and belongs to the technical field of space-based re-entry vehicles. The method comprises: S100, first, using second-order cone programming to solve the transformed de-orbit braking convex optimization problem to obtain a single de-orbit optimal trajectory; S200, secondly, traversing and calculating the boundary of the reachable domain of the de-orbit segment; S300, then, converting the re-entry segment trajectory optimization problem into a convex optimization problem to obtain a single re-entry gliding segment optimal trajectory; S400, further changing the performance index to obtain the nearest point, farthest point, left boundary and right boundary re-entry segment reachable domain feature points; S500, finally, merging the reachable domain of the de-orbit segment and the re-entry gliding segment, and combining the analytical prediction method to calculate the reachable domain boundary of the re-entry gliding segment affected by the no-fly zone. Compared with the existing method, the present invention obtains a wider range of reachable domain, and at the same time, the speed of solving the reachable domain is faster, and has broad engineering application prospects.

Description

Technical Field

[0001] The present invention relates to a method for quickly calculating a reachable area of ​​a space-based re-entry vehicle taking into account no-fly zone constraints, and belongs to the technical field of space-based re-entry vehicles. Background Art

[0002] The reachable area of ​​a space-based re-entry glide vehicle refers to the distribution of landing points that the vehicle can reach. The reachable area allows analysis of the vehicle's flight capabilities and maneuverability during the vehicle's early design phase, providing key information for mission planning. It also serves as an important basis for selecting alternative landing points should the vehicle experience a sudden malfunction during flight and be unable to maintain nominal flight conditions. Calculation of the reachable area is crucial for analyzing the vehicle's flight capabilities and maneuverability during the vehicle's early design phase, providing key information for mission planning. It also serves as an important basis for selecting alternative landing points should the vehicle experience a sudden malfunction during flight and be unable to maintain nominal flight conditions. Analysis of a spacecraft's reachable area is crucial for completing various missions and provides a theoretical basis for various orbital optimization designs.

[0003] Given a set of initial conditions, a vehicle's reachable region is the set of feasible target points. A target point is feasible if the vehicle has sufficient capability to reach it without violating any constraints during flight. The reachable region provides important information for mission planning, so its calculation is crucial to developing more effective mission strategies during the mission planning phase. Traditional calculations of a vehicle's reachable region typically ignore maneuverability during the deorbit phase and consider only the reachable region within the atmosphere. Few studies consider the Earth's shape as an ellipsoid and the shape of the reachable region under different initial launches. Instead, the reachable region under different initial launches is approximated as the boundary of a fixed launch. Considering the impact of no-fly zones on the reachable region helps achieve more feasible mission planning results. Traditional reachable region calculation methods for boost-glide vehicles are no longer adaptable to the requirements of new missions. Therefore, to account for the impact of no-fly zones during mission planning, a fast calculation method for the reachable region of a space-based reentry vehicle that considers no-fly zone constraints is proposed, which has strong engineering practical significance. Summary of the Invention

[0004] The present invention proposes a method for rapidly calculating the reachable domain of a space-based re-entry vehicle taking into account the constraints of a no-fly zone. The reachable domains of the deorbit braking phase and the re-entry gliding phase are solved based on a convex optimization method. Based on the elliptical approximation method, a method for rapidly solving the reachable domain of the re-entry gliding phase is proposed. Furthermore, based on the analytical prediction range method, a method for correcting the reachable domain under the influence of a no-fly zone is given, in order to solve the problems existing in the prior art.

[0005] A method for rapidly calculating a reachable area of ​​a space-based reentry vehicle taking into account no-fly zone constraints is provided. The method comprises the following steps:

[0006] S100. First, a de-orbit trajectory planning method based on convex optimization is proposed. The de-orbit braking problem is transformed into a convex optimization problem, and the second-order cone programming method is used to solve it to obtain the optimal trajectory of a single de-orbit segment.

[0007] S200. Secondly, based on the optimal trajectory of a single off-track segment, the reachable domain is calculated by traversing the boundary. The boundary of the reachable domain is divided into N equal parts on one side. A three-dimensional off-track trajectory planning and reachable domain calculation process is designed to calculate the nearest point, farthest point, left boundary, and right boundary of the reachable domain.

[0008] S300, then, a re-entry glide trajectory planning method based on convex optimization is proposed, which transforms the re-entry trajectory optimization problem into a convex optimization problem, and obtains the optimal trajectory of a single re-entry glide segment;

[0009] S400, further changing the performance index, and obtaining the reentry segment reachable domain by traversing the boundary;

[0010] S500. Finally, the reachable domains of the deorbit segment and the re-entry glide segment are combined, and the reachable capability of the re-entry glide segment affected by the no-fly zone is calculated by combining the analytical prediction method, and finally the reachable domain boundary of the space-based re-entry glide vehicle considering the no-fly zone is obtained.

[0011] Furthermore, in S100, the following steps are specifically included:

[0012] S110, normalization of state equation:

[0013] definition is the thrust-acceleration vector, where T is the ratio of thrust to the initial weight of the engine, m is the ratio of the current mass of the aircraft to the initial mass, and z = lnm is defined. r and v are respectively expressed as R0 and Normalize and convert time t to Normalize and get the changed state equation,

[0014]

[0015] Where σ is a slack variable representing thrust-acceleration, and the state variable is x = (r, v, z) T , the control quantity is u=(τ,σ) T , write the state equation in matrix form,

[0016]

[0017] Among them, v ex =I sp g0, r = ‖r‖, the upper bound of engine thrust acceleration and the definition of relaxation relationship constraints are,

[0018] 0≤σ≤T max e -z (3)

[0019] ||τ||≤σ (4)

[0020] When using the iterative method to find the optimal solution, Equation (4) is regarded as a linear equation;

[0021] S120, solve the initial trajectory:

[0022] First, an initial trajectory is needed, and optimization is performed based on it. Here, the optimal fuel solution for reaching a fixed target point and specifying the time for the deorbit process is selected.

[0023] t f =t0 (5)

[0024] The deorbit process is divided into N segments, each of which is Solving N+1 state variables and control variables, the initial position, velocity, terminal position, and velocity constraints of the aircraft can be expressed as:

[0025]

[0026]

[0027] The deorbit process is constrained by equations (2) to (7). Since the distance from the center of the earth r is a normalized variable, it is approximately 1 during the deorbit process. When the accuracy of the result is not high, r in equation (2) is approximately 1, which becomes a linear equation. The fuel index function is:

[0028]

[0029] Since all constraints and objective functions are convex, this problem is a standard convex optimization problem;

[0030] S130, iteratively solving the time-optimal derailment braking problem:

[0031] Introducing the time variable, both Equation (2) and Equation (7) are nonlinear equality constraints. At the same time, multiple iterations are required based on the initial solution.

[0032] When the convex optimization problem contains nonlinear equality constraints h i When (y)=0,i=1,...,q, it is assumed that the current iteration result is y [k] , where y=(x,u,tf ), solve for y [k+1] The process is as follows:

[0033] In y [k] Taylor expansion of the function is performed at h i (y)=0 can be written as

[0034]

[0035] in, is a second-order correction that needs to be evaluated.

[0036] At the same time, there is an expansion type.

[0037]

[0038] in, For high-order corrections, only approximate solutions are needed It is approximately

[0039] The constraints of the convex optimization problem are,

[0040]

[0041] ||yy [k] ||≤ρ (12)

[0042] The optimization result is y p ,Will In y [k] First-order expansion

[0043]

[0044] Multiply both sides of the equation by (y p -y [k] ) T / 2, we get

[0045]

[0046] The left side of formula (13) is taken as The new constraint is,

[0047]

[0048] ||yy [k] ||≤ρ (15)

[0049] The non-convex constraint h i (y) = 0 is converted to a convex constraint,

[0050] In the time-optimal derailment braking problem, in addition to considering equations (2), (3), and (4), the remaining constraints are:

[0051]

[0052] z(N+1)>z f (17)

[0053] The indicator function is:

[0054] J=t (18)

[0055] The discrete form of formula (2) and formula (16) are nonlinear constraints. According to the above, we can obtain the constraint conditions of the convex function and add the following constraints:

[0056]

[0057] That is to find y [k+1] and u [k+1] , iterates multiple times until convergence, and obtains the time-optimal deorbit braking result.

[0058] Furthermore, in S200, specifically, the reachable domain is calculated by traversing the boundary. Assuming that the reachable domain boundary is divided into N equal parts on one side, the three-dimensional off-track trajectory planning and reachable domain calculation process is designed, including the following steps:

[0059] S210, calculate the closest point in the reachable domain: take the minimum range in the orbital plane as the optimization index, solve the trajectory with the shortest range, and record the minimum range R min ;

[0060] S220, calculate the farthest point in the reachable domain: take the maximum range in the orbital plane as the optimization index, solve the trajectory with the longest range, and record the maximum range R max ;

[0061] S230, calculate the left boundary: first simulate based on the finite thrust guidance algorithm to obtain the remaining fuel mass m of the aircraft r , m r Divide into N equal parts, and assume that the mass of the aircraft when the fuel is exhausted is m0. For the i-th (i=1,2...N-1) boundary point, the nearest boundary point is numbered i=1. It is approximately assumed that there are i·m r / N of fuel is used to change the flight direction and generate lateral motion. The aircraft is set to first perform a lateral braking angle component β = 90° with a constant positive sign. After this part of the fuel is consumed, longitudinal deorbit braking is performed. The longitudinal deorbit braking algorithm is the same as the algorithm in the previous section.

[0062] S240, calculate the right boundary: adopt the method of S230, set the sign of the lateral braking angle component β to be always negative, and obtain the right boundary trajectory.

[0063] Furthermore, in S300, specifically,

[0064] Considering the influence of the earth's rotation, the following reentry vehicle motion equation is established:

[0065]

[0066] During the re-entry process of the spacecraft, it is expressed as,

[0067]

[0068] in, q and n are heat flux, dynamic pressure and overload respectively, k Q is the heat flux coefficient, q max and n max are peak heat flux, peak dynamic pressure and overload peak, respectively.

[0069] So far, the optimal control problem is as follows:

[0070]

[0071] For nonlinear systems,

[0072]

[0073] There are bounds on state variables and control variables,

[0074]

[0075] Boundary constraints,

[0076]

[0077] Among them, ψ ol , ψ ou are the upper and lower bounds of ψ at the initial moment, ψ fl , ψ fu are the upper and lower bounds of the terminal time ψ,

[0078] Path constraints:

[0079] C l ≤C[x(τ),u(τ),τ]≤C u (26)

[0080] Among them, C l with C u are the upper and lower bounds of the path constraint,

[0081] Adopting a sequential convex optimization method based on trust region constraints, we gradually approach the solution of the original problem by solving a series of convex optimization sub-problems.

[0082] After linearization, we get

[0083]

[0084] The equation of state is converted to x k The linearized dynamic equation is

[0085]

[0086] The basic constraints in Equation (26) can also be processed by the above convexification method. The iterative optimization process of sequential convex optimization is linearized around a given reference trajectory. In order to ensure that the subproblem approaches the original problem, the trust region constraint || xx needs to be added. k ||≤Δ,

[0087] Adding a virtual control variable in equation (28), we get:

[0088]

[0089] During the iteration process, the range of the virtual control amount in the objective function is controlled so that it is small in the later stage of the iteration, thereby ensuring the accuracy of the problem, that is,

[0090] J′=J+δ||Cw||2 (30)

[0091] To solve the numerical solution of the infinite-dimensional optimal control problem, discretize it and divide the time N into equal parts. The discretized dynamic equation is:

[0092]

[0093] The superscript k represents the kth iteration, and the subscript i represents the t i time, Rearrange formula (31) to get:

[0094]

[0095] in

[0096]

[0097]

[0098] The above method transforms the reentry trajectory optimization problem into a convex optimization problem.

[0099] Furthermore, in S400, the following steps are specifically included:

[0100] S410, calculating the closest point of the axis range: using the closest vertical range of the initial shot as the performance index, the trajectory is calculated and recorded as L1;

[0101] S420, calculating the farthest point of the axis range: using the farthest range as a performance indicator, calculating the trajectory, which is recorded as L2;

[0102] S430, calculate the left and right boundaries: take the terminal ranges of L1 and L2 as the two endpoints, divide them into n equal parts, and evenly select n-1 points from small to large as the terminal range constraints. Optimize with the maximum and minimum horizontal ranges as performance indicators, and obtain n-1 left and right boundaries. The left boundary trajectory is recorded as L3~L n+1 , the right boundary trajectory is L n+2 ~L 2×n ;

[0103] S440, axis calculation: take the terminal ranges of L1 and L2 as the two endpoints, divide them into n equal parts, evenly select n-1 points from small to large as the terminal range constraints, and set the terminal lateral range constraint to 0. Optimize with flight time as the performance indicator, and you can get n-1 trajectories on the axis, which are L 2×n+1 ~L 3×n -1;

[0104] S450, near-boundary calculation: take the initial upward and leftward range as the performance index and calculate the trajectory L 3×n , taking the initial upward range to the right as the performance index, the trajectory L is calculated 3×n+1 , respectively L 3×n and L1 terminal range, L 3×n+1 The terminal range of L1 is divided into two endpoints, and m points are evenly divided from small to large. m-1 points are evenly selected as the terminal range constraints. The maximum and minimum values ​​of the horizontal range under the fixed vertical range on one side of the initial upward direction are calculated respectively to obtain the trajectory L 3×n+2 ~L 3×n+(m-1)×4 ;

[0105] An ellipse is used to approximate the boundaries of the reentry vehicle's reachable area. The calculation process of the four characteristic points of the ellipse includes:

[0106] S451. Using the minimum vertical range of the initial shot as the performance indicator, calculate the trajectory, which is recorded as L1.

[0107] J=min(S r ) (35)

[0108] S452. Using the maximum vertical range of the initial shot as the performance indicator, calculate the trajectory, which is recorded as L2.

[0109] J=min(-S r ) (36)

[0110] The formula for calculating the fixed range S by the maximum and minimum range points is as follows:

[0111] S=(S m +(S f -S m ) / 2)·c (37)

[0112] Where c is the empirical coefficient, and its value range is 0≤c≤1;

[0113] S453, fixed range S, calculate the trajectory with the maximum horizontal range as the performance index under range S, this trajectory is recorded as L3,

[0114] J=min(-C r ) (38)

[0115] -0.001≤D r -R≤0.001 (39)

[0116] S454, fixed range S, calculate the trajectory with the minimum horizontal range as the performance indicator under range S, this trajectory is recorded as L4,

[0117] J=min(C r ) (40)

[0118] -0.001≤D r -R≤0.001 (41) According to the ellipse calculation formula as shown in formula (42),

[0119]

[0120] The center e, major axis a, and minor axis b1 and b2 of the ellipse are obtained respectively.

[0121]

[0122] Furthermore, in S500, specifically:

[0123] In order to measure the reachability of the aircraft during the re-entry gliding phase under the influence of the no-fly zone, a quantitative calculation of the reachable area under the influence of the no-fly zone is performed based on the analytical method.

[0124] The Dubins curve is used to calculate the no-fly zone avoidance trajectory, assuming that the aircraft is in a balanced flight with no sideslip during the turning process and the speed angle is approximately zero:

[0125]

[0126] The component of lift produces centripetal force, assuming the turning radius is r z , we can get:

[0127]

[0128] Combining equations (44) and (45), ignoring the centrifugal force in the horizontal plane, we obtain an approximate analytical solution for the turning radius of the aircraft:

[0129]

[0130] Then maintain the turning radius r z The bank angle for balanced flight is:

[0131]

[0132] Assuming that the flight altitude before the turn is h1 and the atmospheric density is ρ1, the values ​​after the turn are h2 and ρ2, which are estimated by the atmospheric density fitting formula.

[0133]

[0134] In order to maintain balanced flight after turning, the atmospheric density needs to be satisfied.

[0135]

[0136] Get the relationship between the height change and the roll angle,

[0137]

[0138] The energy-based range prediction is obtained as,

[0139]

[0140] The range reduction caused by turning is calculated analytically using equation (51).

[0141] Beneficial effects of the present invention: The present invention proposes a method for rapidly calculating the reachable domain of a space-based re-entry vehicle taking into account the constraints of a no-fly zone, and solves the reachable domain of the three-dimensional deorbit braking segment and the re-entry gliding segment based on the convex optimization method. Compared with the existing reachable domain calculation method, the lateral maneuverability of the deorbit segment is taken into account, and the range of the reachable domain is broadened. A method for rapidly solving the reachable domain of the re-entry gliding segment is proposed based on the elliptical approximation method, and the range is predicted through an analytical method, and a method for correcting the reachable domain under the influence of a no-fly zone is given, ultimately achieving the multi-stage reachable domain rapid calculation requirement taking into account the constraints of a no-fly zone. Therefore, compared with the existing trajectory planning-based method, the present invention obtains a wider range of reachable domains, and at the same time, compared with the existing trajectory optimization-based method, the reachable domain solution speed is faster, and has broad prospects for engineering applications. BRIEF DESCRIPTION OF THE DRAWINGS

[0142] Figure 1 Schematic diagram of three-dimensional reachable domain;

[0143] Figure 2 Schematic diagram of the ellipse approximate reachable domain;

[0144] Figure 3 is the calculation result of the reachable domain of the derailment braking segment, where: Figure 3 (a) is the height change curve, Figure 3 (b) is the speed change curve, Figure 3 (c) is the velocity angle variation curve, Figure 3 (d) is the velocity angle variation curve, Figure 3 (e) is the three-dimensional trajectory, Figure 3 (f) Mass change curve, Figure 3 (g) is the longitudinal and transverse range change curve, Figure 3 (h) is the reachable domain of the off-orbit segment;

[0145] Figure 4 is the calculation result of the reachable area of ​​the reentry glide phase, where Figure 4 (a) is the three-dimensional trajectory, Figure 4 (b) is the speed change curve, Figure 4 (c) is the longitude and latitude variation curve, Figure 4 (d) is the ground track change curve, Figure 4 (e) is the latitude and longitude of the reachable domain boundary, Figure 4 (f) vertical and horizontal ranges that can reach the domain boundary;

[0146] Figure 5 Comparison chart of reachable domain solution results;

[0147] Figure 6 Calculated for the reachable area without no-fly zones;

[0148] Figure 7 is the relationship between turning distance and loss of range;

[0149] Figure 8 It is the boundary of the accessible area under the influence of the no-fly zone. DETAILED DESCRIPTION

[0150] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. All other embodiments obtained by ordinary technicians in this field based on the embodiments of the present invention without making any creative efforts shall fall within the scope of protection of the present invention.

[0151] Reference Figures 1-8As shown, the present invention proposes a method for quickly calculating the reachable domain of a space-based re-entry vehicle considering the constraints of a no-fly zone. The reachable domains of the deorbit braking phase and the re-entry gliding phase are solved based on a convex optimization method. Based on the elliptical approximation method, a method for quickly solving the reachable domain of the re-entry gliding phase is proposed. Based on the analytical prediction range method, a method for correcting the reachable domain under the influence of a no-fly zone is given to solve the problems existing in the prior art.

[0152] A method for rapidly calculating a reachable area of ​​a space-based reentry vehicle taking into account no-fly zone constraints is provided. The method comprises the following steps:

[0153] S100. First, a de-orbit trajectory planning method based on convex optimization is proposed. The de-orbit braking problem is transformed into a convex optimization problem, and the second-order cone programming method is used to solve it to obtain the optimal trajectory of a single de-orbit segment.

[0154] S200. Secondly, based on the optimal trajectory of a single off-track segment, the reachable domain is calculated by traversing the boundary. The boundary of the reachable domain is divided into N equal parts on one side. A three-dimensional off-track trajectory planning and reachable domain calculation process is designed to calculate the nearest point, farthest point, left boundary, and right boundary of the reachable domain.

[0155] S300, then, a re-entry glide trajectory planning method based on convex optimization is proposed, which transforms the re-entry trajectory optimization problem into a convex optimization problem, and obtains the optimal trajectory of a single re-entry glide segment;

[0156] S400, further changing the performance index, and obtaining the reentry segment reachable domain by traversing the boundary;

[0157] S500. Finally, the reachable domains of the deorbit segment and the re-entry glide segment are combined, and the reachable capability of the re-entry glide segment affected by the no-fly zone is calculated by combining the analytical prediction method, and finally the reachable domain boundary of the space-based re-entry glide vehicle considering the no-fly zone is obtained.

[0158] Furthermore, in S100, in the time-optimal deorbit braking problem, given the altitude and speed at the starting time, the aircraft is required to reach a certain altitude in the shortest time and meet the terminal speed and inclination requirements. The steps for solving the problem using convex optimization are as follows:

[0159] S110, normalization of state equation:

[0160] definition is the thrust-acceleration vector, where T is the ratio of thrust to the initial weight of the engine, m is the ratio of the current mass of the aircraft to the initial mass, and z = lnm is defined. r and v are respectively expressed as R0 and Normalize and convert time t to Normalize and get the changed state equation,

[0161]

[0162] Where σ is a slack variable representing thrust-acceleration, and the state variable is x = (r, v, z) T , the control quantity is u=(τ,σ) T , write the state equation in matrix form,

[0163]

[0164] Among them, v ex =I sp g0, r = ‖r‖, the upper bound of engine thrust acceleration and the definition of relaxation relationship constraints are,

[0165] 0≤σ≤T max e -z (3)

[0166] ||τ||≤σ (4)

[0167] When using the iterative method to find the optimal solution, Equation (2) is regarded as a linear equation;

[0168] S120, solve the initial trajectory:

[0169] When using the iterative method to solve the time-optimal deorbit braking problem, we first need an initial trajectory, based on which we optimize the solution. Here, we select the optimal fuel solution that reaches a fixed target point and specifies the deorbit process time.

[0170] t f =t0 (5)

[0171] To solve the fuel optimal deorbit problem, the deorbit process needs to be divided into N segments, each of which is Solving N+1 state variables and control variables, the initial position, velocity, terminal position, and velocity constraints of the aircraft can be expressed as:

[0172]

[0173]

[0174] The deorbit process is constrained by equations (2) to (7). Since the distance from the center of the earth r is a normalized variable, it is approximately 1 during the deorbit process. When the accuracy of the result is not high, r in equation (2) is approximately 1, which becomes a linear equation. The fuel index function is:

[0175]

[0176] Since all constraints and objective functions are convex, this problem is a standard convex optimization problem;

[0177] S130, iteratively solving the time-optimal derailment braking problem:

[0178] Solving the time-optimal deorbit braking problem does not specify a fixed reentry point or deorbit time. It only requires that the specified velocity and inclination angle be met when reaching a fixed reentry altitude. Therefore, it is necessary to introduce a time variable. In this case, both Equations (2) and (7) are nonlinear equality constraints. At the same time, to ensure the accuracy of the results, multiple iterative solutions are required based on the initial solution.

[0179] When the convex optimization problem contains nonlinear equality constraints h i When (y)=0,i=1,...,q, it is assumed that the current iteration result is y [k] , where y=(x,u,t f ), solve for y [k+1] The process is as follows:

[0180] In y [k] Taylor expansion of the function is performed at h i (y)=0 can be written as

[0181]

[0182] in, is a second-order correction that needs to be evaluated.

[0183] At the same time, there is an expansion type.

[0184]

[0185] in, For high-order corrections, only approximate solutions are needed It is approximately The constraints of the convex optimization problem are,

[0186]

[0187] ||yy [k] ||≤ρ (12)

[0188] The optimization result is y p ,Will In y [k] First-order expansion

[0189]

[0190] Multiply both sides of the equation by (y p -y [k] ) T / 2, we get

[0191]

[0192] The left side of formula (13) is taken as The new constraint is,

[0193]

[0194] ||yy [k] ||≤ρ (15)

[0195] Therefore, the non-convex constraint h i (y)=0 converts to a convex constraint.

[0196] In the time-optimal derailment braking problem, in addition to considering equations (2), (3), and (4), the remaining constraints are:

[0197]

[0198] z(N+1)>z f (17)

[0199] The indicator function is:

[0200] J=t (18)

[0201] The discrete form of formula (2) and formula (16) are nonlinear constraints. According to the above, we can obtain the constraint conditions of the convex function and add the following constraints:

[0202]

[0203] That is to find y [k+1] and u [k+1] , iterates multiple times until convergence, and obtains the time-optimal deorbit braking result.

[0204] Furthermore, in S200 , the remaining fuel is specifically used to accelerate the vehicle initially and decelerate later, thus meeting the required speed. Movement in the longitudinal plane dissipates the vehicle's excess mass. Similarly, lateral maneuvers can be performed, thereby expanding the vehicle's reach during the deorbit braking phase. In three-dimensional space, the braking motors no longer generate only longitudinal braking force, but can also generate lateral braking forces through lateral swing angles.

[0205] The reachable domain is calculated by traversing the boundary. Assuming that the reachable domain boundary is divided into N equal parts on one side, the three-dimensional off-orbit trajectory planning and reachable domain calculation process are designed as follows:

[0206] S210, calculate the closest point in the reachable domain: take the minimum range in the orbital plane as the optimization index, solve the trajectory with the shortest range, and record the minimum range Rmin ;

[0207] S220, calculate the farthest point in the reachable domain: take the maximum range in the orbital plane as the optimization index, solve the trajectory with the longest range, and record the maximum range R max ;

[0208] S230, calculate the left boundary: first simulate based on the finite thrust guidance algorithm to obtain the remaining fuel mass m of the aircraft r , m r Divide into N equal parts, and assume that the mass of the aircraft when the fuel is exhausted is m0. For the i-th (i=1,2...N-1) boundary point, the nearest boundary point is numbered i=1. It is approximately assumed that there are i·m r / N of fuel is used to change the flight direction and generate lateral motion. The aircraft is set to first perform a lateral braking angle component β = 90° with a constant positive sign. After this part of the fuel is consumed, longitudinal deorbit braking is performed. The longitudinal deorbit braking algorithm is the same as the algorithm in the previous section.

[0209] S240, calculate the right boundary: adopt the method of S230, set the sign of the lateral braking angle component β to be always negative, and obtain the right boundary trajectory.

[0210] Furthermore, in S300, specifically,

[0211] Considering the influence of the earth's rotation, the following reentry vehicle motion equation is established:

[0212]

[0213] During the reentry process of the spacecraft, in order to ensure the success of the flight mission, the path constraints such as heat flow, dynamic pressure, and overload must be met, which can be expressed as:

[0214]

[0215] in, q and n are heat flux, dynamic pressure and overload respectively, k Q is the heat flux coefficient, q max and n max are peak heat flux, peak dynamic pressure and overload peak, respectively.

[0216] So far, the optimal control problem is as follows:

[0217]

[0218] For nonlinear systems,

[0219]

[0220] There are bounds on state variables and control variables,

[0221]

[0222] Boundary constraints,

[0223]

[0224] Among them, ψol, ψ ou are the upper and lower bounds of ψ at the initial moment, ψ fl , ψ fu are the upper and lower bounds of the terminal time ψ,

[0225] Path constraints:

[0226] C l ≤C[x(τ),u(τ),τ]≤C u (26)

[0227] Among them, C l with C u are the upper and lower bounds of the path constraint,

[0228] Complex multi-constrained optimal control problems are typically solved through linearization. Here, we employ a sequential convex optimization method based on trust region constraints. By solving a series of convex optimization subproblems, we gradually approach the solution to the original problem.

[0229] After linearization, we get

[0230]

[0231] The equation of state is converted to x k The linearized dynamic equation is

[0232]

[0233] The basic constraints in Equation (26) can also be processed by the above convexification method. The iterative optimization process of sequential convex optimization is linearized around a given reference trajectory. In order to ensure that the subproblem approaches the original problem, the trust region constraint || xx needs to be added. k ||≤Δ,

[0234] In order to avoid the problem of “pseudo-infeasibility” in the early stage of iteration, reduce the dependence on the initial value, and compensate for the error caused by linearization, a virtual control variable is added to equation (28), and we can get:

[0235]

[0236] During the iteration process, the range of the virtual control amount in the objective function is controlled so that it is small in the later stage of the iteration, thereby ensuring the accuracy of the problem, that is,

[0237] J′=J+δ||Cw||2 (30)

[0238] In order to find the numerical solution to the infinite-dimensional optimal control problem, it must be discretized. For the sake of simplicity, the state variables are discretized into equidistant numerical points as usual, and the control quantity at each numerical point is expressed by the expression of the B-spline control point. The number of control points can be adjusted in the actual simulation, and the integral terms in the optimization process are all expressed by trapezoidal integrals. The simulation results show that this method can effectively suppress the jagged phenomenon of the control quantity in the numerical solution, is more conducive to iterative convergence, and reduces the number of variables that need to be solved. Dividing the time N into equal parts, the discretized dynamic equation can be obtained as follows:

[0239]

[0240] The superscript k represents the kth iteration, and the subscript i represents the t i time, Rearrange formula (31) to get:

[0241]

[0242] in

[0243]

[0244]

[0245] After converting the reentry trajectory optimization problem into a convex optimization problem, a single reentry glide trajectory can be solved offline.

[0246] Furthermore, in S400, the calculation of the reachable region requires obtaining the closed boundaries of its left boundary, right boundary, farthest point, and near boundary, so as to obtain the reachable area. The steps for solving the reachable region of the reentry segment are as follows:

[0247] S410, calculating the closest point of the axis range: using the closest vertical range of the initial shot as the performance index, the trajectory is calculated and recorded as L1;

[0248] S420, calculating the farthest point of the axis range: using the farthest range as a performance indicator, calculating the trajectory, which is recorded as L2;

[0249] S430, calculate the left and right boundaries: take the terminal ranges of L1 and L2 as the two endpoints, divide them into n equal parts, and evenly select n-1 points from small to large as the terminal range constraints. Optimize the performance indicators with the maximum and minimum horizontal ranges (including the sign, positive on the left and negative on the right), and obtain n-1 left and right boundaries. The left boundary trajectory is recorded as L3~L n+1 , the right boundary trajectory is L n+2 ~L 2×n ;

[0250] S440, axis calculation: take the terminal ranges of L1 and L2 as the two endpoints, divide them into n equal parts, evenly select n-1 points from small to large as the terminal range constraints, and set the terminal lateral range constraint to 0. Optimize with flight time as the performance indicator, and you can get n-1 trajectories on the axis, which are L 2×n+1 ~L 3×n-1 ;

[0251] S450, near-boundary calculation: take the initial upward and leftward range as the performance index and calculate the trajectory L 3×n , taking the initial upward range to the right as the performance index, the trajectory L is calculated 3×n+1 , respectively L 3×n and L1 terminal range, L 3×n+1 The terminal range of L1 is divided into two endpoints, and m points are evenly divided from small to large. m-1 points are evenly selected as the terminal range constraints. The maximum and minimum values ​​of the horizontal range under the fixed vertical range on one side of the initial upward direction are calculated respectively to obtain the trajectory L 3×n+2 ~L 3×n+(m-1)×4 ;

[0252] To quickly determine the reentry reachable region, an ellipse can be used to approximate the reachable region boundary of the reentry vehicle. This fast reachable region calculation method requires calculating four feature points to fit two ellipse boundaries, approximating them to the reachable region boundary of the reentry vehicle. The calculation process for the four feature points of the ellipse is as follows:

[0253] S451. Using the minimum vertical range of the initial shot as the performance indicator, calculate the trajectory, which is recorded as L1.

[0254] J=min(S r ) (35)

[0255] S452, taking the maximum vertical range of the initial shot as the performance index, calculate the trajectory, which is recorded as L2 (and the minimum range is recorded as S f , let the maximum longitudinal range be R max ),

[0256] J=min(-S r ) (36)

[0257] The formula for calculating the fixed range S by the maximum and minimum range points is as follows:

[0258] S=(S m +(S f -S m ) / 2)·c (37)

[0259] Where c is the empirical coefficient, and its value range is 0≤c≤1;

[0260] S453, fixed range S, calculate the trajectory with the maximum horizontal range as the performance index under range S, the trajectory is recorded as L3 (fixed vertical range is recorded as R, maximum horizontal range is recorded as Z max ),

[0261] J=min(-C r ) (38)

[0262] -0.001≤D r -R≤0.001 (39)

[0263] S454, fixed range S, calculate the trajectory with the minimum horizontal range as the performance indicator under range S, the trajectory is recorded as L4 (record the fixed vertical range as R, record the minimum horizontal range as Z min ),

[0264] J=min(C r ) (40)

[0265] -0.001≤D r -R≤0.001 (41)

[0266] According to the ellipse calculation formula as shown in formula (42),

[0267]

[0268] The center e, major axis a, and minor axis b1 and b2 of the ellipse are obtained respectively.

[0269]

[0270] Furthermore, in S500, specifically:

[0271] In order to measure the reachability of the aircraft during the re-entry glide phase under the influence of the no-fly zone, the reachable area under the influence of the no-fly zone was quantitatively calculated based on the analytical method.

[0272] Given the maximum flight capability of an aircraft, its range-energy relationship is determined. When a no-fly zone exists, lateral maneuvers are necessary to circumvent it, which inevitably increases the ground track length. If the ground track length caused by the no-fly zone exceeds the maximum flight capability, the mission will inevitably fail, and the target point will be unreachable. If the ground track length caused by the no-fly zone exceeds the difference between the maximum flight capability and the longitudinal capability loss caused by lateral movement, the mission will also fail, and the target point will be unreachable. If the ground track length caused by the no-fly zone is less than this difference, the mission will succeed, and the target point will be reached. Therefore, the key issue is how to quantitatively analyze the longitudinal capability loss caused by lateral movement. Through extensive simulation analysis or analytical calculations, the impact of a certain lateral variable on energy / range can be determined, and then calculated and compensated.

[0273] The Dubins curve is used to calculate the no-fly zone avoidance trajectory, assuming that the aircraft is in a balanced flight with no sideslip during the turning process and the speed angle is approximately zero:

[0274]

[0275] The component of lift produces centripetal force, assuming the turning radius is r z , we can get:

[0276]

[0277] Combining equations (44) and (45), ignoring the centrifugal force in the horizontal plane, we obtain an approximate analytical solution for the turning radius of the aircraft:

[0278]

[0279] Then maintain the turning radius r z The bank angle for balanced flight is:

[0280]

[0281] Assuming that the flight altitude before the turn is h1 and the atmospheric density is ρ1, the values ​​after the turn are h2 and ρ2, which are estimated by the atmospheric density fitting formula.

[0282]

[0283] In order to maintain balanced flight after turning, the atmospheric density needs to be satisfied.

[0284]

[0285] Get the relationship between the height change and the roll angle,

[0286]

[0287] The energy-based range prediction is obtained as,

[0288]

[0289] The range reduction caused by turning is calculated analytically using equation (51).

[0290] The effectiveness of the present invention is illustrated by taking the calculation of the reachable area of ​​a space-based re-entry vehicle as an example.

[0291] (1) Simulation of accessible area during derailment braking

[0292] According to the calculation process of the reachable domain of the deorbit braking segment, the maximum simulation time of the deorbit braking segment is set to 1500s, and the expected re-entry angle is θ f =-3.5°, the other simulation conditions remain unchanged, the simulation results are as follows Figure 3 shown.

[0293] Under these simulation conditions, the lateral maneuverability during the deorbit braking phase is approximately 400 km. The vehicle's fuel is completely depleted, and the resulting lateral maneuvering can alter the vehicle's original trajectory. Depending on the range of the lateral maneuver, the trajectory can vary by 1° to 3°, thereby broadening the lateral flight capabilities during the reentry glide phase. Typically, the ground projection of a vehicle's reachable region during the deorbit braking phase is approximately elliptical. In this example, due to the maximum deorbit time constraint during the deorbit braking phase, the elliptical shape is not present at the far boundary of the reachable region. If the deorbit time constraint is removed or set to a larger value, the elliptical shape will become more pronounced.

[0294] (2) Simulation of the reachable area in the re-entry glide phase without a no-fly zone

[0295] By selecting different performance indicators for the reentry trajectory optimization problem, the reachable region under specific flight conditions can be calculated. The head heat flux density constraint during the flight process is Total heating capacity is not more than 5000MJ / m 2 , dynamic pressure constraint Terminal pneumatic overload No more than 15g. Consider the control amplitude constraint as the maximum value of the attack angle and sideslip angle in the aerodynamic data table (using interpolation), that is, α∈[-16°,16°], β∈[-16°,16°], γ c ∈[-90°,90°]. In actual systems, the control quantity is subject not only to amplitude constraints but also to rate constraints. Since the control quantity rate constraint is not required yet, it is set as follows:

[0296] The reachable domain is the set of spatial positions of the end points of multiple trajectory curves. Connecting the terminal positions of each curve gives the reachable domain of the aircraft under this set of initial conditions.

[0297] In order to further illustrate the advantages of ellipse fitting reachable domain in terms of accuracy and speed, trajectory planning method, convex optimization method and ellipse approximation fast calculation method are used for simulation respectively. The simulation results are shown in the figure below. Figure 5 and as shown in Table 1.

[0298]

[0299] Table 1 Comparison of calculation time

[0300] It can be seen that the fast calculation method of the reachable domain based on ellipse approximation has higher calculation accuracy. Compared with the reachable domain calculation method based on convex optimization, its boundary is almost the same, indicating that the ellipse fitting method has higher accuracy. In addition, compared with the reachable domain calculation method based on trajectory planning, the boundary is significantly larger. Figure 5 It can also be seen that the fast calculation method of the reachable domain proposed in this section significantly shortens the calculation time. Although it is still longer than the trajectory planning-based method, it has the ability to perform online calculations.

[0301] (3) Simulation of the reachable area of ​​the no-fly zone constraint during the re-entry glide phase

[0302] The remaining range is simulated using analytical prediction. Assuming that the aircraft turns according to the Dubins curve after passing the closest point, the calculation results of the reachable area with / without the no-fly zone are as follows: Figure 6-Figure 8 shown.

[0303] Depend on Figure 7 It can be seen that as the turning distance increases, the range loss increases, and the two are positively correlated. The range loss effect is more obvious as the turning distance increases. This is because the turning causes the flight altitude to decrease, the resistance to increase, and the aircraft energy loss is too fast. When the no-fly zone is located on the upper side of the accessible area, it has no effect on the lower side. The no-fly zone only affects the area behind the no-fly zone. A target point is randomly set behind the no-fly zone, and the aircraft is guided to fly towards the target point. Figure 8 It can be seen that due to the influence of the no-fly zone, the boundary of the reachable region "collapses" behind the no-fly zone. The method based on analytical prediction of flight range can quickly determine the boundary of the reachable region under the influence of the no-fly zone and whether the target point is reachable.

[0304] It can be verified from the above examples that the method for quickly calculating the reachable area of ​​a space-based reentry vehicle considering the no-fly zone constraints proposed in the present invention is reasonable and effective.

Claims

1. A method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints, characterized by: The method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints comprises the following steps: S100. First, a derailment trajectory planning method based on convex optimization is proposed. The derailment braking problem is transformed into a convex optimization problem, and a second-order cone programming method is used to solve it to obtain the optimal trajectory of a single derailment segment. In S100, the following steps are specifically included: S110, normalization of state equation: definition is the thrust-acceleration vector, where is the ratio of thrust to initial weight of the engine, m is the ratio of current mass of the aircraft to initial mass, and is defined as ,Will , respectively and Normalize and convert time right Normalize and get the changed state equation, (1) in, is a slack variable, representing thrust-acceleration, and the state quantity is , the control quantity is , write the state equation in matrix form, (2) in, , , the upper bound of engine thrust acceleration and the definition of relaxation relationship constraints are, (3) (4) When using the iterative method to find the optimal solution, Equation (2) is regarded as a linear equation; S120, solve the initial trajectory: First, an initial trajectory is needed, and optimization is performed based on it. Here, the optimal fuel solution for reaching a fixed target point and specifying the time for the deorbit process is selected. (5) The deorbit process is divided into Each time interval is , solve The state variables and control variables, the initial position, velocity and terminal position and velocity constraints of the aircraft can be expressed as: (6) (7) The deorbit process is constrained by equations (2) to (7). is a normalized variable, which is approximately 1 during the deorbit process. When the accuracy of the result is not high, the Approximately 1, it becomes a linear equation, and the fuel index function is: (8) Since all constraints and objective functions are convex, this problem is a standard convex optimization problem; S200. Secondly, based on the optimal trajectory of a single off-track segment, the reachable domain is calculated by traversing the boundary. The boundary of the reachable domain is divided into N equal parts on one side. A three-dimensional off-track trajectory planning and reachable domain calculation process is designed to calculate the nearest point, farthest point, left boundary, and right boundary of the reachable domain. S300, then, a re-entry glide trajectory planning method based on convex optimization is proposed, which transforms the re-entry trajectory optimization problem into a convex optimization problem, and obtains the optimal trajectory of a single re-entry glide segment; S400, further changing the performance index, and obtaining the reentry segment reachable domain by traversing the boundary; S500. Finally, the reachable domains of the deorbit segment and the re-entry glide segment are combined, and the reachable capability of the re-entry glide segment affected by the no-fly zone is calculated by combining the analytical prediction method, and finally the reachable domain boundary of the space-based re-entry glide vehicle considering the no-fly zone is obtained.

2. The method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints according to claim 1 is characterized in that: In S100, the following steps are also included: S130, iteratively solving the time-optimal derailment braking problem: Introducing the time variable, both Equation (2) and Equation (7) are nonlinear equality constraints. At the same time, multiple iterations are required based on the initial solution. When the convex optimization problem contains nonlinear equality constraints When , assuming the current iteration result is ,in , solve The process is as follows: exist Taylor expansion of the function is performed at Written as (9) in, is a second-order correction that needs to be evaluated. At the same time, there is an expansion type. (10) in, For high-order corrections, only approximate solutions are needed , which is approximately , The constraints of the convex optimization problem are, (11) (12) The optimization result is ,Will exist First-order expansion Multiply both sides of the equation ,get, (13) The left side of formula (13) is taken as The new constraint is, (14) (15) Non-convex constraints Convert to a convex constraint, In the time-optimal derailment braking problem, in addition to considering equations (2), (3), and (4), the remaining constraints are: (16) (17) The indicator function is: (18) The discrete form of formula (2) and formula (16) are nonlinear constraints. According to the above, we can obtain the constraint conditions of the convex function and add the following constraints: (19) That is to obtain as well as , iterates multiple times until convergence, and obtains the time-optimal deorbit braking result.

3. The method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints according to claim 1, characterized in that: Specifically, in S200, the reachable domain is calculated by traversing the boundary. Assuming that the reachable domain boundary is divided into N equal parts on one side, a three-dimensional off-track trajectory planning and reachable domain calculation process is designed, including the following steps: S210, calculate the closest point in the reachable domain: take the minimum range in the orbital plane as the optimization index, solve the trajectory with the shortest range, and record the minimum range R min ; S220, calculate the farthest point in the reachable domain: take the maximum range in the orbital plane as the optimization index, solve the trajectory with the longest range, and record the maximum range R max ; S230, calculate the left boundary: first simulate based on the finite thrust guidance algorithm to obtain the remaining fuel mass of the aircraft ,Will Divide into Assume that the mass of the aircraft when the fuel is exhausted is , for the Boundary points, , the nearest boundary point number is , it is approximately assumed that the remaining fuel The fuel is used to change the flight direction and generate lateral motion, setting the aircraft to first perform lateral braking angle component , the sign is always positive. When this part of the fuel is consumed, the longitudinal derailment braking is performed. The longitudinal derailment braking algorithm is the same as the algorithm in the previous section. S240, calculate the right boundary: use the method of S230 to make the lateral braking angle component The sign is always negative, which means the right boundary trajectory is obtained.

4. The method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints according to claim 1 is characterized in that: In S300, specifically, Considering the influence of the earth's rotation, the following reentry vehicle motion equation is established: (20) During the re-entry process of the spacecraft, it is expressed as, (21) in, , and are heat flux, dynamic pressure and overload respectively, is the heat flux coefficient, , and are peak heat flux, peak dynamic pressure and overload peak, respectively. So far, the optimal control problem is as follows: (22) For nonlinear systems, (23) There are bounds on state variables and control variables, (24) Boundary constraints, (25) in, 、 The initial moments The upper and lower bounds of , Terminal time The upper and lower bounds of Path constraints: (26) in, and are the upper and lower bounds of the path constraint, Adopting a sequential convex optimization method based on trust region constraints, we gradually approach the solution of the original problem by solving a series of convex optimization sub-problems. After linearization, we get (27) The state equation (23) is The linearized dynamic equation is (28) The basic constraints in Equation (26) can also be handled by the above convex optimization method. The iterative optimization process of sequential convex optimization is linearized around a given reference trajectory. In order to ensure that the subproblem approaches the original problem, a trust region constraint needs to be added. , Adding a virtual control variable in equation (28), we get: (29) During the iteration process, the range of the virtual control amount in the objective function is controlled so that it is small in the later stage of the iteration, thereby ensuring the accuracy of the problem, that is, (30) To solve the numerical solution of infinite-dimensional optimal control problem, we discretize Eq. (29) and transform the time Divide it into equal parts, and the discretized dynamic equation is: (31) Among them, the superscript Indicates the Iterations, subscript express time, , , , rearrange formula (31) to get, (32) in (33) (34) The above method transforms the reentry trajectory optimization problem into a convex optimization problem.

5. The method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints according to claim 1 is characterized in that: In S400, the following steps are specifically included: S410, calculating the closest point of the axis range: using the closest vertical range of the initial shot as the performance indicator, calculating the trajectory, which is recorded as L1; S420, calculating the farthest point of the axis range: using the farthest range as a performance indicator, calculating the trajectory, which is recorded as L2; S430, calculate the left and right boundaries: take the terminal ranges of L1 and L2 as the two endpoints, divide them into n equal parts, and evenly select n-1 points from small to large as the terminal range constraints. Optimize with the maximum and minimum horizontal ranges as performance indicators, and obtain n-1 left and right boundaries. The left boundary trajectory is recorded as L3~L n+1 , the right boundary trajectory is L n+2 ~L 2×n ; S440, axis calculation: take the terminal ranges of L1 and L2 as the two endpoints, divide them into n equal parts, evenly select n-1 points from small to large as the terminal range constraints, and set the terminal lateral range constraint to 0. Optimize with flight time as the performance indicator, and you can get n-1 trajectories on the axis, which are L 2×n+1 ~ L 3×n-1 ; S450, near-boundary calculation: take the initial upward and leftward range as the performance index and calculate the trajectory L 3×n , taking the initial upward range to the right as the performance index, the trajectory L is calculated 3×n+1 , respectively L 3×n and L1 terminal range, L 3×n+1 The terminal range of L1 is divided into two endpoints, and m points are evenly divided from small to large. m-1 points are evenly selected as the terminal range constraints. The maximum and minimum values ​​of the horizontal range under the fixed vertical range on one side of the initial upward direction are calculated respectively to obtain the trajectory L 3×n+2 ~ L 3×n+(m-1)×4 ; An ellipse is used to approximate the boundaries of the reentry vehicle's reachable area. The calculation process of the four characteristic points of the ellipse includes: S451. Using the minimum vertical range of the initial shot as the performance indicator, calculate the trajectory, which is recorded as L1. (35) S452. Using the maximum vertical range of the initial shot as the performance indicator, calculate the trajectory, which is recorded as L2. (36) The formula for calculating the fixed range S by the maximum and minimum range points is as follows: (37) In the formula is the empirical coefficient, and its value range is ; S453, fixed range S, calculate the trajectory with the maximum horizontal range as the performance index under range S, this trajectory is recorded as L3, (38) (39) S454, fixed range S, calculate the trajectory with the minimum horizontal range as the performance indicator under range S, this trajectory is recorded as L4, (40) (41) According to the ellipse calculation formula as shown in formula (42), (42) Get the center of the ellipse , semimajor axis , semi-minor axis They are (43)。 6. The method for rapidly calculating the reachable area of ​​a space-based reentry vehicle considering no-fly zone constraints according to claim 1, characterized in that: In S500, specifically: In order to measure the reachability of the aircraft during the re-entry gliding phase under the influence of the no-fly zone, a quantitative calculation of the reachable area under the influence of the no-fly zone is performed based on the analytical method. The Dubins curve is used to calculate the no-fly zone avoidance trajectory, assuming that the aircraft is in a balanced flight with no sideslip during the turning process and the speed angle is approximately zero: (44) The components of lift produce centripetal force, assuming the turning radius is , we can get: (45) Combining equations (44) and (45), ignoring the centrifugal force in the horizontal plane, we obtain an approximate analytical solution for the turning radius of the aircraft: (46) Maintain turning radius The bank angle for balanced flight is: (47) Assumed flight altitude before turning and atmospheric density , the value after turning is 、 , estimated by the atmospheric density fitting formula, (48) In order to maintain balanced flight after turning, the atmospheric density needs to be satisfied. (49) Get the relationship between the height change and the roll angle, (50) The energy-based range prediction is obtained as, (51) The range reduction caused by turning is calculated analytically using Equation (51).

Images