A conflict-minimized trajectory collaborative planning method considering the time-varying nature of upper-level winds

Abstract

The present invention discloses a conflict - minimized trajectory collaborative planning method considering the time - variability of upper - air winds. By introducing a two - layer planning method, the trajectory planning problem considering the time - variability of upper - air winds in the strategic stage is studied. Considering the influence of flight time and conflicts on flight operation efficiency and safety, taking the flight departure time and altitude layer as decision variables, and taking the shortest flight time and conflict - free flight of flights in the flexible airway airspace as the objective function, a multi - flight strategic trajectory collaborative planning model is constructed. A two - layer planning algorithm adapted to the multi - flight strategic trajectory collaborative planning model is designed according to the characteristics of the strategic trajectory planning problem. The present invention provides a fast implementation method for the conflict - minimized trajectory planning method considering the time - variability of upper - air winds, and provides technical support for reasonably arranging flight plans, ensuring flight safety and reasonable allocation of airspace resources.

Description

Technical Field

[0001] The present invention belongs to the field of air traffic planning, and particularly relates to a conflict-minimizing trajectory collaborative planning method considering the time-varying nature of upper-level winds. Background Art

[0002] Currently, with the continuous development of the civil aviation transportation industry, the contradiction between traffic demand and traffic supply has become increasingly prominent, and the traditional air traffic management method has gradually shown its backwardness. At present, the operation based on four-dimensional trajectories provides a systematic solution to break through this bottleneck. Existing research on trajectory planning only considers the impact of upper-level wind uncertainty on trajectory planning. However, for trajectory planning in the strategic stage, due to the long lead time and large flight time span of flights, the changes of upper-level winds over time and space cannot be ignored. The optimal trajectory planned only considering uncertainty may not be able to improve the operation efficiency over time. In addition, existing research often separates horizontal route planning from conflict detection and resolution, that is, a two-stage trajectory planning model is established. However, during conflict resolution, with the adjustment of flight departure times and altitude levels, the two-stage planning model cannot adjust the optimal route according to specific situations. Summary of the Invention

[0003] Object of the Invention: Aiming at the above disadvantages, the present invention discloses a conflict-minimizing trajectory collaborative planning method considering the time-varying nature of upper-level winds, and establishes a two-layer planning model for trajectory planning, comprehensively considering flight operation efficiency and safety.

[0004] Technical Solution: The present invention provides a conflict-minimizing trajectory collaborative planning method considering the time-varying nature of upper-level winds, which specifically includes the following steps:

[0005] (1) Pre-acquire the ensemble forecast data of time-varying upper-level winds, flight historical flight plans, airspace information, and flight safety interval standards;

[0006] (2) Obtain the entry and exit points of flights in the flexible route airspace according to the flight historical flight plans;

[0007] (3) Perform grid processing on the flexible route airspace according to the airspace information, and establish connection rules between nodes;

[0008] (4) According to the ensemble forecast data of time-varying upper-level winds, calculate the expected arrival time of each flight at the route point, determine the upper-level wind conditions affecting the flight on the current flight segment, obtain the set of flight times of the flight on the flight segment, and fit the flight time distribution of the flight on the current flight segment; use the expected value of the flight segment flight time as the flight segment flight cost and plan the optimal route with the goal of minimizing the flight time, and use the variance as a constraint condition to ensure the predictability of the trajectory;

[0009] (5) Based on the planned optimal route, use the equal-interval sampling method to form a set of track sampling points; fit the time distribution passing through the points according to the set of flight arrival sampling point times; use the grid conflict detection method to conduct a preliminary screening of flight conflicts; further identify the conflicts between flights according to the flight safety interval standard and conflict probability to obtain flights with conflicts with each other.

[0010] (6) According to the dynamic grouping rules, put the flights with conflicts with each other into the same group, and allocate different population sizes to each group according to the number of conflicts and the number of flights within the group. Use the variable population size operator to dynamically design the sizes of various population sizes.

[0011] (7) For each group of flights, adjust the takeoff times and flight altitude layers of each flight within the group to avoid conflicts, and initialize each population; then use a fast genetic algorithm to optimize each population to obtain the optimal solution of each population; obtain the 4D track with the minimum multi-flight conflict through the co-evolution between each sub-population.

[0012] (8) Input the takeoff times and altitude layers of each flight corresponding to the 4D track with the minimum conflict into step (4), and consider the time-varying nature of the upper-air wind to plan the optimal route under the current takeoff time and altitude layer conditions to obtain the current optimal 4D track.

