import math
from dataclasses import dataclass, field
from plot import animate, animate_comparison, animate_input, plot_input
from report import report_resolution
KNOT = 1 / 60
@dataclass
class Aircraft:
start: tuple = None
speed: float = None
heading: float = None
min_scale: float = None
max_scale: float = None
min_turn: float = None
max_turn: float = None
@property
def start_x(self):
return self.start[0]
@property
def start_y(self):
return self.start[1]
@dataclass
class Control:
v: float
theta: float
q: float
omega: float
@property
def velocity_x(self):
return self.v * math.cos(self.theta)
@property
def velocity_y(self):
return self.v * math.sin(self.theta)
class Resolution:
def __init__(self, deviation, v, theta, q, omega):
self.deviation = deviation
self.controls = [
Control(v_i, theta_i, q_i, omega_i)
for v_i, theta_i, q_i, omega_i in zip(v, theta, q, omega, strict=True)
]
@dataclass
class Airspace:
num_aircraft: int
fleet: list = field(init=False)
name: str = None
min_separation: float = None
ray_length: float = None
xlim: tuple = None
ylim: tuple = None
def __post_init__(self):
self.fleet = [Aircraft() for _ in range(self.num_aircraft)]Aircraft 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.
import random
SKEW = 1 / 100
def _fleet_bounds(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 = -math.pi / 16, math.pi / 16
def load_circle(num_aircraft=8, ray_length=100, min_separation=5, extent=None):
airspace = Airspace(num_aircraft)
airspace.name = "Circle air traffic"
airspace.min_separation = min_separation
airspace.ray_length = ray_length
radius = ray_length / 2
reach = radius + 2 * min_separation if extent is None else extent
airspace.xlim = airspace.ylim = (-reach, reach)
_fleet_bounds(airspace)
for i, aircraft in enumerate(airspace.fleet):
bearing = 2 * math.pi / num_aircraft * i
skewed = bearing + SKEW
aircraft.start = (radius * math.cos(skewed), radius * math.sin(skewed))
aircraft.heading = math.pi + bearing
return airspace
def load_sector(
num_aircraft=8,
ray_length=120,
min_separation=5,
sector_angle=math.pi / 2,
extent=None,
):
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 = radius + 2 * min_separation if extent is None else extent
airspace.xlim = airspace.ylim = (-reach, reach)
_fleet_bounds(airspace)
for i, aircraft in enumerate(airspace.fleet):
bearing = sector_angle / (num_aircraft + 1) * (i + 1)
aircraft.start = (radius * math.cos(bearing), radius * math.sin(bearing))
aircraft.heading = math.pi + bearing
return airspace
def load_lines(line_size=4, ray_length=150, min_separation=20, spacing=25, standoff=65):
airspace = Airspace(2 * line_size)
airspace.name = f"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):
southbound = i < line_size
if southbound:
aircraft.start = (0, standoff + spacing * i)
aircraft.heading = 3 * math.pi / 2
else:
aircraft.start = (-standoff - spacing * (i - line_size), 0)
aircraft.heading = 0
return airspace
def load_random(num_aircraft=10, ray_length=70, min_separation=5, spread=50, seed=9):
airspace = Airspace(num_aircraft)
airspace.name = f"{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)
draw = random.Random(seed)
for aircraft in airspace.fleet:
x, y = draw.randint(-spread, spread), draw.randint(-spread, spread)
if x == 0:
x = 1
aircraft.start = (x, y)
heading = math.atan(y / x) + draw.randint(-180, 180) / 180 * math.pi / 4
aircraft.heading = heading + math.pi if x > 0 else heading
return airspaceEach 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
import 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.
def minimize_deviation(
airspace, objective="speed+angle", *, multistart=True, max_nodes=500
):
prob = knitro.Problem()
prob.set_param(knitro.KN_PARAM_MIP_MULTISTART, int(multistart))
prob.set_param(knitro.KN_PARAM_MIP_MAXNODES, max_nodes)
n = airspace.num_aircraft
fleet = airspace.fleet
pairs = [(i, j) for i in range(n) for j in range(n) if i < j]
v = [prob.add_variable() for _ in range(n)]
theta = [
prob.add_variable(
lb=fleet[i].heading + fleet[i].min_turn,
ub=fleet[i].heading + fleet[i].max_turn,
)
for i in range(n)
]
q = [
prob.add_variable(lb=fleet[i].min_scale, ub=fleet[i].max_scale)
for i in range(n)
]
omega = [
prob.add_variable(lb=fleet[i].min_turn, ub=fleet[i].max_turn) for i in range(n)
]
k = [prob.add_variable(vtype=knitro.KN_VARTYPE_BINARY) for _ in pairs]
v_x = [prob.add_variable() for _ in pairs]
v_y = [prob.add_variable() for _ in pairs]
if objective == "speed":
prob.add_objective(prob.nsum((q[i] - 1) ** 2 for i in range(n)))
elif objective == "angle":
prob.add_objective(prob.nsum(omega[i] ** 2 for i in range(n)))
else:
prob.add_objective(prob.nsum((q[i] - 1) ** 2 + omega[i] ** 2 for i in range(n)))
for i in range(n):
prob.add_constraint(v[i] == q[i] * fleet[i].speed)
prob.add_constraint(theta[i] == fleet[i].heading + omega[i])
for ij, (i, j) in enumerate(pairs):
fi, fj = fleet[i], fleet[j]
dx = fi.start_x - fj.start_x
dy = fi.start_y - fj.start_y
prob.add_constraint(
v_x[ij] == v[i] * prob.cos(theta[i]) - v[j] * prob.cos(theta[j])
)
prob.add_constraint(
v_y[ij] == v[i] * prob.sin(theta[i]) - v[j] * prob.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
prob.add_constraint(dot * (2 * k[ij] - 1) <= 0)
prob.add_constraint(k[ij] * (speed_sq * gap_sq - dot**2) >= 0)
prob.solve()
deviation = prob.get_attr(knitro.KN_ATTR_OBJ_VALUE)
return Resolution(
deviation,
[v[i].value for i in range(n)],
[theta[i].value for i in range(n)],
[q[i].value for i in range(n)],
[omega[i].value for i in range(n)],
)Use 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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 762 | 475
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 144 | 8 144
quadratic: 16 84 72 | 16 84 72
nonlinear: 0 308 156 | 0 308 156
total: 16 480 172 | 16 480 172
Knitro using Branch and Bound method with 8 threads.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [1e+00, 1e+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: [0e+00, 0e+00] | [3e+00, 9e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 6.72188e-02 3.54063 0.257175 0.029
1 0.100182 3.29569 6.05678e-02 0.030
2 3.16667 1.14224 0.738079 0.031
3 3.20274 1.12502 0.592681 0.033
4 3.24531 1.11208 0.592666 0.033
5 3.24729 1.11087 0.592664 0.033
6 3.24913 1.11042 0.592663 0.034
7 3.25972 1.10892 0.592661 0.034
8 3.26034 1.10884 0.592661 0.035
9 3.26278 1.10862 0.592661 0.035
10 3.26739 1.10847 0.592661 0.035
11 3.55745 1.01148 1.99015 0.037
12 3.57405 1.00875 1.96784 0.038
13 0.505348 2.28376 161.503 0.041
14 0.147928 2.43499 103.588 0.042
15 7.23331e-02 2.49905 40.5840 0.043
16 7.50926e-02 2.29084 28.9442 0.044
17 5.31074e-02 2.88087 84.2770 0.045
18 3.93734e-02 4.11115 83.2938 0.045
19 3.66123e-02 5.49047 68.0372 0.046
20 2.95464e-02 4.50378 42.6850 0.046
30 5.35761e-03 4.22427 0.321597 0.051
40 0.120275 3.03375 0.202689 0.055
50 0.122674 3.00132 0.262596 0.058
60 0.137777 3.00008 0.323085 0.065
70 0.169631 3.00070 0.310018 0.072
80 0.202107 3.00047 0.347616 0.076
90 0.197809 3.00002 1.00001 0.082
100 0.173090 3.09597 2.39840 0.088
110 0.233551 3.21205 0.490312 0.094
120 0.273110 4.93286 0.382495 0.099
130 4.50843e-02 5.68413 5.88422 0.105
140 0.260762 5.64449 0.223960 0.109
150 0.269207 5.64258 0.239608 0.113
160 0.266407 5.64258 0.239612 0.118
170 0.271061 5.64258 0.239613 0.122
180 0.156611 248.864 1.54511 0.135
190 7.31709e-02 1.04003 2.81568 0.141
200 0.130884 0.957756 0.232943 0.145
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=3'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.199
Knitro deduced that the problem is non-convex.
4 3 0.141507 FCRD -inf 0.339
7 3 0.141507 DURD -inf 0.387
108 11 0.117188 FCRD -inf 2.288
371 13 0.114315 FCRD -inf 6.333
491 13 0.114315 -inf 8.230
501 13 0.114315 -inf 8.305
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.14315207739929647e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 501 (35.798s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 62344 (3.045s)
# of gradient evaluations = 54313 (3.466s)
# of hessian evaluations = 51476 (4.597s)
# of hessian-vector evaluations = 0
# of subproblems processed = 510 (35.901s)
Total program time (secs) = 8.30568 (36.205 CPU time)
Time spent in evaluations (secs) = 11.10807
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 / 1 / 0.013s
Rounding heuristic = 4 / 3 / 0.043s
MPEC heuristic = 1 / 1 / 0.019s
Local search heuristic = 9 / 1 / 0.009s
===========================================================================
Show the full outputHide the full output
Output visualization
report_resolution(circle, circle_resolution)Circle air traffic, deviation 0.1143
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 198kt 0.991 -7.31° 0.0001 0.0163
2 197kt 0.986 -7.18° 0.0002 0.0157
3 220kt 1.100 +4.77° 0.0100 0.0069
4 188kt 0.939 -7.65° 0.0037 0.0178
5 212kt 1.061 +2.35° 0.0038 0.0017
6 196kt 0.978 +2.63° 0.0005 0.0021
7 180kt 0.900 +1.81° 0.0100 0.0010
8 202kt 1.009 -8.97° 0.0001 0.0245
Speed changes cost 0.0283 and turns cost 0.0860animate(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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 780 | 493
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 144 | 8 144
quadratic: 16 84 72 | 16 84 72
nonlinear: 0 308 156 | 0 308 156
total: 16 480 172 | 16 480 172
Knitro using Branch and Bound method with 8 threads.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [1e+00, 1e+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: [0e+00, 0e+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.002
1 0.112848 0.563537 4.89476e-02 0.003
2 0.135652 0.524107 0.115358 0.003
3 0.128534 0.495786 8.17917e-02 0.004
4 0.128951 0.475124 3.45872e-02 0.005
5 0.159908 0.428554 0.187057 0.005
6 0.166266 0.411026 0.101763 0.006
7 0.178823 0.401457 0.101709 0.006
8 0.200214 0.376631 0.101846 0.006
9 0.217983 0.364569 0.110744 0.007
10 0.221502 0.362311 0.121038 0.007
11 0.223613 0.359498 0.123736 0.008
12 0.224689 0.357433 0.127529 0.008
13 0.225041 0.357268 0.128006 0.009
14 0.229491 0.354223 0.135211 0.009
15 0.230189 0.353360 0.136228 0.009
16 0.231143 0.352618 0.137816 0.010
17 0.232692 0.351500 0.199597 0.010
18 0.209986 0.421718 306.029 0.013
19 0.166668 1.38873 353.354 0.014
20 0.148991 1.25681 112.621 0.015
30 9.63115e-03 6.69676 1.03050 0.020
40 0.173667 0.349744 0.205175 0.023
50 0.218319 0.250042 0.342974 0.028
60 0.236815 0.274006 14.9072 0.034
70 0.237494 0.265571 0.297164 0.038
80 0.236609 0.265576 0.299647 0.042
90 0.230066 0.265576 0.299652 0.046
100 4.68447e-02 0.445037 34.0932 0.054
110 0.167349 0.257489 6.89561e-02 0.058
120 0.223163 6.54410e-03 0.305318 0.063
130 0.229910 3.85154e-03 2.15274 0.067
140 0.230338 2.20418e-03 16.5621 0.073
150 0.217982 0.134912 4.00488 0.080
160 0.221364 2.35282e-02 21.1957 0.086
170 0.230054 7.77642e-04 2.97888 0.092
180 0.213906 3.30441e-02 0.267353 0.098
190 0.226976 5.67616e-03 0.380729 0.101
200 0.229365 5.60710e-03 0.383910 0.105
300 0.228835 5.61418e-03 0.383828 0.154
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=2'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.168
Knitro deduced that the problem is non-convex.
4 3 0.189932 LS -inf 0.353
7 3 0.189932 DURD -inf 0.440
13 4 0.164539 FCRD -inf 0.610
377 4 0.164539 -inf 9.850
381 5 0.147586 FCRD -inf 9.930
501 5 0.147586 -inf 12.641
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.47586032148000257e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 501 (38.988s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 85749 (2.530s)
# of gradient evaluations = 66287 (3.090s)
# of hessian evaluations = 63696 (4.576s)
# of hessian-vector evaluations = 0
# of subproblems processed = 509 (39.071s)
Total program time (secs) = 12.64074 (39.258 CPU time)
Time spent in evaluations (secs) = 10.19650
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 / 1 / 0.008s
Rounding heuristic = 3 / 2 / 0.034s
MPEC heuristic = 1 / 1 / 0.010s
Local search heuristic = 9 / 1 / 0.007s
===========================================================================
Show the full outputHide the full output
report_resolution(sector, sector_resolution)Sector circle air traffic, deviation 0.1476
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 199kt 0.997 +9.84° 0.0000 0.0295
2 209kt 1.045 +3.64° 0.0021 0.0040
3 220kt 1.100 -1.26° 0.0100 0.0005
4 188kt 0.938 +7.83° 0.0038 0.0187
5 198kt 0.988 +2.30° 0.0001 0.0016
6 180kt 0.902 +3.07° 0.0097 0.0029
7 180kt 0.900 -6.75° 0.0100 0.0139
8 220kt 1.100 -10.06° 0.0100 0.0308
Speed changes cost 0.0457 and turns cost 0.1019animate(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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 689 | 412
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 132 | 8 132
quadratic: 16 72 60 | 16 72 60
nonlinear: 0 308 144 | 0 308 144
total: 16 468 160 | 16 468 160
Knitro using Branch and Bound method with 8 threads.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [1e+00, 1e+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: [0e+00, 0e+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.003
1 7.75496e-02 3.35474 0.231757 0.004
2 8.19513e-02 3.21520 0.301676 0.005
3 9.88299e-02 3.17358 7.99092e-03 0.006
4 0.108494 3.15188 6.43567e-02 0.007
5 0.114470 3.13906 6.43567e-02 0.008
6 0.176191 3.01892 6.43567e-02 0.009
7 0.178041 3.01599 6.43567e-02 0.009
8 0.182582 3.00905 6.43567e-02 0.010
9 0.184868 3.00564 6.43567e-02 0.010
10 0.185920 3.00409 6.43567e-02 0.011
11 0.186381 3.00342 6.43567e-02 0.011
12 0.186776 3.00291 6.43567e-02 0.012
13 0.192102 3.00091 7.86649e-02 0.012
14 0.191914 3.00090 8.10459e-02 0.013
15 1.72474 1.79412 0.393306 0.014
16 1.72453 1.79344 1.20179 0.015
17 0.265935 2.70494 573.628 0.018
18 8.96656e-02 2.91996 737.124 0.020
19 8.94461e-02 5.32969 4316.46 0.021
20 7.58514e-02 3.65647 413.807 0.023
30 6.34394e-02 4.15919 0.169907 0.030
40 0.202403 3.05014 0.284432 0.036
50 0.242281 3.01614 0.319866 0.044
60 0.256630 3.00067 0.489144 0.050
70 0.299582 3.00025 0.370737 0.054
80 0.361187 3.00013 0.409929 0.059
90 0.348626 3.00068 0.391890 0.065
100 0.377181 3.06702 6.87037 0.072
110 0.351702 4.51520 3.45055 0.077
120 0.352367 4.51521 1.00314 0.084
130 0.352367 4.51522 0.307620 0.087
140 7.66138e-02 5.06326 2.65783 0.096
150 0.158067 17.9656 0.304637 0.101
160 0.324940 3.94006 0.178837 0.106
170 0.310830 4.48604 0.705698 0.111
180 0.312426 3.77274 0.279723 0.115
190 0.281655 3.76977 0.412917 0.121
200 0.285380 3.71071 0.579143 0.128
300 0.271346 0.948146 0.338830 0.184
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=3'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.220
Knitro deduced that the problem is non-convex.
3 2 0.145690 FCRD -inf 0.494
5 4 0.127936 FCRD -inf 0.525
9 5 0.127803 FCRD -inf 0.576
464 5 0.127803 -inf 10.816
504 5 0.127803 -inf 11.757
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.27803216471344072e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 504 (44.503s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 81892 (2.590s)
# of gradient evaluations = 74370 (3.490s)
# of hessian evaluations = 71348 (4.945s)
# of hessian-vector evaluations = 0
# of subproblems processed = 512 (44.632s)
Total program time (secs) = 11.75669 (44.813 CPU time)
Time spent in evaluations (secs) = 11.02561
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 / 1 / 0.009s
Rounding heuristic = 4 / 4 / 0.035s
MPEC heuristic = 1 / 1 / 0.017s
Local search heuristic = 8 / 0 / 0.007s
===========================================================================
Show the full outputHide the full output
report_resolution(lines, lines_resolution)Two converging lines of 4 aircraft, deviation 0.1278
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 184kt 0.918 -1.74° 0.0067 0.0009
2 207kt 1.036 +6.65° 0.0013 0.0135
3 183kt 0.915 -5.80° 0.0073 0.0103
4 197kt 0.987 +0.57° 0.0002 0.0001
5 220kt 1.100 -11.25° 0.0100 0.0386
6 183kt 0.913 +3.89° 0.0075 0.0046
7 201kt 1.003 -3.30° 0.0000 0.0033
8 184kt 0.919 +7.50° 0.0065 0.0171
Speed changes cost 0.0395 and turns cost 0.0883animate(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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 1194 | 766
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 220 | 10 220
quadratic: 20 135 110 | 20 135 110
nonlinear: 0 495 245 | 0 495 245
total: 20 760 265 | 20 760 265
Knitro using Branch and Bound method with 8 threads.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [1e+00, 1e+00] | [2e+00, 2e+00]
linear constraints: [1e+00, 9e+01] | [1e+00, 9e+01]
quadratic objective: [1e+00, 1e+00] | [1e+00, 1e+00]
quadratic constraints: [4e+00, 2e+02] | [4e+00, 2e+02]
variable bounds: [4e-03, 5e+00] | [4e-03, 5e+00]
constraint bounds: [0e+00, 0e+00] | [2e-01, 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.103605 3.22868 0.107332 0.004
2 4.89015 0.854109 0.218278 0.004
3 4.93462 0.841633 0.282926 0.005
4 4.95527 0.835522 0.195734 0.006
5 4.96749 0.830258 0.257397 0.007
6 5.08504 0.804497 0.892809 0.008
7 5.07960 0.802849 1.12228 0.009
8 1.38788 27.3932 161.666 0.012
9 0.271000 2.46839 150.680 0.014
10 0.132914 2.53831 258.220 0.015
11 8.81280e-02 2.55683 227.242 0.016
12 9.51359e-02 2.85589 148.502 0.018
13 8.67860e-02 2.89291 535.882 0.018
14 8.17720e-02 2.91036 101.875 0.019
15 6.68842e-02 3.08143 7.66917 0.021
16 6.33241e-02 4.90986 17.2804 0.022
17 6.18649e-02 5.61051 11.0144 0.022
18 5.14276e-02 6.47105 5.62500 0.023
19 2.95353e-02 7.93384 1.69989 0.024
20 1.92739e-02 8.42419 7.73381 0.024
30 0.119953 3.38907 0.146853 0.030
40 0.242642 3.03895 13.5696 0.040
50 0.255793 3.00799 0.268851 0.047
60 0.282533 3.00340 0.277509 0.055
70 0.289978 3.00321 12.3969 0.064
80 0.435936 3.13060 1.36347 0.071
90 0.394385 3.21974 6.17229 0.081
100 0.414259 3.40951 2.84111 0.089
110 0.418626 3.55325 4158.20 0.097
120 0.290272 4.46978 4.19519e-02 0.104
130 0.423950 3.55418 0.357555 0.110
140 0.425715 3.52916 0.374923 0.117
150 0.425715 3.52916 0.374924 0.124
160 0.176747 128.214 208.417 0.142
170 0.239290 63.9966 159.189 0.154
180 6.65064e-02 1.31120 0.163962 0.162
190 0.344037 1.08168 0.268068 0.168
200 0.348019 1.10923 0.278267 0.175
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=2'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.238
Knitro deduced that the problem is non-convex.
2 2 5.64034e-02 FCRD -inf 0.365
31 11 4.63675e-02 FCRD -inf 1.223
49 16 3.18116e-02 FCRD -inf 1.559
335 19 3.07516e-02 FCRD -inf 8.465
453 20 3.04827e-02 FCRD -inf 11.162
505 20 3.04827e-02 -inf 12.361
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 3.04827461593148996e-02
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 505 (71.002s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 80181 (9.219s)
# of gradient evaluations = 67843 (9.318s)
# of hessian evaluations = 64509 (11.487s)
# of hessian-vector evaluations = 0
# of subproblems processed = 515 (71.203s)
Total program time (secs) = 12.36088 (73.336 CPU time)
Time spent in evaluations (secs) = 30.02369
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.013s
Rounding heuristic = 6 / 6 / 0.091s
MPEC heuristic = 1 / 0 / 0.022s
Local search heuristic = 10 / 0 / 0.011s
===========================================================================
Show the full outputHide the full output
report_resolution(scattered, scattered_resolution)10 random aircraft, deviation 0.0305
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 197kt 0.984 -2.47° 0.0003 0.0019
2 204kt 1.020 -5.62° 0.0004 0.0096
3 200kt 1.000 -0.00° 0.0000 0.0000
4 195kt 0.975 -3.90° 0.0006 0.0046
5 202kt 1.011 +1.06° 0.0001 0.0003
6 200kt 1.000 -0.00° 0.0000 0.0000
7 186kt 0.931 +0.26° 0.0048 0.0000
8 210kt 1.051 -2.14° 0.0026 0.0014
9 199kt 0.996 +1.46° 0.0000 0.0007
10 205kt 1.024 -2.91° 0.0006 0.0026
Speed changes cost 0.0094 and turns cost 0.0211animate(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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 746 | 475
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 144 | 8 142
quadratic: 8 84 64 | 8 84 64
nonlinear: 0 308 156 | 0 308 156
total: 8 480 164 | 8 478 164
Knitro using Branch and Bound method with 8 threads.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [1e+00, 1e+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: [0e+00, 0e+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.003
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.004
4 3.23076 1.11165 0.228238 0.005
5 3.23619 1.11009 0.228220 0.005
6 3.24073 1.10873 0.228205 0.005
7 3.24076 1.10871 0.228204 0.006
8 3.25752 1.10209 0.416178 0.007
9 3.25752 1.10204 0.417866 0.007
10 3.32814 1.07582 0.966317 0.008
11 3.32847 1.07540 0.965880 0.008
12 3.32992 1.07459 0.965038 0.008
13 3.33029 1.07438 0.964821 0.009
14 3.33162 1.07373 0.964143 0.009
15 3.33203 1.07333 0.963720 0.009
16 3.33218 1.07302 0.963399 0.009
17 3.33252 1.07275 0.963119 0.010
18 3.33270 1.07261 0.962970 0.010
19 3.33353 1.07173 0.962060 0.010
20 3.33393 1.07127 0.961583 0.011
30 3.40037e-03 4.24762 22.5805 0.020
40 2.05340e-02 4.15865 0.199309 0.023
50 7.77868e-02 3.02745 7.48154e-02 0.027
60 7.99006e-02 3.00424 2.46743 0.031
70 7.99803e-02 3.00182 24.0718 0.037
80 7.99798e-02 3.29840 0.454522 0.044
90 7.99788e-02 3.46497 3.63568 0.049
100 7.99778e-02 3.81528 1.10936 0.053
110 7.99786e-02 4.45580 1.22996 0.060
120 7.99791e-02 5.09529 1.18738 0.066
130 7.99924e-02 5.82648 0.874951 0.071
140 7.99949e-02 6.11611 0.957559 0.077
150 7.99953e-02 6.14284 0.621196 0.082
160 7.99939e-02 6.03943 0.116986 0.085
170 7.99946e-02 6.03926 7.46980e-06 0.089
180 3.20679e-02 5.82187 0.402848 0.095
190 7.94032e-02 6.01879 1.91322e-03 0.099
200 8.00000e-02 6.03949 8.91146e-08 0.103
300 6.69471e-02 0.738192 0.140974 0.167
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=2'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.179
Knitro deduced that the problem is non-convex.
4 3 3.87627e-02 LS -inf 0.294
14 5 2.53137e-02 FCRD -inf 0.519
38 6 -3.55271e-15 FCRD -inf 1.082
419 6 -3.55271e-15 -inf 8.215
500 6 -3.55271e-15 -inf 9.806
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 (35.494s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 63175 (2.632s)
# of gradient evaluations = 55305 (3.251s)
# of hessian evaluations = 52579 (4.536s)
# of hessian-vector evaluations = 0
# of subproblems processed = 508 (35.570s)
Total program time (secs) = 9.81039 (35.966 CPU time)
Time spent in evaluations (secs) = 10.41903
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 / 1 / 0.006s
Rounding heuristic = 3 / 2 / 0.035s
MPEC heuristic = 1 / 1 / 0.009s
Local search heuristic = 8 / 1 / 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.12° 0.0000 0.0253
2 200kt 1.000 -9.14° 0.0000 0.0255
3 200kt 1.000 -9.06° 0.0000 0.0250
4 200kt 1.000 -9.12° 0.0000 0.0253
5 200kt 1.000 -9.15° 0.0000 0.0255
6 200kt 1.000 -9.08° 0.0000 0.0251
7 200kt 1.000 -9.37° 0.0000 0.0267
8 200kt 1.000 -10.15° 0.0000 0.0314
Speed changes cost 0.0000 and turns cost 0.2099animate(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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 738 | 475
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 144 | 0 144
quadratic: 8 84 64 | 8 84 64
nonlinear: 0 308 156 | 0 308 156
total: 8 480 164 | 8 480 164
Knitro using Branch and Bound method with 8 threads.
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: [0e+00, 0e+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.009
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.1589 9.93056 0.013
13 3.62005e-02 4.06117 3.86961 0.014
14 3.27233e-02 3.51610 15.4852 0.014
15 2.78475e-02 4.32539 12.4341 0.014
16 1.88123e-02 5.79541 2.88415 0.015
17 9.71776e-03 7.12891 0.966853 0.015
18 5.50561e-03 7.24259 2.34103 0.016
19 2.37570e-05 7.53239 1.05143 0.016
20 1.69350e-05 7.10993 1.17332 0.016
30 3.87644e-02 3.00907 0.154736 0.020
40 5.70515e-02 3.37510 0.323245 0.024
50 0.110739 3.00407 0.340834 0.028
60 9.71554e-02 3.00490 0.316685 0.032
70 0.213902 3.06732 0.962980 0.038
80 0.195995 4.19685 3.65054 0.044
90 0.195744 5.08596 0.926176 0.048
100 0.173951 5.64284 3.07867 0.054
110 0.144876 5.94350 1.21999 0.059
120 0.106285 6.15267 1.26553 0.063
130 0.100614 6.12617 0.295422 0.067
140 1.37641e-02 6.16627 0.846082 0.072
150 8.79448e-02 6.09330 0.273867 0.076
160 0.100663 6.12621 0.295346 0.081
170 0.101521 6.12621 0.295344 0.085
180 7.84601e-02 333.073 20.0131 0.094
190 8.63582e-02 206.828 2593.80 0.103
200 7.89572e-03 1.15934 5.37403 0.109
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=2'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.159
Knitro deduced that the problem is non-convex.
2 2 0.135908 FCRD -inf 0.235
6 2 0.135908 MPEC -inf 0.318
15 5 9.96958e-02 FCRD -inf 0.583
71 6 9.81345e-02 FCRD -inf 1.837
296 13 8.74859e-02 FCRD -inf 5.690
452 13 8.74859e-02 -inf 7.985
500 13 8.74859e-02 -inf 8.747
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 8.74859078490169834e-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 (35.235s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 60715 (2.973s)
# of gradient evaluations = 52827 (3.481s)
# of hessian evaluations = 50041 (4.641s)
# of hessian-vector evaluations = 0
# of subproblems processed = 509 (35.311s)
Total program time (secs) = 8.74983 (36.301 CPU time)
Time spent in evaluations (secs) = 11.09564
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.009s
Rounding heuristic = 4 / 4 / 0.024s
MPEC heuristic = 1 / 1 / 0.009s
Local search heuristic = 10 / 0 / 0.008s
===========================================================================
Show the full outputHide the full output
report_resolution(circle, circle_angle)Circle air traffic, deviation 0.0875
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 220kt 1.100 -7.57° 0.0100 0.0175
2 180kt 0.900 +2.55° 0.0100 0.0020
3 196kt 0.981 -8.27° 0.0004 0.0209
4 220kt 1.100 +2.02° 0.0100 0.0012
5 180kt 0.900 -9.78° 0.0100 0.0292
6 208kt 1.041 +4.21° 0.0017 0.0054
7 196kt 0.979 +3.33° 0.0004 0.0034
8 185kt 0.924 +5.13° 0.0058 0.0080
Speed changes cost 0.0483 and turns cost 0.0875animate(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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 673 | 412
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 132 | 8 124
quadratic: 8 72 52 | 8 72 52
nonlinear: 0 308 144 | 0 308 144
total: 8 468 152 | 8 460 152
Knitro using Branch and Bound method with 8 threads.
Initial points
--------------
No initial point provided for the root node relaxation.
No primal point provided for the MIP.
Coefficient range:
linear objective: [1e+00, 1e+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: [0e+00, 0e+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.002
1 3.73919e-02 3.35474 7.03742e-02 0.003
2 8.10047e-02 3.21524 4.82246e-02 0.003
3 9.95160e-02 3.17349 2.46878e-02 0.004
4 0.115326 3.14041 2.46607e-02 0.004
5 0.116868 3.13410 2.46065e-02 0.005
6 6.06150e-02 3.10141 501.595 0.007
7 6.57510e-02 3.04217 241.022 0.007
8 5.87079e-02 3.05790 65.8467 0.008
9 5.70787e-02 3.07258 25.9031 0.009
10 5.43009e-02 3.12400 9.46803 0.010
11 5.45696e-02 3.12398 2.28203 0.010
12 5.33813e-02 20.1568 39.4931 0.011
13 5.29813e-02 3.17803 9.72065 0.012
14 6.25045e-02 3.16128 11.5426 0.012
15 5.86843e-02 3.17503 4.64306 0.013
16 5.58975e-02 3.20476 4.15142 0.013
17 5.37845e-02 3.44702 3.28764 0.013
18 5.14899e-02 4.05638 2.82405 0.014
19 5.00609e-02 4.33822 7.37330 0.015
20 6.25095e-02 4.26311 0.234088 0.015
30 7.90515e-02 4.03205 7.80277e-02 0.019
40 7.96474e-02 4.12306 0.159775 0.025
50 7.99288e-02 4.12203 1.73602e-03 0.029
60 7.99972e-02 4.12182 9.08496e-04 0.033
70 8.00000e-02 4.12181 9.03037e-04 0.036
80 0.221698 3.58790 3.57967 0.046
90 7.99020e-02 4.12083 0.191752 0.051
100 7.98585e-02 4.12204 0.139016 0.056
110 7.98076e-02 3.88751 0.418695 0.064
120 7.95656e-02 4.17355 39.0574 0.068
130 7.27659e-02 3.98620 0.103809 0.072
140 7.97913e-02 3.73646 0.377700 0.080
150 7.99097e-02 3.73173 0.138244 0.084
160 7.99929e-02 3.73140 0.147125 0.090
170 7.99986e-02 3.73138 1.02438e-02 0.096
180 7.99995e-02 3.73137 0.183991 0.105
190 3.62969e-02 546.096 317.613 0.113
200 3.42794e-02 25.7051 233.519 0.121
300 5.55842e-02 0.612244 0.103360 0.177
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=2'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.209
Knitro deduced that the problem is non-convex.
2 2 8.88468e-03 FCRD -inf 0.422
192 2 8.88468e-03 -inf 10.179
392 2 8.88468e-03 -inf 20.435
EXIT: Satisfactory solution found.
Final Statistics for MIP
------------------------
Final objective value = 8.88467672725745672e-03
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 402 (35.230s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 66565 (1.197s)
# of gradient evaluations = 59857 (2.138s)
# of hessian evaluations = 57408 (3.403s)
# of hessian-vector evaluations = 0
# of subproblems processed = 407 (35.342s)
Total program time (secs) = 20.95801 (35.511 CPU time)
Time spent in evaluations (secs) = 6.73812
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 / 1 / 0.005s
MPEC heuristic = 1 / 0 / 0.047s
Local search heuristic = 6 / 0 / 0.004s
===========================================================================
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.
concurrent_evals 0
feastol 1e-06
feastol_abs 1e-06
mip_maxnodes 500
mip_multistart 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 expressions: 665 | 412
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 132 | 0 132
quadratic: 8 72 52 | 8 72 52
nonlinear: 0 308 144 | 0 308 144
total: 8 468 152 | 8 468 152
Knitro using Branch and Bound method with 8 threads.
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: [0e+00, 0e+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.002
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.004
4 4.39379e-05 3.14062 8.63099e-03 0.005
5 4.38107e-05 3.13227 8.62887e-03 0.005
6 4.39476e-05 3.12933 9.04633e-03 0.006
7 4.65794e-05 3.12641 1.07590e-02 0.006
8 1.50120e-03 3.11605 1.14353e-02 0.007
9 1.51575e-03 3.11580 1.18066e-02 0.007
10 1.93462e-03 3.11502 2.18571e-02 0.007
11 1.75872e-02 3.03448 4.19092e-02 0.008
12 2.41674e-02 3.03407 7.19739e-02 0.009
13 2.68664e-02 3.03317 5.61009e-02 0.009
14 3.37145e-02 2.70804 6.62414e-02 0.010
15 4.28439e-02 2.70730 9.28552e-02 0.010
16 6.32686e-02 2.70618 0.146580 0.011
17 6.69071e-02 2.70545 0.144305 0.012
18 2.09821e-02 3.04561 162.410 0.014
19 1.82101e-02 3.02317 149.903 0.015
20 4.12031e-02 2.83444 81.9839 0.016
30 2.85318e-03 3.86900 12.0619 0.022
40 8.48628e-02 4.62622 3.34222e-02 0.026
50 0.292459 4.29087 1.86350 0.033
60 0.241959 4.37091 4.50241 0.041
70 0.240272 4.50547 1.15762 0.047
80 0.234951 4.51439 0.637594 0.051
90 0.252216 4.51532 0.368095 0.055
100 0.252474 4.51522 0.369515 0.059
110 0.252484 4.51521 0.369566 0.062
120 9.93774e-02 4.21148 41.8915 0.071
130 2.88283e-02 4.12156 0.407726 0.076
140 2.92181e-02 3.34937 0.220527 0.083
150 0.296078 3.70620 0.284164 0.087
160 0.305516 3.71849 0.375416 0.093
170 0.303875 3.71856 0.367957 0.098
180 0.303382 3.71894 7.54356 0.104
190 0.251627 3.71860 2.62683 0.109
200 0.244254 3.71904 0.283273 0.114
300 0.208516 0.556798 5.38064 0.172
WARNING: Evaluations seem expensive and concurrent evaluations are disabled,
please consider reducing the number of threads for the branch-and-bound
using option 'mip_numthreads'.
Recommended value: 'mip_numthreads=2'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.191
Knitro deduced that the problem is non-convex.
2 2 8.79343e-02 FCRD -inf 0.407
192 2 8.79343e-02 -inf 10.299
245 3 8.28300e-02 FCRD -inf 13.118
464 3 8.28300e-02 -inf 21.032
501 2 8.28300e-02 -inf 22.672
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 8.28299605301864722e-02
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 501 (44.843s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 82180 (1.601s)
# of gradient evaluations = 75117 (2.767s)
# of hessian evaluations = 72084 (4.305s)
# of hessian-vector evaluations = 0
# of subproblems processed = 507 (45.046s)
Total program time (secs) = 22.67240 (45.266 CPU time)
Time spent in evaluations (secs) = 8.67248
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.013s
Rounding heuristic = 2 / 2 / 0.060s
MPEC heuristic = 1 / 0 / 0.021s
Local search heuristic = 7 / 0 / 0.005s
===========================================================================
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)