const KNOT = 1 / 60
mutable struct Aircraft
start::Tuple{Float64,Float64}
speed::Float64
heading::Float64
min_scale::Float64
max_scale::Float64
min_turn::Float64
max_turn::Float64
Aircraft() = new((0.0, 0.0), 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)
end
start_x(aircraft::Aircraft) = aircraft.start[1]
start_y(aircraft::Aircraft) = aircraft.start[2]
struct Control
v::Float64
theta::Float64
q::Float64
omega::Float64
end
velocity_x(control::Control) = control.v * cos(control.theta)
velocity_y(control::Control) = control.v * sin(control.theta)
struct Resolution
deviation::Float64
controls::Vector{Control}
end
function Resolution(deviation, v, theta, q, omega)
return Resolution(
deviation,
[
Control(v_i, theta_i, q_i, omega_i) for
(v_i, theta_i, q_i, omega_i) in zip(v, theta, q, omega)
],
)
end
mutable struct Airspace
num_aircraft::Int
fleet::Vector{Aircraft}
name::String
min_separation::Float64
ray_length::Float64
xlim::Tuple{Float64,Float64}
ylim::Tuple{Float64,Float64}
end
function Airspace(num_aircraft::Int)
return Airspace(
num_aircraft,
[Aircraft() for _ in 1:num_aircraft],
"",
0.0,
0.0,
(0.0, 1.0),
(0.0, 1.0),
)
endAircraft Deconfliction
Adjust aircraft headings and speeds by the smallest amount that keeps every pair of trajectories conflict-free.
Introduction
We consider the problem of avoiding collisions between aircraft in the sky.
Given a set of aircraft with their current speed and straight direction at the current instant, the controller needs to give each aircraft a new heading and speed that ensure no collision will occur. This problem is known as conflict resolution or tactical deconfliction, as opposed to strategic deconfliction, which happens hours before flight, and collision avoidance, which addresses imminent conflicts.
The trajectories live in the Euclidean 2D plane, since we assume every aircraft flies at the same flight level. We also assume there is no conflict at the initial instant.
Here are the trajectories of eight aircraft, first flown on their nominal headings until they collide, then as the solver reroutes them:
Problem description
Input
- A set of aircraft \(A\); for each aircraft \(i\in A\):
- \(\overrightarrow{P_i^0} = (X_i, Y_i)\) the initial position of aircraft \(i\)
- \(\overrightarrow{V_i} = (V_i\cos \Theta _i, V_i\sin \Theta _i)\) the nominal vector of velocity of aircraft \(i\), where \(V_i\) is the magnitude of nominal speed of aircraft \(i\), in NM/min, and \(\Theta _i\) is the nominal heading angle of aircraft \(i\), in radians
- \(0\le Q_i^{\text{min}}\le Q_i^{\text{max}}\) the minimum and maximum scale factors applied to nominal speed of aircraft \(i\)
- \(\Omega _i^{\text{min}}\le \Omega _i^{\text{max}}\), \(\left[\Omega _i^{\text{min}},~\Omega _i^{\text{max}}\right] \subset \left[-\dfrac{\pi}2,~ \dfrac{\pi}2\right]\), the minimum and maximum variation of the heading angle of aircraft \(i\), in radians
- A safety distance \(D\) between aircraft, in NM (nautical miles)
Time is counted in minutes throughout, so a nominal speed of 200 knots is 200 * KNOT with KNOT one nautical mile per hour expressed in NM/min.
Problem:
Choose the new angles and speeds of the aircraft in given ranges such that the resulting trajectories are conflict-free.
Objective:
Minimize the sum of the squares of the differences in speeds and angles of the aircraft.
Modeling the separation constraint
We denote here for each aircraft \(i\in A\):
- \(v_i\) the magnitude of modified speed of aircraft \(i\)
- \(\theta _i\) the modified heading angle of aircraft \(i\)
- \(q_i\) the scale factor applied to nominal speed of aircraft \(i\)
- \(\omega _i\) the variation of the heading angle of aircraft \(i\)
Given the set of aircraft \(A\) with initial positions \(\overrightarrow{P_i^0}\) and nominal vectors of velocity \(\overrightarrow{V_i}\), the aircraft deconfliction problem consists in finding new vectors \(\overrightarrow{v_i}\) such that the resulting trajectories are conflict-free, that is, such that \[ \forall i,j \in A,~i<j,~\forall t\ge 0, \qquad \lVert \overrightarrow{p_i}(t) - \overrightarrow{p_j}(t) \rVert \ge D \qquad (1) \]
where \(\overrightarrow{v_i} = (v_i\cos \theta _i, v_i\sin \theta _i)\) is the modified vector of velocity of aircraft \(i\) and \(\overrightarrow{p_i}(t) := \overrightarrow{P_i^0} + t\overrightarrow{v_i}\) is the new position of aircraft \(i\) at instant \(t\ge 0\).
Such a condition cannot be directly plugged into a mathematical programming solver, mainly because of the fact that one of its quantifiers is \(t\), which ranges on a continuous domain. There are several ways of turning \((1)\) into mathematical programming constraints (see “Aircraft deconfliction via Mathematical Programming: Review and insights”, Pelegrín and D’Ambrosio, 2020).
We will use here the “minimum-distance time equations” method which consists of finding, for each pair of aircraft, the critical time instant at which this constraint has to be satisfied. Since distances are positive, we base our analysis on the squared version of \((1)\).
Given a pair of aircraft \(i,j\in A\), the squared distance between \(i\) and \(j\) at time \(t \ge 0\) is
\[ f_{ij} (t)=\lVert \overrightarrow{p_i}(t) - \overrightarrow{p_j}(t) \rVert^2 =\lVert \overrightarrow{v_{ij}} \rVert^2 t^2 + 2 t \langle \overrightarrow{P^0_{ij}} , \overrightarrow{v_{ij}} \rangle + \lVert \overrightarrow{P^0_{ij}} \rVert^2 \qquad \qquad (2) \]
where \(~~\overrightarrow{P^0_{ij}} =\overrightarrow{P^0_i} - \overrightarrow{P^0_j}~~\) and \(~~\overrightarrow{v_{ij}} = \overrightarrow{v_i} - \overrightarrow{v_j}\)
By differentiating \(f_{ij}\), we obtain that the minimum separation between aircraft \(i\) and \(j\) is attained at time \(\qquad t^{\text{min}}_{ij} = \dfrac{- \langle \overrightarrow{P^0_{ij}} , \overrightarrow{v_{ij}} \rangle}{\lVert \overrightarrow{v_{ij}} \rVert^2}\)
If \(t^{\text{min}}_{ij} <0\), the minimum separation between \(i\) and \(j\) was attained in the past, and their trajectories are diverging in the given time horizon. Otherwise, by substituting \(t^{\text{min}}_{ij}\) in \((2)\), we obtain the minimum squared distance between \(i\) and \(j\) during their observed trajectories: \[ f^{\text{min}}_{ij} = \dfrac{- \langle \overrightarrow{P^0_{ij}} , \overrightarrow{v_{ij}} \rangle^2}{\lVert \overrightarrow{v_{ij}} \rVert^2} + \lVert \overrightarrow{P^0_{ij}} \rVert^2 \]
A new separation condition, also equivalent to \((1)\), can be thus stated as follows: \[ \forall i, j \in A,~i<j, \qquad t^{\text{min}}_{ij}<0 \quad \text{ or } \quad f^{\text{min}}_{ij} \ge D^2 \]
\[ \forall i, j \in A,~i<j, \qquad t^{\text{min}}_{ij}<0 \quad \text{ or } \quad \lVert \overrightarrow{v_{ij}} \rVert^2 (\lVert \overrightarrow{P^0_{ij}} \rVert^2 -D^2) - \langle \overrightarrow{P^0_{ij}} , \overrightarrow{v_{ij}} \rangle^2 \ge 0 \qquad (3) \]
We introduce a new variable \(k_{ij}=\left\{ \begin{array}{ll} 1 \text{ if } t^{\text{min}}_{ij} \ge 0\\[.2cm] 0 \text{ otherwise} \end{array} \right.\) so we can rewrite \((3)\) as:
\[ \forall i, j \in A,~i<j, \qquad k_{ij}\left(\lVert \overrightarrow{v_{ij}} \rVert^2 (\lVert \overrightarrow{P^0_{ij}} \rVert^2 -D^2) - \langle \overrightarrow{P^0_{ij}} , \overrightarrow{v_{ij}} \rangle^2 \right) \ge 0 \]
with the constraint checking the sign of \(t^{\text{min}}_{ij}\):
\[ \forall i, j \in A,~i<j, \qquad t^{\text{min}}_{ij} (2k_{ij} -1) \ge 0 \]
If we denote \(\overrightarrow{v_{ij}}=(v_{x,ij},v_{y,ij})\) with \(v_{x,ij}\) and \(v_{y,ij}\) two new variables for each \(i,j \in A\), the previous constraints become:
- Sign of \(t_{ij}^{\text{min}}\):
\[ \forall i, j \in A,~i<j, \qquad ((X_i - X_j) \cdot v_{x,ij} + (Y_i - Y_j) \cdot v_{y,ij}) \cdot (2 k_{ij} - 1) \le 0 \]
- Pairwise aircraft separation:
\[ \forall i, j \in A,~i<j, \qquad k_{ij}\left((v_{x,ij} ^ 2 + v_{y,ij} ^ 2) \cdot ((X_i - X_j)^2 + (Y_i - Y_j)^2 -D^2) - ((X_i - X_j) \cdot v_{x,ij} + (Y_i - Y_j) \cdot v_{y,ij})^2 \right) \ge 0 \]
Nonlinear model
Now that we have modeled the separation constraint in the nonlinear programming formalism, we can write the complete model.
Variables
- For each aircraft \(i\in A\):
- \(v_i \in \mathbb{R}\) the magnitude of modified speed of aircraft \(i\), in NM/min
- \(\theta _i \in \mathbb{R}\) the modified heading angle of aircraft \(i\), in radians
- \(q_i \in [Q_i^{\text{min}},~ Q_i^{\text{max}}]\) the scale factor applied to nominal speed of aircraft \(i\)
- \(\omega _i \in [\Omega _i^{\text{min}},~ \Omega _i^{\text{max}}]\) the variation of the heading angle of aircraft \(i\)
- For each pair of aircraft \(i, j\in A\):
- \(k_{ij} \in \{0, 1\}\) a binary variable which is \(1\) if minimum separation between aircraft \(i\) and \(j\) is attained in the future, else \(0\)
- \(v_{x,ij},~v_{y,ij} \in \mathbb{R}\) the coordinates of the difference between the modified vectors of velocity of aircraft \(i\) and \(j\) (see above)
Objective: minimize the sum of the squared deviations of the aircraft speeds and heading angles
\[ \min \sum_{i \in A} (q_i-1)^2 + \omega _i ^ 2 \]
Constraints
- Definition of the speed scale factor \(q_i\)
\[ \forall i \in A, \qquad v_i=q_i\cdot V_i \]
- Definition of the variation of the heading angle \(\omega _i\)
\[ \forall i \in A, \qquad \theta _i=\Theta _i + \omega _i \]
- Definition of \(v_x\) and \(v_y\)
\[ \forall i,j \in A,~i<j, \qquad v_{x,ij} = v_i \cdot \cos(\theta _i) - v_j \cdot \cos(\theta _j) \] \[ \forall i,j \in A,~i<j, \qquad v_{y,ij} = v_i \cdot \sin(\theta _i) - v_j \cdot \sin(\theta _j) \]
- Sign of \(t_{ij}^{\text{min}}\)
\[ \forall i, j \in A,~i<j, \qquad ((X_i - X_j) \cdot v_{x,ij} + (Y_i - Y_j) \cdot v_{y,ij}) \cdot (2 k_{ij} - 1) \le 0 \]
- Pairwise aircraft separation
\[ \forall i, j \in A,~i<j, \qquad k_{ij}\left((v_{x,ij} ^ 2 + v_{y,ij} ^ 2) \cdot ((X_i - X_j)^2 + (Y_i - Y_j)^2 -D^2) - ((X_i - X_j) \cdot v_{x,ij} + (Y_i - Y_j) \cdot v_{y,ij})^2 \right) \ge 0 \]
Input data
An Aircraft is where it starts, how fast it nominally flies and on what heading, plus how far its speed and heading may be pushed. A Control holds what the solver decides for one aircraft, a speed and a heading together with the two deviations that reach them. An Airspace is a fleet and the min_separation every pair has to keep, and a Resolution is one Control per aircraft together with the deviation the fleet paid.
Speeds are in nautical miles per minute, so KNOT converts from knots, and distances are in nautical miles throughout.
We will consider 4 use cases for the problem. The circle and the sector are stress tests rather than pictures of real traffic, with every pair in conflict and every conflict happening at once. The two converging lines and the ten random aircraft are closer to what a controller sees.
using Random
const SKEW = 1 / 100
function fleet_bounds!(airspace::Airspace)
for aircraft in airspace.fleet
aircraft.speed = 200 * KNOT
aircraft.min_scale, aircraft.max_scale = 0.9, 1.1
aircraft.min_turn, aircraft.max_turn = -pi / 16, pi / 16
end
return airspace
end
function load_circle(; num_aircraft=8, ray_length=100, min_separation=5, extent=nothing)
airspace = Airspace(num_aircraft)
airspace.name = "Circle air traffic"
airspace.min_separation = min_separation
airspace.ray_length = ray_length
radius = ray_length / 2
reach = extent === nothing ? radius + 2 * min_separation : extent
airspace.xlim = airspace.ylim = (-reach, reach)
fleet_bounds!(airspace)
for (i, aircraft) in enumerate(airspace.fleet)
bearing = 2 * pi / num_aircraft * (i - 1)
aircraft.start = (radius * cos(bearing + SKEW), radius * sin(bearing + SKEW))
aircraft.heading = pi + bearing
end
return airspace
end
function load_sector(;
num_aircraft=8,
ray_length=120,
min_separation=5,
sector_angle=pi / 2,
extent=nothing,
)
airspace = Airspace(num_aircraft)
airspace.name = "Sector circle air traffic"
airspace.min_separation = min_separation
airspace.ray_length = ray_length
radius = ray_length / 2
reach = extent === nothing ? radius + 2 * min_separation : extent
airspace.xlim = airspace.ylim = (-reach, reach)
fleet_bounds!(airspace)
for (i, aircraft) in enumerate(airspace.fleet)
bearing = sector_angle / (num_aircraft + 1) * i
aircraft.start = (radius * cos(bearing), radius * sin(bearing))
aircraft.heading = pi + bearing
end
return airspace
end
function load_lines(;
line_size=4, ray_length=150, min_separation=20, spacing=25, standoff=65
)
airspace = Airspace(2 * line_size)
airspace.name = "Two converging lines of $line_size aircraft"
airspace.min_separation = min_separation
airspace.ray_length = ray_length
tail = standoff + spacing * (line_size - 1)
airspace.xlim = (-tail - 3 * min_separation, ray_length)
airspace.ylim = (-ray_length, tail + 3 * min_separation)
fleet_bounds!(airspace)
for (i, aircraft) in enumerate(airspace.fleet)
if i <= line_size
aircraft.start = (0.0, standoff + spacing * (i - 1))
aircraft.heading = 3 * pi / 2
else
aircraft.start = (-standoff - spacing * (i - line_size - 1), 0.0)
aircraft.heading = 0.0
end
end
return airspace
end
function load_random(;
num_aircraft=10, ray_length=70, min_separation=5, spread=50, seed=9
)
airspace = Airspace(num_aircraft)
airspace.name = "$num_aircraft random aircraft"
airspace.min_separation = min_separation
airspace.ray_length = ray_length
reach = spread + 2 * min_separation
airspace.xlim = airspace.ylim = (-reach, reach)
fleet_bounds!(airspace)
rng = MersenneTwister(seed)
for aircraft in airspace.fleet
x, y = float(rand(rng, (-spread):spread)), float(rand(rng, (-spread):spread))
x == 0 && (x = 1.0)
aircraft.start = (x, y)
heading = atan(y / x) + rand(rng, -180:180) / 180 * pi / 4
aircraft.heading = x > 0 ? heading + pi : heading
end
return airspace
endEach aircraft carries a bubble of radius \(D/2\), so two of them are in conflict exactly when their bubbles overlap. The faint dashed ray behind a track is the heading it would have flown had nothing been asked of it, and the figures play the motion out, so the closest approach can be watched rather than inferred.
The circle is eight aircraft spread evenly round its rim, all converging on the center at the same speed. Left alone they arrive together and every bubble overlaps every other.
circle = load_circle()
animate_input(circle; show_circle=true)Model implementation
using JuMP
using KNITROminimize_deviation(airspace) builds the model, solves it with Knitro, and returns a Resolution. The objective keyword picks what is paid for. The default charges for both the speed and the heading deviation, and the two later sections charge for one of them alone.
The separation between a pair is nonconvex, and the binary k decides whether the closest approach lies in the future or is already behind them. That makes the model a mixed-integer nonlinear program. Knitro’s MIP multistart is switched on and a node limit keeps the solve time bounded, so what comes back is the best point found rather than a proven optimum.
function minimize_deviation(
airspace::Airspace, objective="speed+angle"; multistart=true, max_nodes=500
)
model = Model(KNITRO.Optimizer)
set_attribute(model, "mip_multistart", Int(multistart))
set_attribute(model, "mip_maxnodes", max_nodes)
n = airspace.num_aircraft
fleet = airspace.fleet
pairs = [(i, j) for i in 1:n for j in 1:n if i < j]
@variable(model, v[1:n])
@variable(
model,
fleet[i].heading + fleet[i].min_turn <=
theta[i in 1:n] <=
fleet[i].heading + fleet[i].max_turn
)
@variable(model, fleet[i].min_scale <= q[i in 1:n] <= fleet[i].max_scale)
@variable(model, fleet[i].min_turn <= omega[i in 1:n] <= fleet[i].max_turn)
@variable(model, k[1:length(pairs)], Bin)
@variable(model, v_x[1:length(pairs)])
@variable(model, v_y[1:length(pairs)])
if objective == "speed"
@objective(model, Min, sum((q[i] - 1)^2 for i in 1:n))
elseif objective == "angle"
@objective(model, Min, sum(omega[i]^2 for i in 1:n))
else
@objective(model, Min, sum((q[i] - 1)^2 + omega[i]^2 for i in 1:n))
end
for i in 1:n
aircraft = fleet[i]
@constraint(model, v[i] == aircraft.speed * q[i])
@constraint(model, theta[i] == aircraft.heading + omega[i])
end
for (ij, (i, j)) in enumerate(pairs)
fi, fj = fleet[i], fleet[j]
dx = start_x(fi) - start_x(fj)
dy = start_y(fi) - start_y(fj)
@constraint(model, v_x[ij] == v[i] * cos(theta[i]) - v[j] * cos(theta[j]))
@constraint(model, v_y[ij] == v[i] * sin(theta[i]) - v[j] * sin(theta[j]))
dot = dx * v_x[ij] + dy * v_y[ij]
speed_sq = v_x[ij]^2 + v_y[ij]^2
gap_sq = dx^2 + dy^2 - airspace.min_separation^2
@constraint(model, dot * (2 * k[ij] - 1) <= 0)
@constraint(model, k[ij] * (speed_sq * gap_sq - dot^2) >= 0)
end
optimize!(model)
deviation = objective_value(model)
return Resolution(deviation, value.(v), value.(theta), value.(q), value.(omega))
endUse cases
Circle air traffic
The circle above is unrealistic on purpose. It has served to stress existing resolution methods in order to explore their limits, since every pair of aircraft is in conflict and all of the conflicts happen at the same instant.
circle_resolution = minimize_deviation(circle);=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 0 variables (0%) and 0 constraints (0%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 116
bounds: lower upper range | lower upper range
0 0 52 | 0 0 52
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 88 28 0
Number of constraints: 128 | 128
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 16 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 8 88 | 8 88
quadratic: 16 84 72 | 16 84 72
nonlinear: 0 364 192 | 0 364 192
total: 16 480 208 | 16 480 208
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [2e+00, 2e+00] | [2e+00, 2e+00]
linear constraints: [7e-01, 1e+02] | [7e-01, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [1e+00, 2e+02] | [1e+00, 2e+02]
variable bounds: [2e-01, 9e+00] | [2e-01, 9e+00]
constraint bounds: [3e+00, 9e+00] | [3e+00, 9e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 6.72188e-02 3.54063 0.257175 2.668
1 0.100182 3.29569 6.05678e-02 12.015
2 3.16667 1.14224 0.738079 12.016
3 3.20274 1.12502 0.592681 12.016
4 3.24531 1.11208 0.592666 12.017
5 3.24729 1.11087 0.592664 12.017
6 3.24913 1.11042 0.592663 12.018
7 3.25972 1.10892 0.592661 12.018
8 3.26034 1.10884 0.592661 12.018
9 3.26278 1.10862 0.592661 12.019
10 3.26739 1.10847 0.592661 12.019
11 3.55745 1.01148 1.99015 12.020
12 3.57405 1.00875 1.96784 12.021
13 0.505348 2.28376 161.503 12.025
14 0.147928 2.43499 103.588 12.026
15 7.23331e-02 2.49905 40.5840 12.027
16 7.50926e-02 2.29084 28.9442 12.028
17 5.31074e-02 2.88087 84.2770 12.029
18 3.93734e-02 4.11115 83.2938 12.029
19 3.66123e-02 5.49047 68.0372 12.030
20 2.95465e-02 4.50380 42.6850 12.030
30 5.35760e-03 4.22427 0.321597 12.035
40 0.120275 3.03375 0.202689 12.040
50 0.122674 3.00132 0.262596 12.044
60 0.137777 3.00008 0.323085 12.053
70 0.169631 3.00070 0.310018 12.060
80 0.202107 3.00047 0.347616 12.065
90 0.197809 3.00002 1.00001 12.072
100 0.173090 3.09597 2.39834 12.079
110 0.233777 3.21218 0.511213 12.087
120 0.253769 4.27526 0.845448 12.093
130 4.19282e-02 5.51236 29.6873 12.101
140 5.18564e-02 1.40888 8.70441e-02 12.106
150 0.189232 1.35194 0.284203 12.112
160 0.192742 1.35943 0.295708 12.117
170 0.109371 273.666 10.0615 12.127
180 0.112658 219.237 3092.82 12.137
190 4.67609e-02 1.42147 7.80755e-02 12.143
200 0.137433 0.962878 0.310872 12.149
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 12.207
Knitro deduced that the problem is non-convex.
2 2 0.159118 DURD -inf 12.296
19 10 0.117188 FCRD -inf 12.712
100 10 0.117188 -inf 17.522
113 13 0.111862 FCRD -inf 18.125
156 14 0.105362 FCRD -inf 20.372
200 14 0.105362 -inf 22.826
300 14 0.105362 -inf 28.157
346 17 0.104537 FCRD -inf 30.450
400 17 0.104537 -inf 33.218
500 17 0.104537 -inf 38.093
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.04537083624540728e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 500 (35.187s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 50442 (3.413s)
# of gradient evaluations = 46475 (0.672s)
# of hessian evaluations = 43776 (13.092s)
# of hessian-vector evaluations = 0
# of subproblems processed = 510 (35.269s)
Total program time (secs) = 38.09299 (38.088 CPU time)
Time spent in evaluations (secs) = 17.17750
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.004s
Rounding heuristic = 5 / 4 / 0.045s
MPEC heuristic = 1 / 0 / 0.006s
Local search heuristic = 10 / 0 / 0.008s
===========================================================================
Show the full outputHide the full output
Output visualization
report_resolution(circle, circle_resolution)Circle air traffic, deviation 0.1045
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 194kt 0.970 -6.74° 0.0009 0.0138
2 196kt 0.981 -5.63° 0.0003 0.0096
3 199kt 0.996 -6.18° 0.0000 0.0116
4 201kt 1.007 -6.10° 0.0001 0.0114
5 200kt 1.002 -8.64° 0.0000 0.0227
6 209kt 1.043 -0.54° 0.0019 0.0001
7 217kt 1.087 -6.71° 0.0076 0.0137
8 180kt 0.900 +1.67° 0.0100 0.0008
Speed changes cost 0.0207 and turns cost 0.0838animate(circle, circle_resolution; show_circle=true, show_nominal=true)An intuitive solution would follow a “roundabout” pattern, whereby every aircraft deviates its course by the same angle. Knitro finds something cheaper and far less symmetric. Every aircraft turns the same way, but by very different amounts, and one of them pays with a noticeable change of speed instead. The model is a non-convex mixed-integer nonlinear program and the search stops at the node limit set above, so what you see is the best point found rather than a proven optimum.
The animation shows the speeds as well as the angles. It runs until the last aircraft has left the circle, so the order in which they cross the boundary is the order of their speeds. The ones that accelerated go first and the ones that slowed down are the last to disappear. The faint dashed rays are the nominal headings, kept underneath the solid tracks so that each deviation can be read off directly.
Sector of a circle
This use case is a variant of the circle, with 8 aircraft uniformly distributed on a quarter of the circumference of the circle, all converging on the center of the circle at the same speed.
sector = load_sector()plot_input(sector; show_circle=true)Let’s solve it and visualize the solution.
sector_resolution = minimize_deviation(sector);=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 0 variables (0%) and 0 constraints (0%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 116
bounds: lower upper range | lower upper range
0 0 52 | 0 0 52
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 88 28 0
Number of constraints: 128 | 128
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 16 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 8 88 | 8 88
quadratic: 16 84 72 | 16 84 72
nonlinear: 0 364 192 | 0 364 192
total: 16 480 208 | 16 480 208
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [2e+00, 2e+00] | [2e+00, 2e+00]
linear constraints: [1e+00, 5e+01] | [1e+00, 5e+01]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [5e+00, 1e+02] | [5e+00, 1e+02]
variable bounds: [2e-01, 5e+00] | [2e-01, 5e+00]
constraint bounds: [3e+00, 5e+00] | [3e+00, 5e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 6.72188e-02 3.12288 0.257175 0.004
1 0.112848 0.563537 4.89476e-02 0.005
2 0.135652 0.524107 0.115359 0.005
3 0.128514 0.495787 8.17697e-02 0.006
4 0.129066 0.475118 3.39921e-02 0.007
5 0.160060 0.429708 0.191115 0.007
6 0.164256 0.416835 0.130177 0.008
7 0.166564 0.406372 9.28360e-02 0.009
8 0.174601 0.397475 9.28707e-02 0.009
9 0.186263 0.384774 9.28515e-02 0.010
10 0.190239 0.380894 9.28722e-02 0.010
11 0.208287 0.367969 9.30022e-02 0.011
12 0.209154 0.366536 9.28548e-02 0.011
13 0.211561 0.364863 9.30568e-02 0.012
14 0.212107 0.364503 9.30783e-02 0.012
15 0.215467 0.362703 9.30857e-02 0.013
16 0.203050 0.395386 992.552 0.016
17 0.161823 0.599180 657.264 0.017
18 0.161165 0.752737 291.627 0.017
19 0.131711 1.14674 39.4429 0.018
20 0.129515 1.14534 100.280 0.019
30 1.59200e-02 4.06067 0.409956 0.024
40 0.181950 0.335293 0.245706 0.029
50 0.241832 0.247544 0.525840 0.035
60 0.271814 0.269771 0.354073 0.041
70 0.272839 0.268718 0.360512 0.047
80 0.271512 0.268763 0.362162 0.052
90 0.273205 0.268763 0.362165 0.057
100 4.88301e-02 0.681941 1.76959 0.067
110 0.172026 0.263045 6.24875e-02 0.073
120 0.253561 6.96614e-02 0.338845 0.078
130 0.260396 6.57391e-02 0.368594 0.084
140 0.260397 6.57383e-02 0.368601 0.089
150 0.245166 9.26053 5.79789 0.099
160 4.52224e-02 1.60279 0.686233 0.106
170 0.229788 0.120672 0.225016 0.112
180 0.259760 6.61621e-02 0.365953 0.117
190 0.260397 6.57385e-02 0.368599 0.123
200 0.260397 6.57383e-02 0.368601 0.128
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.177
Knitro deduced that the problem is non-convex.
2 2 0.170358 DURD -inf 0.440
100 2 0.170358 -inf 7.315
200 2 0.170358 -inf 15.460
300 1 0.170358 -inf 23.748
400 1 0.170358 -inf 31.756
EXIT: Satisfactory solution found.
Final Statistics for MIP
------------------------
Final objective value = 1.70357831420728800e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 402 (31.712s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 68599 (4.694s)
# of gradient evaluations = 51422 (0.753s)
# of hessian evaluations = 49266 (4.653s)
# of hessian-vector evaluations = 0
# of subproblems processed = 408 (31.786s)
Total program time (secs) = 31.92325 (31.923 CPU time)
Time spent in evaluations (secs) = 10.09964
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.006s
Rounding heuristic = 1 / 0 / 0.025s
MPEC heuristic = 1 / 0 / 0.010s
Local search heuristic = 6 / 0 / 0.004s
===========================================================================
Show the full outputHide the full output
report_resolution(sector, sector_resolution)Sector circle air traffic, deviation 0.1704
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 180kt 0.900 +8.00° 0.0100 0.0195
2 180kt 0.901 -1.85° 0.0097 0.0010
3 194kt 0.970 -3.23° 0.0009 0.0032
4 186kt 0.929 -9.70° 0.0050 0.0286
5 207kt 1.037 -3.96° 0.0014 0.0048
6 198kt 0.988 -10.09° 0.0001 0.0310
7 220kt 1.100 -4.47° 0.0100 0.0061
8 210kt 1.051 -10.92° 0.0026 0.0363
Speed changes cost 0.0398 and turns cost 0.1306animate(sector, sector_resolution; show_circle=true, show_nominal=true)We can see that in this example, changing the headings is not enough, and some aircraft accelerate a lot while others slow down significantly.
Two converging lines
Two lines of four aircraft spaced 25 NM in trail are converging towards each other at the same speed. The safety distance was chosen arbitrarily to be 20 NM. That makes the conflicts between the two lines impossible to resolve by staging alone, without deviating any path. Slotting an aircraft from one line between two consecutive aircraft of the other would need 40 NM of room, and only 25 NM is available.
lines = load_lines()plot_input(lines)lines_resolution = minimize_deviation(lines);=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 0 variables (0%) and 0 constraints (0%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 116
bounds: lower upper range | lower upper range
0 0 52 | 0 0 52
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 88 28 0
Number of constraints: 128 | 128
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 16 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 8 76 | 8 76
quadratic: 16 72 60 | 16 72 60
nonlinear: 0 364 180 | 0 364 180
total: 16 468 196 | 16 468 196
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [2e+00, 2e+00] | [2e+00, 2e+00]
linear constraints: [1e+00, 1e+02] | [1e+00, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [5e+01, 3e+02] | [5e+01, 3e+02]
variable bounds: [2e-01, 5e+00] | [2e-01, 5e+00]
constraint bounds: [5e+00, 5e+00] | [5e+00, 5e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 6.72188e-02 3.54063 0.257175 0.002
1 7.75496e-02 3.35474 0.231757 0.003
2 8.19513e-02 3.21520 0.301676 0.004
3 9.88299e-02 3.17358 7.99092e-03 0.005
4 0.108494 3.15188 6.43567e-02 0.005
5 0.114470 3.13906 6.43567e-02 0.006
6 0.176191 3.01892 6.43567e-02 0.007
7 0.178041 3.01599 6.43567e-02 0.007
8 0.182582 3.00905 6.43567e-02 0.008
9 0.184868 3.00564 6.43567e-02 0.008
10 0.185920 3.00409 6.43567e-02 0.008
11 0.186381 3.00342 6.43567e-02 0.009
12 0.186776 3.00291 6.43567e-02 0.009
13 0.192102 3.00091 7.86649e-02 0.010
14 0.191914 3.00090 8.10459e-02 0.010
15 1.72474 1.79412 0.393306 0.011
16 1.72453 1.79344 1.20179 0.011
17 0.265935 2.70494 573.628 0.014
18 8.96656e-02 2.91996 737.124 0.015
19 8.94462e-02 5.32968 4316.46 0.016
20 7.58514e-02 3.65642 413.807 0.017
30 6.34392e-02 4.15920 0.169906 0.023
40 0.202403 3.05014 0.284439 0.029
50 0.242208 3.01612 0.319779 0.037
60 0.221421 3.00068 0.361088 0.043
70 0.299582 3.00025 0.370737 0.048
80 0.361187 3.00013 0.409929 0.054
90 0.348626 3.00068 0.391890 0.061
100 0.377180 3.06702 6.87023 0.068
110 0.351701 4.51520 3.45051 0.074
120 0.352367 4.51521 1.00317 0.082
130 0.352367 4.51522 0.551276 0.086
140 7.66142e-02 5.06342 2.66001 0.096
150 0.206321 12.1907 0.226657 0.102
160 0.324940 3.94006 0.178846 0.107
170 0.310830 4.48594 0.705709 0.113
180 0.312426 3.77274 0.279723 0.118
190 0.281655 3.76977 0.412917 0.125
200 0.285380 3.71071 0.579143 0.132
300 0.270887 0.948146 0.337764 0.198
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.256
Knitro deduced that the problem is non-convex.
5 2 0.160126 FCRD -inf 1.210
6 3 0.151434 FCRD -inf 1.309
46 5 0.134583 FCRD -inf 4.582
100 5 0.134583 -inf 9.446
200 5 0.134583 -inf 18.816
300 5 0.134583 -inf 28.098
400 4 0.134583 -inf 37.162
485 5 0.131180 FCRD -inf 45.188
500 5 0.131180 -inf 46.380
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.31180492497220413e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 500 (46.100s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 82882 (5.454s)
# of gradient evaluations = 75676 (1.107s)
# of hessian evaluations = 72582 (5.869s)
# of hessian-vector evaluations = 0
# of subproblems processed = 508 (46.203s)
Total program time (secs) = 46.38002 (46.379 CPU time)
Time spent in evaluations (secs) = 12.42971
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.006s
Rounding heuristic = 4 / 4 / 0.024s
MPEC heuristic = 1 / 0 / 0.010s
Local search heuristic = 9 / 0 / 0.007s
===========================================================================
Show the full outputHide the full output
report_resolution(lines, lines_resolution)Two converging lines of 4 aircraft, deviation 0.1312
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 220kt 1.100 +10.10° 0.0100 0.0311
2 184kt 0.919 -5.52° 0.0065 0.0093
3 201kt 1.007 +1.92° 0.0001 0.0011
4 217kt 1.087 +6.88° 0.0075 0.0144
5 181kt 0.905 +2.33° 0.0090 0.0017
6 205kt 1.023 -6.34° 0.0005 0.0122
7 182kt 0.911 +6.97° 0.0079 0.0148
8 190kt 0.952 +2.98° 0.0023 0.0027
Speed changes cost 0.0438 and turns cost 0.0874animate(lines, lines_resolution; show_nominal=true)Ten random aircraft
In this use case, a set of 10 aircraft flying at a nominal speed of 200 knots was positioned and oriented at random. This initial configuration generates conflicts that have to be avoided.
scattered = load_random()plot_input(scattered)scattered_resolution = minimize_deviation(scattered);=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 0 variables (0%) and 0 constraints (0%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 175 | 175
bounds: lower upper range | lower upper range
0 0 75 | 0 0 75
free fixed | free fixed
100 0 | 100 0
cont. binary integer | cont. binary integer
130 45 0 | 130 45 0
Number of constraints: 200 | 200
eq. ineq. range | eq. ineq. range
linear: 20 0 0 | 20 0 0
quadratic: 0 45 0 | 0 45 0
nonlinear: 90 45 0 | 90 45 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 10 130 | 10 130
quadratic: 20 135 110 | 20 135 110
nonlinear: 0 585 300 | 0 585 300
total: 20 760 320 | 20 760 320
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [2e+00, 2e+00] | [2e+00, 2e+00]
linear constraints: [1e+00, 1e+02] | [1e+00, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [4e+00, 2e+02] | [4e+00, 2e+02]
variable bounds: [1e-01, 5e+00] | [1e-01, 5e+00]
constraint bounds: [7e-02, 4e+00] | [7e-02, 4e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 5.78517e-02 3.50966 0.283161 0.003
1 0.108513 3.21734 0.174435 0.004
2 5.20812 0.771571 0.363942 0.005
3 5.24085 0.760382 0.169102 0.006
4 5.26608 0.752906 0.147678 0.007
5 5.29380 0.746346 0.646228 0.008
6 5.31483 0.741860 0.644595 0.009
7 5.32653 0.740448 0.643375 0.010
8 5.35552 0.736227 0.642738 0.011
9 1.32975 1.97011 276.082 0.015
10 0.212147 2.62163 100.223 0.016
11 9.35060e-02 5.45793 190.248 0.018
12 9.31095e-02 2.65913 328.768 0.019
13 0.151613 11.0968 195.821 0.020
14 0.145621 11.0473 165.378 0.022
15 0.113608 170.205 96.6601 0.023
16 0.128919 166.674 170.887 0.024
17 0.129938 151.593 171.490 0.025
18 0.138274 44.2650 157.028 0.026
19 0.130162 2.36871 81.7029 0.027
20 9.92052e-02 2.82132 50.7336 0.028
30 2.71777e-03 15.7965 35.8772 0.034
40 0.210107 3.56185 0.132312 0.041
50 0.312771 3.14269 0.167840 0.049
60 0.369084 3.02549 0.334520 0.056
70 0.436202 3.18051 0.542480 0.064
80 0.461637 3.41303 0.358562 0.071
90 0.481939 3.42142 0.628917 0.080
100 0.485524 3.42222 0.378856 0.087
110 0.486496 3.41256 12.2216 0.098
120 0.105822 4.42356 12.5266 0.111
130 0.290364 3.94042 9.94445e-02 0.119
140 0.378266 3.80940 0.484296 0.128
150 0.375030 5.09472 0.251148 0.139
160 0.381267 5.10549 0.264279 0.147
170 0.381267 5.10549 0.264280 0.156
180 0.244895 206.693 78.2108 0.172
190 0.253812 143.443 696.281 0.188
200 7.24138e-02 1.47544 2.47650 0.200
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.303
Knitro deduced that the problem is non-convex.
2 2 4.57177e-02 FCRD -inf 0.618
19 13 2.23971e-02 FCRD -inf 1.076
76 15 2.10437e-02 FCRD -inf 7.172
100 15 2.10437e-02 -inf 9.870
200 15 2.10437e-02 -inf 21.588
300 15 2.10437e-02 -inf 32.061
400 15 2.10437e-02 -inf 42.481
500 15 2.10437e-02 -inf 53.041
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 2.10437182637335241e-02
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 500 (52.694s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 66591 (7.358s)
# of gradient evaluations = 56572 (1.282s)
# of hessian evaluations = 53640 (7.308s)
# of hessian-vector evaluations = 0
# of subproblems processed = 507 (52.812s)
Total program time (secs) = 53.04152 (53.040 CPU time)
Time spent in evaluations (secs) = 15.94811
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.006s
Rounding heuristic = 3 / 3 / 0.019s
MPEC heuristic = 1 / 0 / 0.014s
Local search heuristic = 8 / 0 / 0.019s
===========================================================================
Show the full outputHide the full output
report_resolution(scattered, scattered_resolution)10 random aircraft, deviation 0.0210
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 200kt 1.000 -0.00° 0.0000 0.0000
2 189kt 0.947 +0.20° 0.0028 0.0000
3 208kt 1.040 -1.11° 0.0016 0.0004
4 190kt 0.952 +1.05° 0.0023 0.0003
5 199kt 0.993 -3.31° 0.0001 0.0033
6 200kt 1.000 +0.00° 0.0000 0.0000
7 210kt 1.050 -5.02° 0.0025 0.0077
8 200kt 1.000 +0.00° 0.0000 0.0000
9 202kt 1.008 +0.17° 0.0001 0.0000
10 200kt 1.000 +0.00° 0.0000 0.0000
Speed changes cost 0.0093 and turns cost 0.0117animate(scattered, scattered_resolution; show_nominal=true)With ten aircraft the faint dashed headings earn their keep. The ones that never conflicted with anyone are easy to pick out, since they did not change speed or direction and their solid track lies straight on top of the dashed one.
Model extensions
The goal of this section is to show how the solutions change when we change the objective function.
Objective only taking into account the changes in speed
The new objective function is to minimize the sum of the squared differences in the aircraft speeds (rather than the heading angles).
\[ \min \sum_{i \in A} (q_i-1)^2 \]
The heading angles can still be adjusted; they are simply not considered in the objective function.
circle_speed = minimize_deviation(circle, "speed");=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 1 variable (1%) and 1 constraint (1%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 115
bounds: lower upper range | lower upper range
0 0 52 | 1 0 50
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 87 28 0
Number of constraints: 128 | 127
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 15 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 8 88 | 8 86
quadratic: 8 84 64 | 8 84 64
nonlinear: 0 364 192 | 0 364 192
total: 8 480 200 | 8 478 200
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [2e+00, 2e+00] | [2e+00, 2e+00]
linear constraints: [7e-01, 1e+02] | [7e-01, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [1e+00, 2e+02] | [1e+00, 2e+02]
variable bounds: [2e-01, 9e+00] | [2e-01, 9e+00]
constraint bounds: [3e+00, 9e+00] | [3e+00, 8e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 1.45418e-02 3.54063 3.10810e-02 0.002
1 5.19648e-02 3.29569 1.03427e-02 0.003
2 3.14434 1.14224 0.228421 0.004
3 3.18763 1.12475 0.228383 0.005
4 3.23076 1.11165 0.228238 0.005
5 3.23619 1.11009 0.228220 0.006
6 3.24073 1.10873 0.228205 0.006
7 3.24076 1.10871 0.228204 0.007
8 3.25752 1.10209 0.416178 0.008
9 3.25752 1.10204 0.417866 0.008
10 3.32814 1.07582 0.966317 0.009
11 3.32847 1.07540 0.965880 0.010
12 3.32991 1.07460 0.965042 0.010
13 3.33028 1.07439 0.964827 0.011
14 3.33161 1.07374 0.964153 0.011
15 3.33202 1.07334 0.963735 0.012
16 3.33217 1.07304 0.963421 0.012
17 3.33249 1.07277 0.963143 0.013
18 3.33267 1.07263 0.962994 0.013
19 3.33352 1.07175 0.962078 0.014
20 3.33392 1.07129 0.961599 0.014
30 2.37628e-03 6.54262 42.0037 0.023
40 2.89851e-03 5.46212 0.532703 0.028
50 7.65595e-02 3.01653 0.115204 0.033
60 7.97516e-02 3.00142 0.581047 0.040
70 7.99820e-02 3.00083 0.343723 0.048
80 7.99821e-02 3.18187 6.81848 0.054
90 7.99786e-02 3.17512 43.0317 0.062
100 7.99784e-02 3.41540 0.920719 0.070
110 7.99777e-02 4.22112 5.15029 0.077
120 7.99471e-02 4.65585 3.29023 0.082
130 0.143272 5.25325 5654.57 0.089
140 2.91638e-02 5.69935 7.93333 0.095
150 6.36936e-02 5.62310 5.63813e-02 0.101
160 7.23134e-02 5.55429 9.62196e-02 0.107
170 7.23134e-02 5.55429 9.61947e-02 0.113
180 4.87612e-02 368.490 19.0094 0.124
190 6.65134e-02 135.121 1156.61 0.135
200 2.96405e-02 1.56080 9.99740e-02 0.141
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.202
Knitro deduced that the problem is non-convex.
2 2 3.79830e-02 DURD -inf 0.337
6 4 3.17609e-02 FCRD -inf 0.487
9 5 3.55271e-15 FCRD -inf 0.575
100 5 3.55271e-15 -inf 7.188
200 5 3.55271e-15 -inf 14.875
300 4 3.55271e-15 -inf 21.324
400 4 3.55271e-15 -inf 27.296
500 3 3.55271e-15 -inf 34.709
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 3.55271367880050093e-15
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 500 (34.449s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 58420 (3.988s)
# of gradient evaluations = 51917 (4.007s)
# of hessian evaluations = 49236 (4.140s)
# of hessian-vector evaluations = 0
# of subproblems processed = 508 (34.534s)
Total program time (secs) = 34.70938 (34.708 CPU time)
Time spent in evaluations (secs) = 12.13419
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.005s
Rounding heuristic = 3 / 2 / 0.040s
MPEC heuristic = 1 / 0 / 0.009s
Local search heuristic = 8 / 0 / 0.008s
===========================================================================
Show the full outputHide the full output
report_resolution(circle, circle_speed)Circle air traffic, deviation 0.0000
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 200kt 1.000 -9.44° 0.0000 0.0271
2 200kt 1.000 -9.37° 0.0000 0.0268
3 200kt 1.000 -9.81° 0.0000 0.0293
4 200kt 1.000 -9.39° 0.0000 0.0268
5 200kt 1.000 -9.43° 0.0000 0.0271
6 200kt 1.000 -9.63° 0.0000 0.0283
7 200kt 1.000 -9.51° 0.0000 0.0276
8 200kt 1.000 -10.15° 0.0000 0.0314
Speed changes cost 0.0000 and turns cost 0.2244animate(circle, circle_speed; show_circle=true, show_nominal=true)Objective only taking into account the changes in heading angle
The new objective function is to minimize the sum of the squared differences in the aircraft heading angles (rather than the speeds).
\[ \min \sum_{i \in A} \omega _i ^ 2 \]
The speeds can still be adjusted; they are simply not considered in the objective function.
circle_angle = minimize_deviation(circle, "angle");=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 0 variables (0%) and 0 constraints (0%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 116
bounds: lower upper range | lower upper range
0 0 52 | 0 0 52
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 88 28 0
Number of constraints: 128 | 128
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 16 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 0 88 | 0 88
quadratic: 8 84 64 | 8 84 64
nonlinear: 0 364 192 | 0 364 192
total: 8 480 200 | 8 480 200
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [0e+00, 0e+00] | [0e+00, 0e+00]
linear constraints: [7e-01, 1e+02] | [7e-01, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [1e+00, 2e+02] | [1e+00, 2e+02]
variable bounds: [2e-01, 9e+00] | [2e-01, 9e+00]
constraint bounds: [3e+00, 9e+00] | [3e+00, 9e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 5.26770e-02 3.54063 0.257175 0.002
1 4.82175e-02 3.29569 5.96041e-02 0.003
2 2.10935e-02 1.14212 0.702056 0.004
3 1.23964e-02 1.12496 0.183246 0.004
4 9.92801e-03 1.11194 0.134726 0.005
5 7.29465e-03 1.10251 0.564974 0.006
6 7.46382e-03 1.10199 0.564708 0.006
7 6.60931e-03 1.09900 0.563273 0.007
8 2.20889e-02 2.28377 51.7832 0.010
9 6.81079e-02 8.76824 59.9695 0.010
10 6.05663e-02 2.61709 49.1633 0.011
11 2.03799e-02 9.41288 26.9507 0.012
12 7.50253e-02 13.1590 9.93056 0.013
13 3.62005e-02 4.06117 3.86961 0.014
14 3.27233e-02 3.51610 15.4852 0.015
15 2.78475e-02 4.32539 12.4341 0.015
16 1.88123e-02 5.79541 2.88415 0.016
17 9.71779e-03 7.12891 0.966854 0.016
18 5.50564e-03 7.24259 2.34104 0.016
19 2.37572e-05 7.53239 1.05143 0.017
20 1.69350e-05 7.10993 1.17333 0.017
30 3.87644e-02 3.00907 0.154736 0.022
40 5.70515e-02 3.37510 0.323245 0.026
50 0.110739 3.00407 0.340834 0.031
60 9.71554e-02 3.00490 0.316685 0.036
70 0.213902 3.06732 0.962980 0.043
80 0.195995 4.19685 3.65054 0.049
90 0.195744 5.08596 0.926179 0.055
100 0.173951 5.64284 3.07860 0.061
110 0.144693 5.94055 1.20528 0.068
120 0.106139 6.15259 1.36421 0.073
130 0.100639 6.12617 0.295429 0.077
140 1.42762e-02 6.12565 0.475445 0.084
150 9.37876e-02 6.10896 0.289775 0.089
160 0.100794 6.12621 0.295344 0.094
170 0.101446 6.12621 0.295343 0.100
180 7.84576e-02 332.315 20.0168 0.112
190 8.70612e-02 202.158 1679.94 0.121
200 1.21461e-02 0.973926 2.98301 0.128
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.189
Knitro deduced that the problem is non-convex.
2 2 9.59858e-02 FCRD -inf 0.341
35 12 8.74859e-02 FCRD -inf 1.549
78 13 7.52755e-02 FCRD -inf 4.241
100 13 7.52755e-02 -inf 5.546
200 13 7.52755e-02 -inf 11.303
300 13 7.52755e-02 -inf 17.274
400 12 7.52755e-02 -inf 23.259
500 12 7.52755e-02 -inf 30.421
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 7.52754803790306876e-02
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 500 (30.191s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 57403 (3.907s)
# of gradient evaluations = 51014 (0.783s)
# of hessian evaluations = 48315 (4.001s)
# of hessian-vector evaluations = 0
# of subproblems processed = 508 (30.243s)
Total program time (secs) = 30.42140 (30.420 CPU time)
Time spent in evaluations (secs) = 8.69053
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.005s
Rounding heuristic = 3 / 3 / 0.016s
MPEC heuristic = 1 / 0 / 0.008s
Local search heuristic = 8 / 0 / 0.005s
===========================================================================
Show the full outputHide the full output
report_resolution(circle, circle_angle)Circle air traffic, deviation 0.0753
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 200kt 1.002 +1.18° 0.0000 0.0004
2 180kt 0.900 -0.34° 0.0100 0.0000
3 200kt 0.998 -8.41° 0.0000 0.0216
4 198kt 0.989 -6.60° 0.0001 0.0133
5 200kt 1.000 -5.83° 0.0000 0.0104
6 202kt 1.009 -6.83° 0.0001 0.0142
7 220kt 1.100 +2.60° 0.0100 0.0021
8 180kt 0.900 -6.63° 0.0100 0.0134
Speed changes cost 0.0302 and turns cost 0.0753animate(circle, circle_angle; show_circle=true, show_nominal=true)We can now compare the solutions obtained with the new objectives to the previous ones.
animate_comparison(
circle,
[circle_resolution, circle_speed, circle_angle];
show_circle=true,
show_nominal=true,
)From left to right, the initial objective, the minimization of speed, and the minimization of angle.
We can see that if we only minimize the speeds, the aircraft keep their nominal speed exactly, for an objective of 0, and pay for it with heading changes close to the largest ones allowed, leaving the fleet almost symmetric. If we minimize only the heading angles, the deviations shrink and the aircraft adjust their speeds instead.
Here is the comparison on another example.
lines_speed = minimize_deviation(lines, "speed");=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 4 variables (3%) and 4 constraints (3%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 112
bounds: lower upper range | lower upper range
0 0 52 | 0 0 48
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 84 28 0
Number of constraints: 128 | 124
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 12 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 8 76 | 8 68
quadratic: 8 72 52 | 8 72 52
nonlinear: 0 364 180 | 0 364 180
total: 8 468 188 | 8 460 188
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [2e+00, 2e+00] | [2e+00, 2e+00]
linear constraints: [1e+00, 1e+02] | [1e+00, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [5e+01, 3e+02] | [5e+01, 3e+02]
variable bounds: [2e-01, 5e+00] | [2e-01, 5e+00]
constraint bounds: [5e+00, 5e+00] | [5e+00, 5e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 1.45418e-02 3.54063 4.13890e-02 0.005
1 3.73919e-02 3.35474 7.03742e-02 0.006
2 8.10047e-02 3.21524 4.82246e-02 0.007
3 9.95160e-02 3.17349 2.46878e-02 0.007
4 0.115326 3.14041 2.46607e-02 0.008
5 0.116868 3.13410 2.46065e-02 0.009
6 6.06150e-02 3.10141 501.595 0.012
7 6.57510e-02 3.04217 241.022 0.012
8 5.87079e-02 3.05790 65.8467 0.013
9 5.70787e-02 3.07258 25.9031 0.014
10 5.43009e-02 3.12400 9.46803 0.015
11 5.45696e-02 3.12398 2.28203 0.016
12 5.33813e-02 20.1568 39.4827 0.017
13 5.29788e-02 3.17809 9.66822 0.018
14 6.24974e-02 3.16136 11.4455 0.019
15 5.86858e-02 3.17503 4.64860 0.019
16 5.58987e-02 3.20478 4.13894 0.020
17 5.37912e-02 3.44643 3.27600 0.020
18 5.14921e-02 4.05703 2.81915 0.021
19 5.00565e-02 4.33924 7.74166 0.022
20 6.24937e-02 4.26343 0.234903 0.022
30 7.90515e-02 4.03205 7.80277e-02 0.027
40 7.96474e-02 4.12306 0.159775 0.035
50 7.99288e-02 4.12203 1.73602e-03 0.041
60 7.99972e-02 4.12182 9.08496e-04 0.046
70 8.00000e-02 4.12181 9.03046e-04 0.051
80 0.239645 3.51401 4.41885 0.061
90 8.01765e-02 4.12041 0.171894 0.069
100 7.99988e-02 4.12182 4.91948e-03 0.076
110 8.00000e-02 4.12181 1.06906e-03 0.082
120 8.00000e-02 4.12181 1.40439e-02 0.088
130 3.78380e-02 1005.84 823.201 0.102
140 3.90262e-02 461.144 2151.03 0.113
150 2.87227e-02 0.904555 1947.72 0.125
160 3.59581e-02 0.960793 8.44314e-02 0.131
170 7.89498e-02 0.774691 3.58371e-02 0.136
180 7.98760e-02 0.898706 0.172859 0.143
190 7.99559e-02 0.909398 0.298114 0.153
200 7.83407e-02 0.835443 0.646513 0.162
300 7.00305e-02 0.667219 1.10519e-02 0.234
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.248
Knitro deduced that the problem is non-convex.
12 2 8.92649e-03 FCRD -inf 3.245
100 2 8.92649e-03 -inf 12.122
200 2 8.92649e-03 -inf 22.222
300 1 8.92649e-03 -inf 33.049
400 1 8.92649e-03 -inf 44.370
EXIT: Satisfactory solution found.
Final Statistics for MIP
------------------------
Final objective value = 8.92648763319581207e-03
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 412 (45.436s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 73982 (5.029s)
# of gradient evaluations = 66165 (4.944s)
# of hessian evaluations = 63170 (5.138s)
# of hessian-vector evaluations = 0
# of subproblems processed = 417 (45.506s)
Total program time (secs) = 45.66663 (45.665 CPU time)
Time spent in evaluations (secs) = 15.11145
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.005s
Rounding heuristic = 1 / 1 / 0.012s
MPEC heuristic = 1 / 0 / 0.010s
Local search heuristic = 6 / 0 / 0.011s
===========================================================================
Show the full outputHide the full output
lines_angle = minimize_deviation(lines, "angle");=======================================
Commercial License
Artelys Knitro 16.0.0
=======================================
Knitro changing mip_method from AUTO to 1.
No start point provided -- Knitro computing one.
Knitro presolve eliminated 0 variables (0%) and 0 constraints (0%) in 0.00s.
datacheck 0
feastol 1e-06
feastol_abs 1e-06
hessian_no_f 1
mip_maxnodes 500
mip_multistart 1
mip_numthreads 1
ms_numthreads 1
numthreads 1
opttol 1e-06
opttol_abs 0.001
Knitro changing mip_root_nlpalg from AUTO to 1.
Knitro changing mip_node_nlpalg from AUTO to 1.
Knitro changing mip_branchrule from AUTO to 2.
Knitro changing mip_selectrule from AUTO to 2.
Knitro changing mip_mir from AUTO to 2.
Knitro changing mip_clique from AUTO to 0.
Knitro changing mip_zerohalf from AUTO to 0.
Knitro changing mip_liftproject from AUTO to 0.
Knitro changing mip_knapsack from AUTO to 2.
Knitro changing mip_gomory from AUTO to 0.
Knitro changing mip_cut_flowcover from AUTO to 2.
Knitro changing mip_cut_probing from AUTO to 1.
Knitro changing mip_rounding from AUTO to 3.
Knitro changing mip_heuristic_strategy from AUTO to 1.
Knitro changing mip_heuristic_feaspump from AUTO to 1.
Knitro changing mip_heuristic_misqp from AUTO to 0.
Knitro changing mip_heuristic_mpec from AUTO to 1.
Knitro changing mip_heuristic_diving from AUTO to 1926.
Knitro changing mip_heuristic_fixpropagate from AUTO to 62.
Knitro changing mip_heuristic_lns from AUTO to 0.
Knitro changing mip_heuristic_localsearch from AUTO to 1.
Knitro changing mip_pseudoinit from AUTO to 1.
WARNING: Problem appears to have nonlinear equalities and be non-convex.
The Knitro mixed integer solver is designed for convex problems.
For non-convex problems it is only a heuristic, and the reported
bounds and optimality claims cannot be verified.
Problem Characteristics | Presolved
-----------------------
Problem type: MINLP
Objective: minimize / quadratic
Number of variables: 116 | 116
bounds: lower upper range | lower upper range
0 0 52 | 0 0 52
free fixed | free fixed
64 0 | 64 0
cont. binary integer | cont. binary integer
88 28 0 | 88 28 0
Number of constraints: 128 | 128
eq. ineq. range | eq. ineq. range
linear: 16 0 0 | 16 0 0
quadratic: 0 28 0 | 0 28 0
nonlinear: 56 28 0 | 56 28 0
Number of nonzeros:
objective Jacobian Hessian | objective Jacobian Hessian
linear: 0 76 | 0 76
quadratic: 8 72 52 | 8 72 52
nonlinear: 0 364 180 | 0 364 180
total: 8 468 188 | 8 468 188
Knitro using Branch and Bound method with 1 thread.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [0e+00, 0e+00] | [0e+00, 0e+00]
linear constraints: [1e+00, 1e+02] | [1e+00, 1e+02]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [5e+01, 3e+02] | [5e+01, 3e+02]
variable bounds: [2e-01, 5e+00] | [2e-01, 5e+00]
constraint bounds: [5e+00, 5e+00] | [5e+00, 5e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 5.26770e-02 3.54063 0.257175 0.003
1 4.01587e-02 3.35474 0.229253 0.003
2 9.20997e-04 3.21527 4.17678e-02 0.004
3 4.93176e-05 3.17365 8.64113e-03 0.005
4 4.39379e-05 3.14062 8.63099e-03 0.005
5 4.38107e-05 3.13227 8.62887e-03 0.006
6 4.39476e-05 3.12933 9.04633e-03 0.006
7 4.65794e-05 3.12641 1.07590e-02 0.006
8 1.50100e-03 3.11605 1.14343e-02 0.007
9 1.51147e-03 3.11587 1.16989e-02 0.008
10 1.95582e-03 3.11502 2.22809e-02 0.008
11 1.77754e-02 3.03440 4.23211e-02 0.009
12 2.41939e-02 3.03401 7.16238e-02 0.009
13 2.67764e-02 3.03310 5.64285e-02 0.010
14 3.34906e-02 2.76246 6.60677e-02 0.011
15 3.45595e-02 2.76236 7.01969e-02 0.011
16 4.21641e-02 2.76184 9.73138e-02 0.011
17 4.69802e-02 2.76150 0.117850 0.012
18 5.27271e-02 2.76104 0.144438 0.012
19 7.18592e-02 2.76004 0.207980 0.013
20 7.22881e-02 2.75949 0.112073 0.013
30 4.80492e-03 3.62229 11.1314 0.023
40 7.95252e-03 4.03828 0.125040 0.027
50 0.245115 3.17919 0.176294 0.031
60 0.303816 3.24752 1.10154 0.039
70 0.306485 3.24905 0.180578 0.044
80 0.308360 3.24929 7.79261e-02 0.048
90 0.308423 3.24931 7.73327e-05 0.053
100 0.308425 3.24931 7.53820e-04 0.057
110 0.308425 3.24931 128.620 0.064
120 0.255160 2.79171 5465.08 0.074
130 2.63550e-02 1.12227 26.6158 0.081
140 0.214042 0.909520 0.174481 0.086
150 0.261350 0.688403 0.284920 0.093
160 0.276936 0.645614 0.444670 0.101
170 0.256760 0.642514 0.314587 0.107
180 0.257425 0.641402 0.341307 0.114
190 0.257612 0.641332 0.418131 0.120
200 0.221759 0.641332 0.267869 0.125
300 0.219575 1787.22 167.493 0.194
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.244
Knitro deduced that the problem is non-convex.
8 2 0.111754 FCRD -inf 1.943
48 3 8.28300e-02 FCRD -inf 5.910
83 4 8.08092e-02 FCRD -inf 9.113
100 4 8.08092e-02 -inf 10.664
200 4 8.08092e-02 -inf 20.047
300 4 8.08092e-02 -inf 29.304
400 4 8.08092e-02 -inf 39.180
500 4 8.08092e-02 -inf 48.275
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 8.08091971167578388e-02
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 500 (47.933s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 85456 (5.678s)
# of gradient evaluations = 77068 (1.164s)
# of hessian evaluations = 73873 (6.080s)
# of hessian-vector evaluations = 0
# of subproblems processed = 507 (48.088s)
Total program time (secs) = 48.27524 (48.274 CPU time)
Time spent in evaluations (secs) = 12.92234
Cuts statistics (gen / add)
---------------------------
Knapsack cuts = 0 / 0
Mixed-integer rounding cuts = 0 / 0
Flow-cover cuts = 0 / 0
Probing cuts = 0 / 0
Heuristics statistics (calls / successes / time)
------------------------------------------------
Feasibility pump = 1 / 0 / 0.006s
Rounding heuristic = 3 / 3 / 0.081s
MPEC heuristic = 1 / 0 / 0.010s
Local search heuristic = 8 / 0 / 0.006s
===========================================================================
Show the full outputHide the full output
animate_comparison(
lines, [lines_resolution, lines_speed, lines_angle]; show_nominal=true
)We do not see in this example a significant change in the flight paths, but we still observe that, depending on what we minimize, the results are different, and the aircraft do not cross the opposite line in the same order.
References
- “Aircraft deconfliction via Mathematical Programming: Review and insights” - Mercedes Pelegrín, Claudia d’Ambrosio (2020)
- “Resolution of Conflicts Involving Many Aircraft via Semidefinite Programming” - E. Frazzoli, Z.-H. Mao, J.-H. Oh, E. Feron (1999)