[0013] (9) Repeat steps (4) to (8) until the maximum number of iterations N is reached, and output the 4D track with the minimum multi-flight conflict.

[0014] Further, step (3) includes the following steps:

[0015] (31) Establish a rectangular search area that covers all entry and exit points of the flexible route airspace. The latitude range of the area is from λ min to λ max , and the longitude range is from φ min to φ max . Define a regular airway point grid inside, with a latitude interval of Δ lat and a longitude interval of Δ lon .

[0016] (32) Construct airway connection rules. With the current airway point as the center, fully connect all airway points within the 7×7 square area. The set of these airway points is defined as connectable airway points; when connecting to the connectable airway points, some airway points within the 7×7 range overlap. Define the overlapping airway points within the 7×7 range as non-connectable airway points; model the airway network as a graph , where represents the airway network nodes, and ε represents the airway segments.

[0017] Further, step (4) includes the following steps:

[0018] (41) Calculate the set of flight times for flight segments:

[0019]

[0020] Among them, is the set of EPS ensemble prediction members. For flight f ∈ F, at a given altitude layer h ∈ H, According to the calculation method of flight segment flight time, when the expected arrival time of f at flight segment e is t under the condition of member k, the flight time of f on flight segment e is the number of EPS members;

[0021] (42) Fit the probability distribution of flight segment flight time:

[0022] For a feasible route of flight f Let the arrival time of f at flight segment e i at the starting point be and the arrival time at the end point be The mathematical expectation μ and variance σ of the flight segment flight time 2 The unbiased estimators of are:

[0023]

[0024] Among them, and are calculated by the following formula:

[0025]

[0026]

[0027] Among them, δ f is the departure time of the flight, and ΔT is the time from takeoff to entering FRA;

[0028] (43) Calculate the mathematical expectation and variance of flight segment flight time:

[0029] There are no time-varying points on the flight segment: When and the wind field conditions do not change. At this time, the flight time of the flight segment Variance

[0030] There are time-varying points on the flight segment: When and the wind field conditions change. When the flight is affected by Wind1, and the probability density function of the flight segment flight time is f1(t); when When the flight is affected by Wind2 and the probability density function of the flight time of the flight segment is f2(t), then:

[0031]

[0032]

[0033]

[0034] (44) According to the expected value and variance of the flight time of the flight segment, the genetic algorithm is used to plan the optimal route, and combined with the take-off time and altitude layer, the optimal flight track with the shortest flight time is obtained.

[0035] Further, the implementation process of the step (5) is as follows:

[0036] (51) The track sampling is carried out by means of equal-interval sampling:

[0037] If there are time-varying points on the flight segment, before and after t vp , the flight is affected by different wind fields. The flight segment is segmented at t vp , regarded as two sub-flight segments, and t vp is used as the sampling point; if there are no time-varying points on the flight segment, the track sampling is directly carried out at equal intervals; calculate the set of passing times of the flight track sampling points, assume that it follows a normal distribution, and fit the probability distribution;

[0038] (52) The four-dimensional coordinates of the track sampling points are discretely mapped to the corresponding grid cell A m.n.l.k for the preliminary screening of flight conflicts of flights by using the grid conflict detection method;

[0039] Establish a 4D spatio-temporal grid corresponding to the flight conflict preliminary screening method: Define the neighborhood of a grid A m.n.l.k as itself and its 3 3 -1 = 26 grid cells in the spatial dimension, denoted as:

[0040] C m,n,l,k =[c m,n,l-1,k c m,n,l,k c m,n,l+1,k

[0041] Among them:

[0042]

[0043]

[0044]

[0045] In the formula, the matrix c m,n,l,k represents the nine-grid neighborhood of this layer; the matrix c​m,n,l-1,k Represents the nine-grid neighborhood of the upper layer of the current grid, c m,n,l+1,k Represents the nine-grid neighborhood of the lower layer of the current grid; a M,N,L,k Represents the grid with coordinates L, M, N, k, where M takes m, m + 1, m - 1; N takes n, n - 1, n + 1; L takes l - 1, l + 1, l;

[0046] Check the 27 grids in the neighborhood of grid A in sequence m.n.l.k If there are track sampling points of other flights in any grid in the neighborhood, it indicates the existence of potential flight conflicts, and then precisely detect the conflicts between flights; otherwise, it indicates the absence of potential flight conflicts;

[0047] Precisely detect the conflicts between flights: For the track sampling points s and s’ of two different flights in the same grid, calculate the vertical distance and horizontal distance between the flights within the set of expected passing times of the flights according to their coordinates. When both the vertical distance and the horizontal distance violate the flight safety interval, calculate the conflict probability; otherwise, it indicates the absence of flight conflicts;

[0048] Calculate the conflict probability: Assume that the expected flight time of flight f on flight segment e i is and the variance is Define the maximum arrival time error:

[0049]

[0050]

[0051] where, is the expected moment of arrival at the sampling point. Assume and Increase linearly according to the flight time;

[0052] When there is no t on the flight segment vp then:

[0053]

[0054]

[0055]

[0056] When there is t on the flight segment vp at t vp before and after, the flight is affected by different wind fields; at this time, the flight segment is segmented at t vp and regarded as two sub-flight segments, and t vp is used as the sampling point, and the parameters of the corresponding sampling point are then calculated according to when there is no t on the flight segment vpTime calculation;

[0057] Consider flight f1 in the airspace with sampling point s1; flight f2 with sampling point s2; then there is:

[0058]

[0059] Since they are independent of each other, their joint probability density function is:

[0060]

[0061] The conflict probability is:

[0062]

[0063]

[0064] If the conflict probability between the sampling points of two flight tracks is less than a predetermined value, it is determined as a flight conflict, and flights with conflicts with each other are obtained.

[0065] Furthermore, step (6) includes the following steps:

[0066] (61) Adopt a dynamic grouping strategy to put flights with conflicts with each other into the same group, that is, satisfy the following correlation:

[0067]

[0068] Among them, F represents the set composed of all flights; group k Represents the Kth group; And Are the ith and i'th flights in group k respectively, C ii' = 1 indicates that there is a conflict between flights i and i', and C ii' = 0 indicates that there is no conflict between flights i and i';

[0069] (62) Allocate different population sizes to each group according to the number of conflicts within the group and the number of flights, and adopt a variable population size operator to dynamically design the sizes of various populations; the variable population size operator is as follows:

[0070]

[0071] In the formula, s k is the population size of the kth group, pop_size is the initial population size, C k represents the number of conflicts within the kth group, L k represents the number of flights within the kth group, C represents the initial number of conflicts; L represents the initial number of flights.

[0072] Further, step (7) includes the following steps:

[0073] (71) Population initialization;

[0074] (72) For each initialized population group, use the genetic algorithm for optimization to obtain the optimal solution for each group; the optimal solution refers to the individual corresponding to the maximum value of the fitness function, and the expression of the fitness function is as follows:

[0075]

[0076]

[0077] where n represents the number of flights; represents the number of conflicts of the i-th flight in the k-th group; δ max represents the maximum δ i ; λ GH represents the ground delay cost coefficient, t i represents the flight departure time, t orig represents the initial flight departure time, T max represents the maximum allowable delay; λ FL represents the altitude layer adjustment cost coefficient, fl i represents the flight altitude layer, fl orig represents the initial flight altitude layer, FL max represents the maximum allowable altitude layer adjustment;

[0078] (73) Combine the optimal solutions among various populations to obtain a 4D trajectory with minimized multi-flight conflicts.

[0079] Advantageous effects: Compared with the prior art, the advantageous effects of the present invention are as follows: 1. The present invention not only considers the uncertainty of the upper-air wind, but also considers the time-varying nature of the upper-air wind, which is more in line with the actual situation of flight operations; it converts the uncertainty and time-varying nature of the upper-air wind into the uncertainty and time-varying nature of the flight passing time; 2. The present invention establishes a two-layer programming model for trajectory planning; after conflict resolution according to the current optimal trajectory, due to the adjustment of the departure time and flight altitude layer, the optimal flight route of the flight will change at this time; the two-layer programming model can optimize the flight route in a timely manner according to the current conditions, ensuring the maximization of flight operation efficiency and at the same time ensuring the safety of flight operations. Brief Description of the Drawings

[0080] Figure 1 is the flowchart of the present invention;

[0081] Figure 2 is the schematic diagram of the connection rule of the flight route network constructed by the present invention;

[0082] Figure 3It is a schematic diagram of the wind field change on both sides of the time-varying point in the present invention;

[0083] Figure 4 It is a schematic diagram of the conflict probability calculation method in the present invention;

[0084] Figure 5 It is a flowchart of the variable population size cooperative co-evolution algorithm based on dynamic grouping in the present invention;

[0085] Figure 6 It is a flowchart of the population initialization in the present invention. Detailed implementation manners

[0086] The present invention will be further described in detail below with reference to the accompanying drawings.

[0087] As Figure 1 shown, a conflict minimization track cooperative planning method considering the time-variation of high-altitude wind in the present invention specifically includes the following steps:

[0088] Step 1: Obtain the time-varying high-altitude wind ensemble forecast data, flight historical flight plans, airspace information, and flight safety interval standards.

[0089] Step 2: Obtain the longitude and latitude coordinates of the entry and exit points of the flight in the flexible route airspace, the time of entering the flexible route airspace, and the flight altitude layer according to the historical flight plan data (FPL) of the flight.

[0090] Step 3: Perform grid processing on the flexible route airspace according to the airspace information, and establish the connection rules between nodes.

[0091] Establish a rectangular search area covering all the entry and exit points of the flexible route airspace. The latitude range of the area is from λ min to λ max , and the longitude range is from φ min to φ max . Define a regular route point grid inside, with a latitude interval of Δ lat and a longitude interval of Δ lon .

[0092] Construct the route connection rules. As Figure 2 shown, where the circle represents the current route point; to ensure the flexibility of the route, the connectable route points are the route points within the 7×7 square area centered on the current route point, which are represented by diamonds in the figure; when connecting to the connectable route points, some route points overlap within the 7×7 range. Therefore, to reduce the complexity, the overlapping route points are defined as non-connectable route points, which are represented by squares in the figure; the route network can be modeled as a graph where represents the route network nodes and ε represents the flight segments.

[0093] Step 4: According to the time-varying upper-air wind ensemble forecast data, calculate the expected arrival time of each flight at the waypoint, determine the upper-air wind conditions affecting the flight on the current flight segment, obtain the set of flight times of the flight on the flight segment, and fit the flight time distribution of the flight on the current flight segment; take the expected value of the flight segment flight time as the flight segment flight cost and plan the optimal route with the goal of minimizing the flight time, and use the variance as a constraint condition to ensure the predictability of the flight track.

[0094] (4.1) Calculate the set of flight times for the flight segment:

[0095] Preprocess the ensemble forecast data of the wind to obtain the wind value at each waypoint. In Figure 2 , each waypoint contains the following information: latitude λ, longitude φ, altitude h, east-west wind component W E (t) and north-south wind component W N (t), and the two-norm of the wind is:

[0096]

[0097] Calculate the length of flight segment e:

[0098]

[0099] where

[0100]

[0101]

[0102] where and are the Cartesian coordinates of the start point w o and the end point w d of the flight segment respectively; assuming the Earth is a sphere and R is its average radius; · represents the dot product. For a given waypoint w, the Cartesian coordinates are calculated as follows:

[0103]

[0104] Calculate the course angle, wind direction angle and wind force values at the start and end points of flight segment e. The direction angle θ e is calculated as follows:

[0105]

[0106] Then the wind o at the start point w d and the end point w wo of flight segment e are:

[0107]

[0108]

[0109] Among them, θ wo (t) is the wind direction angle at the starting point, and θ wd (t) is the wind direction angle at the ending point.

[0110] Calculate the flight time of a flight under a single wind field condition:

[0111] The wind force value TW e (t) on the flight segment e adopts the average value of the wind force values at the starting point and the ending point, assuming that the wind force value on the flight segment is constant:

[0112]

[0113] Then the flight time of the flight on the flight segment e is:

[0114]

[0115] Among them, T a is the true airspeed of the flight, which is a constant value during the cruise phase.

[0116] Obtain the set of flight passing times: Let be the set of EPS ensemble forecast members. For a flight f ∈ F, at a given altitude layer h ∈ H, According to the above flight segment flight time calculation method, when the expected time for f to reach the flight segment e is t under the condition of member k, the flight time of f on the flight segment e Finally, obtain the set of flight times is the number of EPS members.

[0117] (4.2) Fit the probability distribution of the flight segment flight time:

[0118] By assuming that the set of flight times of the flight on the flight segment follows a normal distribution, fit the probability distribution of the flight times of the flight on each flight segment.

[0119] For a feasible route of flight f Let f reach the flight segment e i The time at the starting point is The time to reach the ending point is The mathematical expectation μ and variance σ 2 of the flight segment flight time are the unbiased estimators:

[0120]

[0121] Among them, and are calculated by the following formula:

[0122]

[0123]

[0124] δ f δ is the take-off time of the flight, and ΔT is the time from take-off to entering FRA, assuming it is a fixed value.

[0125] (4.3) Calculate the mathematical expectation and variance of the flight time of the flight segment:

[0126] There is no time-varying point on the flight segment: When and the wind field conditions do not change. At this time, the flight time of the flight segment Variance

[0127] There is a time-varying point on the flight segment: When and the wind field conditions change, as shown in Figure 3 When the flight is affected by Wind1, and the probability density function of the flight time of the flight segment is f1(t); when the flight is affected by Wind2, and the probability density function of the flight time of the flight segment is f2(t). Then:

[0128]

[0129]

[0130]

[0131] (4.4) According to the expected value and variance of the flight time of the flight segment, use the genetic algorithm to plan the optimal route, and combine the take-off time and altitude layer to obtain the optimal flight track with the shortest flight time.

[0132] Step 5: Based on the 4D flight track of the aircraft planned in Step 4, use the equal-interval sampling method to form a set of flight track sampling points; fit the time distribution passing through the points according to the set of times when the aircraft arrives at the sampling points; use the grid conflict detection method to conduct a preliminary screening of the flight conflicts of the aircraft; further identify the conflicts between the aircraft according to the aircraft safety interval standard and the conflict probability to obtain the flights with conflicts with each other.

[0133] (5.1) Flight track sampling: Adopt the equal-interval sampling method.

[0134] There is a time-varying point on the flight segment: Before and after t vp the flight is affected by different wind fields, as shown in Figure 3 At this time, divide the flight segment at t vp and regard it as two sub-flight segments, and regard t vpAs a sampling point. When there is no time-varying point on the flight segment: directly sample at equal intervals. Calculate the set of passing times of the flight track sampling points of the flight, assume that it follows a normal distribution, and fit the probability distribution.

[0135] (5.2) Grid conflict detection method: Discretize and map the four-dimensional coordinates of the flight track sampling points to the corresponding grid cell A m.n.l.k In, the grid conflict detection method is used for the initial screening of flight conflicts.

[0136] S1: Establish the 4D spatio-temporal grid corresponding to the following flight conflict initial screening method: Define a grid A m.n.l.k The neighborhood of is composed of itself and its 3 3 -1 = 26 grid cells in the spatial dimension, denoted as:

[0137] C m,n,l,k =[c m,n,l-1,k c m,n,l,k c m,n,l+1,k

[0138] Where:

[0139]

[0140]

[0141]

[0142] In the formula, the matrix c m,n,l,k Represents the nine-grid neighborhood of this layer; the matrix c m,n,l-1,k Represents the nine-grid neighborhood of the upper layer of the current grid, and c m,n,l+1,k Represents the nine-grid neighborhood of the lower layer of the current grid; a M,N,L,k Represents the grid with coordinates L, M, N, k, where M takes m, m + 1, m - 1; N takes n, n - 1, n + 1; L takes l - 1, l + 1, l.

[0143] S2: Check the 27 grids in the neighborhood of grid A m.n.l.k in turn. If there are flight track sampling points of other flights in any grid in the neighborhood, it indicates that there is a potential flight conflict, and step S3 is executed; otherwise, it indicates that there is no potential flight conflict.

[0144] S3: Perform precise detection of conflicts between flights: For the flight track sampling points s and s' of two different flights in the same grid, calculate the vertical distance and horizontal distance between the flights within the set of expected passing times of the flights according to their coordinates. When both the vertical distance and the horizontal distance violate the flight safety interval, execute S4; otherwise, it indicates that there is no flight conflict.

[0145] S4: Calculate the conflict probability: Let flight f be on flight segment e​i The expected value of the flight time above is The variance is Define the maximum error of the arrival time:

[0146]

[0147]

[0148] is the expected time to reach the sampling point. Let and Increase linearly according to the flight time.

[0149] ① When there is no t on the flight segment vp then:

[0150]

[0151]

[0152]

[0153] ② When there is t on the flight segment vp at this time, before and after t vp the flight is affected by different wind fields, as Figure 3 shown. At this time, the flight segment is segmented at t vp and regarded as two sub-flight segments, and t vp is used as the sampling point, and the parameters of the corresponding sampling point are calculated according to ①.

[0154] Consider flight f1 in the airspace with sampling point s1; flight f2 with sampling point s2. Then there is:

[0155]

[0156] Since they are independent of each other, their joint probability density function is:

[0157]

[0158] The conflict probability is as Figure 4 shown, that is:

[0159]

[0160]

[0161] If the conflict probability of the sampling points of two flight tracks is less than the predetermined value, it is determined as a flight conflict, and the flights with conflicts with each other are obtained.

[0162] Step 6: Put the flights with conflicts among each other into the same group according to the dynamic grouping rule, and assign different population sizes to each group according to the number of conflicts within the group and the number of aircraft. Adopt the variable population size operator to dynamically design the sizes of various populations, specifically as Figure 5 shown.

[0163] (6.1) Adopt the dynamic grouping strategy to put the flights with conflicts among each other into the same group, that is, satisfy the following relevance:

[0164]

[0165] where F represents the set composed of all flights; group k represents the Kth group; and are the ith and i'th flights in group k respectively, and C ii' =1 indicates that there is a conflict between flight i and i', and C ii' =0 indicates that there is no conflict between flight i and i'.

[0166] (6.2) Assign different population sizes to each group according to the number of conflicts within the group and the number of flights, and adopt the variable population size operator to dynamically design the sizes of various populations; the variable population size operator is as follows:

[0167]

[0168] In the formula, s k is the population size of the kth group, pop_size is the initial population size, C k represents the number of conflicts within the kth group, L k represents the number of flights within the kth group, C represents the initial number of conflicts; L represents the initial number of flights.

[0169] Step 7: For each group of flights, adjust the departure times and flight altitude levels of each flight within the group to avoid conflicts, and initialize each population; then adopt a fast genetic algorithm to optimize each population to obtain the optimal solution of each population; obtain the 4D trajectory with the minimum multi-aircraft conflict through the co-evolution among each sub-population.

[0170] Population initialization: As Figure 6 shown, during the population initialization process, the flight plans of each flight in each group constitute the initial feasible solution X initial of the group. For each flight in the initial feasible solution, according to its number of conflicts, the flight with the most conflicts is preferentially selected by means of roulette. When selecting the conflict resolution method, generate a random number r, 0 < r < 1, and let the probability of selecting the ground waiting strategy to resolve the conflict be P GH, when r < P GH , modify the flight departure time; when P GH ≤ r, modify the flight altitude level.

[0171] For each population group after initialization, use the genetic algorithm for optimization to obtain the optimal solution for each group; the optimal solution refers to the individual corresponding to the maximum value of the fitness function, and the expression of the fitness function is as follows:

[0172]

[0173]

[0174] Among them, n represents the number of flights; represents the number of conflicts of the i-th flight in the k-th group; σ max represents the maximum value of σ i ; λ GH represents the ground delay cost coefficient, t i represents the flight departure time, t orig represents the initial flight departure time, T max represents the maximum allowable delay; λ FL represents the altitude level adjustment cost coefficient, fl i represents the flight altitude level, fl orig represents the initial flight altitude level, FL max represents the maximum allowable altitude level adjustment;

[0175] The larger the function fitness value, the better the current solution. Each group iterates independently N times, terminates the optimization within the group, and selects the individual with the maximum fitness as the optimal solution within the group. Combine the optimal solutions among all populations to obtain a 4D trajectory with minimized multi-flight conflicts.

[0176] Step 8: Input the departure times and altitude levels of each aircraft corresponding to the conflict-minimized 4D trajectory into Step 4, and plan the optimal route considering the time-varying nature of the upper-level wind at the current departure time and altitude level to obtain the current optimal 4D trajectory;

[0177] Obtain the set of flight departure times and altitude levels after conflict resolution at the current evolutionary generation; use the set of departure times and altitude levels as inputs to participate in the upper-level planning model, i.e., Step 4. Since the upper-level wind is time-varying and the wind field conditions are different at different altitude levels, based on the current departure time and altitude level, plan the optimal route for each flight; combine the optimal route with the departure time and altitude level to generate the optimal trajectory and input it into the lower-level model, i.e., Step 5, as the initial condition for the next generation.

[0178] Step 9: Repeat Steps 4 to 8 until the maximum number of iterations N is reached; output the 4D trajectory with minimized multi-aircraft conflicts.

[0179] The above are only the preferred embodiments of the present invention, and do not impose any form of limitation on the present invention. Those skilled in the art can make various equivalent changes and improvements based on the above embodiments. Any equivalent changes and modifications made within the scope of the claims shall fall within the protection scope of the present invention.

Claims

1. A conflict - minimizing trajectory collaborative planning method considering the time - variability of upper - air winds, characterized in that, It includes the following steps: (1) Pre-acquire time-varying upper-air wind ensemble forecast data, flight historical flight plans, airspace information, and flight safety interval standards; (2) Obtain the entry and exit points of flights in the flexible route airspace according to the flight historical flight plans; (3) Perform grid processing on the flexible route airspace according to the airspace information and establish connection rules between nodes; (4) According to the time-varying upper-air wind ensemble forecast data, calculate the expected arrival time of each flight at the route point, determine the upper-air wind conditions affecting the flight of the flight on the current flight segment, obtain the flight time set of the flight on the flight segment, and fit the flight time distribution of the flight on the current flight segment; Take the expected value of the flight segment flight time as the flight segment flight cost and plan the optimal route with the goal of minimizing the flight time, and use the variance as a constraint condition to ensure the predictability of the flight track; (5) Based on the planned optimal route, use the equal-interval sampling method to form a flight track sampling point set; Fit the passing time distribution according to the flight arrival sampling point time set; Use the grid conflict detection method to conduct a preliminary screening of flight conflicts; Further identify the conflicts between flights according to the flight safety interval standard and the conflict probability to obtain flights with conflicts with each other; (6) Put the flights with conflicts with each other into the same group according to the dynamic grouping rule, allocate different population sizes to each group according to the number of conflicts and the number of flights in the group, and use the variable population size operator to dynamically design the size of each population; (7) For each group of flights, adjust the takeoff time and flight altitude layer of each flight in the group to avoid conflicts and initialize each population; Then use a fast genetic algorithm to optimize each population to obtain the optimal solution of each population; Obtain the 4D flight track with the minimum multi-flight conflict through the co-evolution between each sub-population; (8) Input the takeoff time and altitude layer of each flight corresponding to the 4D flight track with the minimum conflict into step (4), and plan the optimal route considering the time-variability of the upper-air wind under the current takeoff time and altitude layer conditions to obtain the current optimal 4D flight track; (9) Repeat steps (4) to (8) until the maximum number of iterations N is reached, and output the 4D flight track with the minimum multi-flight conflict; The said step (4) includes the following steps: (41) Calculate the flight time set of the flight segment; where \(K\) is the set of EPS ensemble forecast members, flight \(f\in F\), at a given altitude level \(h\in H\), According to the flight time calculation method for flight segments, when the expected arrival time of \(f\) at flight segment \(e\) is \(t\) under the condition of member \(k\), the flight time of \(f\) on flight segment \(e\) \(\vert K\vert\) is the number of EPS members; (42) Fit the probability distribution of the flight segment flight time; For a feasible route of flight f Let f reach flight segment e i The time at the starting point is The time to reach the end point is The unbiased estimators of the mathematical expectation μ and variance σ of the flight time of the flight segment are 2 as follows: wherein, and are calculated by the following formula: where δ f is the departure time of the flight, and ΔT is the time from the flight's departure to its entry into FRA; (43) Calculate the mathematical expectation and variance of the flight segment flight time; There is no time-varying point on the flight segment: When and the wind field conditions do not change, and at this time, the variance of the flight time of the flight segment variance There are time-varying points on the flight segment: When and occur, the wind field conditions change. When occurs, the flight is affected by Wind1, and the probability density function of the flight time of the flight segment is f1(t); when occurs, the flight is affected by Wind2, and the probability density function of the flight time of the flight segment is f2(t), then: (44) According to the expected value and variance of the flight segment flight time, use the genetic algorithm to plan the optimal route, and combine the takeoff time and altitude layer to obtain the optimal flight track with the shortest flight time; The implementation process of the said step (5) is as follows: (51) Conduct flight track sampling in the way of equal-interval sampling; If there is a time-varying point on the flight segment, before and after t vp , the flight is affected by different wind fields. The flight segment is segmented at t vp , regarded as two sub-flight segments, and t vp is used as the sampling point; if there is no time-varying point on the flight segment, the track is sampled at equal intervals directly; calculate the set of passing times of the flight track sampling points, assume that it follows a normal distribution, and fit the probability distribution; (52) Discretely map the four-dimensional coordinates of the track sampling points to the corresponding grid cell A m.n.l.k Among them, the grid conflict detection method is used for the initial screening of flight conflicts Establish a 4D spatio-temporal grid corresponding to the initial screening method for flight conflicts: Define a grid A m.n.l.k The neighborhood of is itself and its 3 3 -1 = 26 grid cells in the spatial dimension, denoted as: C m,n,l,k = [c m,n,l-1,k c m,n,l,k c m,n,l+1,k ​ Wherein: In the formula, matrix c m,n,l,k represents the nine-grid neighborhood of this layer; matrix c m,n,l-1,k represents the nine-grid neighborhood of the upper layer of the current grid, and c m,n,l+1,k represents the nine-grid neighborhood of the lower layer of the current grid; a M,N,L,k represents the grid with coordinates L, M, N, k, where M takes m, m + 1, m - 1; N takes n, n - 1, n + 1; L takes l - 1, l + 1, l; Check grid A in sequence m.n.l.k For the 27 grids in the neighborhood, if there are track sampling points of other flights in any grid in the neighborhood, it indicates that there is a potential flight conflict, and then the conflict between flights is detected precisely; otherwise, it indicates that there is no potential flight conflict; Precisely detect the conflicts between flights: For the flight track sampling points s and s' of two different flights in the same grid, calculate the vertical distance and horizontal distance between the flights within the flight expected passing point time set, and calculate the conflict probability when both the vertical distance and the horizontal distance violate the flight safety interval; Otherwise, it indicates that there is no flight conflict; Calculate the conflict probability: Assume that the expected flight time of flight f on flight segment e i is and the variance is Define the maximum arrival time error: Among them, To reach the expected time of the sampling point, let and Increase linearly according to the flight time; When t does not exist on the flight segment vp then: When there is a t on the flight segment vp At this time, before and after t vp , the flight is affected by different wind fields; at this time, the flight segment is segmented at t vp , regarded as two sub-flight segments, and t vp is used as the sampling point, and the parameters of the corresponding sampling point are calculated again according to the situation when there is no t on the flight segment vp ; Consider flight f1 in the airspace, with sampling point s1; Flight f2, with sampling point s2; Then there is: Since they are independent of each other, their joint probability density function is as follows: The conflict probability is: If the conflict probability between two track sampling points is less than a predetermined value, it is determined as a flight conflict, and flights with conflicts with each other are obtained; Step (6) includes the following steps: (61) Adopt a dynamic grouping strategy to put flights with conflicts with each other into the same group, that is, satisfy the following correlation: Among them, F represents the set composed of all flights; group k represents the Kth group; and are the ith and i'th flights in group k respectively. C ii′ = 1 indicates that there is a conflict between flights i and i', and C ii′ = 0 indicates that there is no conflict between flights i and i'; (62) Allocate different population sizes to each group according to the number of conflicts and the number of flights within the group, and adopt a variable population size operator to dynamically design the sizes of various populations; the variable population size operator is as follows: where s k is the population size of the k-th group, pop_size is the initial population size, C k represents the number of conflicts within the k-th group, L k represents the number of flights within the k-th group, C represents the initial number of conflicts; L represents the initial number of flights; Step (7) includes the following steps: (71) Population initialization; (72) Use the genetic algorithm to optimize each initialized population group to obtain the optimal solution of each group; the optimal solution refers to the individual corresponding to the maximum value of the fitness function, and the expression of the fitness function is as follows: Among them, n represents the number of flights; represents the number of conflicts of the i-th flight in the k-th group; δ max represents the maximum δ i ; λ GH represents the ground delay cost coefficient, t i represents the departure time of the flight, t orig represents the initial departure time of the flight, T max represents the maximum allowable delay; λ FL represents the altitude layer adjustment cost coefficient, fl i represents the flight altitude layer of the flight, fl orig represents the initial flight altitude layer of the flight, FL max represents the maximum allowable altitude layer adjustment amount; (73) Combine the optimal solutions between each population to obtain a 4D track with minimized multi-flight conflicts.

2. The conflict-minimizing trajectory collaborative planning method considering the time-varying nature of high-altitude winds according to claim 1, characterized in that, Step (3) includes the following steps: (31) Establish a rectangular search area covering all entry and exit points of the flexible airway airspace. The latitude range of the area is from λ min to λ max , and the longitude range is from φ min to φ max . Define a regular airway point grid inside, with a latitude interval of Δ lat and a longitude interval of Δ lon ; (32) Construct a route connection rule. With the current route point as the center, fully connect all route points within the 7×7 square area. This set of route points is defined as connectable route points; when connecting to connectable route points, some route points within the 7×7 range overlap. The overlapping route points within the 7×7 range are defined as non-connectable route points; the route network is modeled as a graph G(W,E), where W represents the route network nodes and E represents the flight segments.

Images