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 pyomo.environ as pyo
SOLVER_NAME = "knitroampl"minimize_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
):
model = pyo.ConcreteModel()
n = airspace.num_aircraft
fleet = airspace.fleet
model.Aircraft = pyo.RangeSet(0, n - 1)
pair_ij = [(i, j) for i in range(n) for j in range(n) if i < j]
model.Pairs = pyo.RangeSet(0, len(pair_ij) - 1)
offsets = [
(fleet[i].start_x - fleet[j].start_x, fleet[i].start_y - fleet[j].start_y)
for i, j in pair_ij
]
model.v = pyo.Var(model.Aircraft)
def theta_bounds_rule(model, i):
heading = fleet[i].heading
return (heading + fleet[i].min_turn, heading + fleet[i].max_turn)
model.theta = pyo.Var(model.Aircraft, bounds=theta_bounds_rule)
def q_bounds_rule(model, i):
return (fleet[i].min_scale, fleet[i].max_scale)
model.q = pyo.Var(model.Aircraft, bounds=q_bounds_rule)
def omega_bounds_rule(model, i):
return (fleet[i].min_turn, fleet[i].max_turn)
model.omega = pyo.Var(model.Aircraft, bounds=omega_bounds_rule)
model.k = pyo.Var(model.Pairs, within=pyo.Binary)
model.v_x = pyo.Var(model.Pairs)
model.v_y = pyo.Var(model.Pairs)
def objective_rule(model):
if objective == "speed":
return sum((model.q[i] - 1) ** 2 for i in model.Aircraft)
if objective == "angle":
return sum((model.omega[i]) ** 2 for i in model.Aircraft)
else:
return sum(
(model.q[i] - 1) ** 2 + (model.omega[i]) ** 2 for i in model.Aircraft
)
model.obj = pyo.Objective(rule=objective_rule, sense=pyo.minimize)
def q_def_rule(model, i):
return model.v[i] == model.q[i] * fleet[i].speed
model.q_def = pyo.Constraint(model.Aircraft, rule=q_def_rule)
def omega_def_rule(model, i):
return model.theta[i] == fleet[i].heading + model.omega[i]
model.omega_def = pyo.Constraint(model.Aircraft, rule=omega_def_rule)
def v_x_def_rule(model, ij):
i, j = pair_ij[ij]
first = model.v[i] * pyo.cos(model.theta[i])
second = model.v[j] * pyo.cos(model.theta[j])
return model.v_x[ij] == first - second
model.v_x_def = pyo.Constraint(model.Pairs, rule=v_x_def_rule)
def v_y_def_rule(model, ij):
i, j = pair_ij[ij]
first = model.v[i] * pyo.sin(model.theta[i])
second = model.v[j] * pyo.sin(model.theta[j])
return model.v_y[ij] == first - second
model.v_y_def = pyo.Constraint(model.Pairs, rule=v_y_def_rule)
def time_sign_rule(model, ij):
dx, dy = offsets[ij]
dot = dx * model.v_x[ij] + dy * model.v_y[ij]
return dot * (2 * model.k[ij] - 1) <= 0
model.time_sign = pyo.Constraint(model.Pairs, rule=time_sign_rule)
def aircraft_separation_rule(model, ij):
dx, dy = offsets[ij]
dot = dx * model.v_x[ij] + dy * model.v_y[ij]
speed_sq = model.v_x[ij] ** 2 + model.v_y[ij] ** 2
gap_sq = dx**2 + dy**2 - airspace.min_separation**2
return model.k[ij] * (speed_sq * gap_sq - dot**2) >= 0
model.aircraft_separation = pyo.Constraint(
model.Pairs, rule=aircraft_separation_rule
)
solver = pyo.SolverFactory(SOLVER_NAME)
solver.options["mip_multistart"] = int(multistart)
solver.options["mip_maxnodes"] = max_nodes
solver.solve(model, tee=True)
deviation = model.obj()
return Resolution(
deviation,
[model.v[i]() for i in model.Aircraft],
[model.theta[i]() for i in model.Aircraft],
[model.q[i]() for i in model.Aircraft],
[model.omega[i]() for i in model.Aircraft],
)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)Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 172 | 0 364 172
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: [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 0.101553 3.55345 0.258152 0.085
1 0.120589 3.28560 0.203434 0.087
2 3.90963 0.856145 0.365245 0.088
3 3.94844 0.843892 0.200204 0.089
4 3.97307 0.834525 0.200209 0.090
5 3.97914 0.832968 0.292517 0.091
6 4.00245 0.827085 0.287503 0.091
7 4.00471 0.826546 0.283536 0.092
8 4.00565 0.826340 0.283572 0.092
9 4.00600 0.826253 0.283709 0.093
10 4.00603 0.826243 0.283687 0.093
11 4.00640 0.826175 0.283469 0.094
12 4.01380 0.822943 0.246145 0.095
13 4.01357 0.822822 0.246145 0.096
14 4.01315 0.822755 0.246145 0.096
15 4.01282 0.822714 0.246145 0.097
16 4.00985 0.822549 0.246145 0.097
17 4.01024 0.822500 0.259152 0.098
18 4.01047 0.822498 0.260308 0.098
19 4.00218 0.822480 0.269215 0.099
20 4.00224 0.822479 0.269759 0.099
30 8.82829e-02 3.31892 39.7586 0.112
40 2.64099e-03 9.89247 0.309954 0.117
50 0.163960 3.24235 0.394019 0.122
60 0.181243 3.01298 0.341351 0.128
70 0.182600 3.00370 0.372163 0.137
80 0.237699 3.00093 0.425602 0.146
90 0.245992 3.00041 0.367414 0.154
100 0.232630 3.00004 0.359843 0.159
110 0.163602 3.00000 0.357909 0.164
120 0.157444 3.00000 0.324168 0.170
130 0.104050 3.16491 2.80922 0.183
140 0.161754 4.41072 0.150028 0.189
150 0.174327 4.75221 0.247245 0.195
160 0.192318 5.18001 0.840974 0.203
170 0.204681 5.91278 4.59367 0.214
180 0.144004 6.06269 2.01077 0.220
190 0.168696 5.87080 0.353786 0.228
200 0.168267 5.87234 0.363636 0.234
300 6.84440e-02 0.908955 0.127244 0.324
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=6'
Tree search
-----------
Nodes Best solution Best bound Gap Time
Expl | Unexpl value value (secs)
--------------- ------------- ---------- --- ------
1 1 -inf 0.397
Knitro deduced that the problem is non-convex.
4 3 0.174626 LS -inf 0.604
7 3 0.174626 DURD -inf 0.683
47 8 0.117188 FCRD -inf 1.736
447 8 0.117188 -inf 8.991
500 8 0.117188 -inf 9.960
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.17187872574968210e-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 (37.508s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 63130 (3.237s)
# of gradient evaluations = 55167 (2.279s)
# of hessian evaluations = 52308 (5.030s)
# of hessian-vector evaluations = 0
# of subproblems processed = 507 (37.600s)
Total program time (secs) = 9.96361 (38.070 CPU time)
Time spent in evaluations (secs) = 10.54587
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 = 2 / 1 / 0.040s
MPEC heuristic = 1 / 1 / 0.013s
Local search heuristic = 8 / 1 / 0.011s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.11718787257496821; optimality
gap Infinity; 500 nodes; 507 subproblem solvesShow the full outputHide the full output
Output visualization
report_resolution(circle, circle_resolution)Circle air traffic, deviation 0.1172
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 200kt 1.000 -6.93° 0.0000 0.0146
2 200kt 1.000 -6.93° 0.0000 0.0146
3 200kt 1.000 -6.93° 0.0000 0.0146
4 200kt 1.000 -6.93° 0.0000 0.0146
5 200kt 1.000 -6.93° 0.0000 0.0146
6 200kt 1.000 -6.93° 0.0000 0.0146
7 200kt 1.000 -6.93° 0.0000 0.0146
8 200kt 1.000 -6.93° 0.0000 0.0146
Speed changes cost 0.0000 and turns cost 0.1172animate(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)Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 172 | 0 364 172
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: [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 0.101553 3.12288 0.258152 0.079
1 8.51674e-02 0.609462 0.117950 0.080
2 7.43031e-02 0.567242 0.120071 0.081
3 0.101819 0.520165 0.244245 0.082
4 0.128662 0.493582 0.335064 0.082
5 0.140952 0.477708 0.350732 0.083
6 0.150174 0.466191 0.342431 0.084
7 0.153166 0.460049 0.337922 0.084
8 0.147267 0.449470 0.330155 0.085
9 0.170828 0.423597 0.335290 0.085
10 0.184764 0.412866 0.377511 0.085
11 0.185858 0.409445 0.367126 0.086
12 0.186160 0.398698 0.366972 0.086
13 0.186975 0.391068 0.359772 0.086
14 0.188108 0.387013 0.350169 0.087
15 0.177155 0.499809 698.966 0.090
16 0.135299 7.37526 603.141 0.090
17 0.128443 4.24286 604.085 0.091
18 0.127609 1.45445 168.033 0.092
19 0.126782 1.73124 215.514 0.093
20 0.129792 1.80596 21.1543 0.094
30 1.17401e-02 5.83856 0.376576 0.099
40 0.186463 0.143220 0.282999 0.103
50 0.222158 0.195982 0.330318 0.108
60 0.228033 0.189407 0.293814 0.112
70 0.228125 0.189191 0.294153 0.116
80 0.226657 0.166926 4.39348 0.124
90 7.52986e-02 0.306311 0.524930 0.132
100 0.172933 0.247637 0.144256 0.137
110 0.205806 0.188877 0.323494 0.144
120 0.206793 0.186447 0.329295 0.148
130 0.206793 0.186446 0.329297 0.153
140 0.160610 18.8356 4.20320 0.162
150 0.134996 2.19102 50.1884 0.176
160 0.168338 0.336008 5.64759e-02 0.183
170 0.204995 0.189128 0.319143 0.189
180 0.206793 0.186446 0.329295 0.196
190 0.206793 0.186446 0.329297 0.202
200 8.16909e-02 52.5912 38.2654 0.215
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=4'
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.
2 2 0.170358 FCRD -inf 0.357
52 5 0.163431 MPEC -inf 2.167
129 7 0.162876 FCRD -inf 3.923
203 8 0.162152 FCRD -inf 5.488
423 8 0.162152 -inf 9.665
441 9 0.155516 FCRD -inf 9.961
505 9 0.155516 -inf 11.051
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.55515927982071034e-01
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 505 (43.787s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 83010 (5.534s)
# of gradient evaluations = 65357 (3.576s)
# of hessian evaluations = 62804 (7.217s)
# of hessian-vector evaluations = 0
# of subproblems processed = 515 (44.122s)
Total program time (secs) = 11.05120 (44.372 CPU time)
Time spent in evaluations (secs) = 16.32761
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.073s
Rounding heuristic = 4 / 4 / 0.020s
MPEC heuristic = 2 / 1 / 0.192s
Local search heuristic = 10 / 0 / 0.011s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.15551592798207103; optimality
gap Infinity; 505 nodes; 515 subproblem solvesShow the full outputHide the full output
report_resolution(sector, sector_resolution)Sector circle air traffic, deviation 0.1555
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 220kt 1.100 +10.29° 0.0100 0.0322
2 180kt 0.900 +1.47° 0.0100 0.0007
3 195kt 0.975 -0.33° 0.0006 0.0000
4 200kt 0.999 -7.26° 0.0000 0.0160
5 184kt 0.921 -11.25° 0.0062 0.0386
6 220kt 1.100 +1.03° 0.0100 0.0003
7 208kt 1.041 -3.50° 0.0017 0.0037
8 197kt 0.987 -9.10° 0.0002 0.0252
Speed changes cost 0.0387 and turns cost 0.1168animate(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)Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 160 | 0 364 160
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: [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 0.101553 3.55345 0.258152 0.035
1 0.101556 3.36308 0.235736 0.036
2 0.115687 3.10786 0.177454 0.036
3 0.160862 3.02546 0.177299 0.037
4 0.161602 3.02396 0.177293 0.037
5 0.162632 3.01907 0.178058 0.038
6 0.163473 3.01750 0.178051 0.038
7 0.164146 3.01611 0.178045 0.038
8 0.164592 3.01523 0.178041 0.039
9 0.166093 3.01459 0.178039 0.039
10 0.171248 3.01434 0.178038 0.039
11 0.456495 2.65308 0.211362 0.040
12 0.457340 2.65268 0.211759 0.040
13 0.463492 2.65122 0.219366 0.041
14 0.496866 2.64862 0.344791 0.041
15 0.504501 2.64838 0.372903 0.041
16 0.505842 2.64825 0.374366 0.042
17 0.510924 2.64788 0.260005 0.042
18 9.22136e-02 3.04474 1171.01 0.045
19 9.47384e-02 3.02195 797.751 0.046
20 9.95374e-02 2.95756 463.043 0.046
30 2.52728e-02 2.74744 0.628846 0.051
40 0.287401 3.32438 0.222237 0.055
50 0.332716 3.30873 1.23227 0.060
60 0.351430 3.16288 0.489387 0.066
70 0.304579 3.01921 0.307253 0.071
80 0.362355 2.97684 0.338238 0.075
90 0.363115 2.99513 0.339760 0.079
100 0.363121 2.99531 0.339764 0.083
110 0.363121 2.99531 0.339763 0.087
120 0.262553 2.59873 7076.71 0.095
130 6.10330e-02 2.67561 7.94473e-02 0.101
140 6.11528e-02 2.85416 0.420556 0.109
150 0.292147 3.05772 0.220846 0.114
160 0.298238 3.08126 0.999274 0.120
170 0.296312 3.40793 0.303474 0.126
180 0.298714 3.40313 0.884142 0.132
190 0.292702 3.45336 0.383485 0.137
200 0.292705 3.45345 0.385178 0.141
300 0.312001 979.232 160.346 0.202
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.251
Knitro deduced that the problem is non-convex.
2 2 0.151563 FCRD -inf 0.471
22 3 0.147864 FCRD -inf 1.602
160 4 0.134583 FCRD -inf 7.008
316 4 0.134583 -inf 11.779
500 4 0.134583 -inf 17.086
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 1.34583377541099836e-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 (48.162s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 82232 (3.238s)
# of gradient evaluations = 75318 (2.090s)
# of hessian evaluations = 72283 (6.487s)
# of hessian-vector evaluations = 0
# of subproblems processed = 507 (48.343s)
Total program time (secs) = 17.08622 (48.470 CPU time)
Time spent in evaluations (secs) = 11.81470
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.021s
Rounding heuristic = 3 / 3 / 0.015s
MPEC heuristic = 1 / 0 / 0.028s
Local search heuristic = 8 / 0 / 0.009s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.13458337754109984; optimality
gap Infinity; 500 nodes; 507 subproblem solvesShow the full outputHide the full output
report_resolution(lines, lines_resolution)Two converging lines of 4 aircraft, deviation 0.1346
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 220kt 1.100 +8.69° 0.0100 0.0230
2 180kt 0.900 -8.54° 0.0100 0.0222
3 199kt 0.994 +0.27° 0.0000 0.0000
4 205kt 1.027 +2.65° 0.0007 0.0021
5 180kt 0.902 +4.21° 0.0097 0.0054
6 200kt 1.002 -3.88° 0.0000 0.0046
7 220kt 1.100 -9.50° 0.0100 0.0275
8 191kt 0.953 +4.82° 0.0022 0.0071
Speed changes cost 0.0427 and turns cost 0.0919animate(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)Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 265 | 0 585 265
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: [2e+00, 2e+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: [2e-01, 4e+00] | [2e-01, 4e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 0.114040 3.51462 0.288168 0.040
1 0.148656 3.25674 0.273705 0.041
2 3.94446 1.13434 0.346533 0.042
3 3.99124 1.11693 0.370207 0.043
4 4.02401 1.10296 0.296914 0.044
5 4.08366 1.08532 0.297044 0.045
6 4.10607 1.07998 0.297049 0.046
7 4.16719 1.06356 0.297081 0.047
8 4.19390 1.05652 0.297145 0.048
9 4.19505 1.05622 0.297145 0.048
10 4.23792 1.04440 0.297200 0.049
11 0.681575 2.26520 255.518 0.053
12 0.116300 8.03140 395.167 0.054
13 9.39040e-02 2.43276 320.733 0.055
14 0.127073 10.9554 60.8691 0.057
15 8.60269e-02 11.2919 21.4560 0.058
16 6.23748e-02 23.2736 8.27536 0.059
17 4.93164e-02 2.96777 3.82036 0.060
18 2.39293e-02 2.90991 2.17163 0.061
19 3.75075e-03 2.79508 1.86581 0.061
20 3.67895e-03 3.16503 0.938192 0.062
30 0.180611 3.64297 0.217721 0.068
40 0.300173 3.82760 0.163957 0.073
50 0.302342 3.82357 0.195880 0.079
60 0.302349 3.82358 0.199803 0.085
70 0.312496 3.49505 3654.84 0.093
80 0.116348 3.58580 0.123129 0.099
90 0.300234 3.82746 0.165917 0.106
100 0.302347 3.82357 0.199150 0.113
110 0.302349 3.82358 0.199802 0.120
120 0.143231 231.496 27.9043 0.136
130 7.61813e-02 0.682389 6.42875 0.149
140 0.313237 0.445423 0.138442 0.155
150 0.370705 0.591155 0.389146 0.162
160 0.373275 0.591153 0.391555 0.169
170 0.370508 0.591153 0.391577 0.176
180 0.370169 189.207 35.4009 0.193
190 6.60263e-02 0.729832 0.187923 0.202
200 0.362898 0.590540 0.313346 0.209
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.232
Knitro deduced that the problem is non-convex.
4 3 0.156338 LS -inf 0.464
7 3 0.144752 FP -inf 0.610
10 5 9.31129e-02 FCRD -inf 0.766
25 7 8.82639e-02 FCRD -inf 1.242
48 11 7.82914e-02 FCRD -inf 2.013
54 14 3.26906e-02 MPEC -inf 2.164
60 16 3.24125e-02 FCRD -inf 2.259
72 17 3.07516e-02 FCRD -inf 2.534
113 18 3.04827e-02 FCRD -inf 3.546
501 18 3.04827e-02 -inf 13.607
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 3.04827454697917943e-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 (74.027s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 77355 (11.311s)
# of gradient evaluations = 63498 (8.112s)
# of hessian evaluations = 60308 (13.180s)
# of hessian-vector evaluations = 0
# of subproblems processed = 514 (74.430s)
Total program time (secs) = 13.61597 (76.494 CPU time)
Time spent in evaluations (secs) = 32.60348
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.024s
Rounding heuristic = 8 / 7 / 0.130s
MPEC heuristic = 2 / 2 / 0.149s
Local search heuristic = 14 / 1 / 0.020s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.030482745469791794;
optimality gap Infinity; 501 nodes; 514 subproblem solvesShow 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")Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 164 | 0 364 164
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: [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.64582e-02 3.55345 3.05972e-02 0.043
1 4.56422e-02 3.28561 2.29852e-02 0.044
2 3.89829 0.856258 0.150859 0.045
3 3.94076 0.843534 0.151471 0.046
4 3.97189 0.833665 0.231361 0.046
5 4.00828 0.823478 0.334034 0.047
6 4.03889 0.813494 0.334034 0.048
7 4.04194 0.811045 0.334034 0.048
8 4.09258 0.794543 0.334034 0.049
9 4.09657 0.793827 0.292402 0.049
10 4.09642 0.793814 0.293034 0.050
11 4.09403 0.793689 0.322451 0.050
12 4.09349 0.793616 0.373147 0.050
13 4.09342 0.793573 0.387557 0.051
14 4.09355 0.793508 0.417314 0.051
15 4.09358 0.793486 0.419745 0.051
16 4.09364 0.793330 0.794787 0.052
17 4.09380 0.793196 0.794803 0.052
18 4.09510 0.791642 0.794753 0.053
19 4.09550 0.791228 0.794739 0.053
20 4.09564 0.791143 0.794737 0.053
30 1.30528e-02 3.53551 7.39230 0.063
40 2.12311e-03 7.83346 0.542139 0.068
50 7.55251e-02 3.81510 0.137004 0.072
60 7.94823e-02 3.36425 4.54556e-03 0.076
70 7.98849e-02 3.88702 1.83275 0.083
80 7.98925e-02 4.02246 5.92809e-02 0.088
90 7.99786e-02 4.09075 3.00728e-04 0.092
100 8.00000e-02 4.09351 8.27981e-08 0.095
110 0.180115 3.75356 2615.99 0.100
120 3.36009e-02 4.54391 0.200057 0.105
130 7.79705e-02 4.05206 3.99581e-03 0.109
140 7.99838e-02 4.09138 3.20093e-05 0.113
150 8.00000e-02 4.09351 5.70013e-05 0.118
160 4.66079e-02 356.678 5.11409 0.127
170 5.29243e-02 180.217 1065.93 0.136
180 2.49393e-02 1.05299 0.925685 0.141
190 5.43843e-02 1.24314 0.132702 0.146
200 6.63786e-02 1.23093 0.181564 0.150
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.195
Knitro deduced that the problem is non-convex.
2 2 1.77636e-15 FCRD -inf 0.296
192 2 1.77636e-15 -inf 8.060
392 2 1.77636e-15 -inf 16.305
EXIT: Satisfactory solution found.
Final Statistics for MIP
------------------------
Final objective value = 1.77635683940025046e-15
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 402 (27.049s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 44655 (1.484s)
# of gradient evaluations = 40401 (0.781s)
# of hessian evaluations = 38289 (3.488s)
# of hessian-vector evaluations = 0
# of subproblems processed = 408 (27.110s)
Total program time (secs) = 16.74518 (27.240 CPU time)
Time spent in evaluations (secs) = 5.75341
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.008s
Rounding heuristic = 1 / 1 / 0.002s
MPEC heuristic = 1 / 1 / 0.017s
Local search heuristic = 6 / 0 / 0.014s
===========================================================================
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.43° 0.0000 0.0271
2 200kt 1.000 -9.20° 0.0000 0.0258
3 200kt 1.000 -9.46° 0.0000 0.0273
4 200kt 1.000 -9.42° 0.0000 0.0270
5 200kt 1.000 -9.39° 0.0000 0.0269
6 200kt 1.000 -9.36° 0.0000 0.0267
7 200kt 1.000 -9.28° 0.0000 0.0262
8 200kt 1.000 -10.42° 0.0000 0.0331
Speed changes cost 0.0000 and turns cost 0.2200animate(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")Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 164 | 0 364 164
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: [3e+00, 9e+00] | [3e+00, 9e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 6.34515e-02 3.55246 0.259314 0.023
1 5.80073e-02 3.28113 0.249650 0.023
2 2.44779e-02 0.826694 0.195719 0.024
3 9.91123e-03 0.814831 9.86229e-02 0.025
4 4.31944e-03 0.807155 7.91697e-02 0.025
5 1.45101e-02 0.798745 0.134488 0.026
6 2.32096e-02 0.796683 0.160600 0.027
7 3.20242e-02 0.795174 0.191622 0.027
8 3.45805e-02 0.794838 0.199533 0.027
9 5.79916e-02 0.792640 0.256708 0.028
10 6.47278e-02 0.792158 0.270308 0.028
11 8.34950e-02 0.780122 0.298607 0.029
12 5.00144e-02 1.84013 370.053 0.032
13 4.55939e-02 2.47006 513.644 0.033
14 4.44929e-02 12.2676 504.862 0.034
15 4.11517e-02 2.66781 467.708 0.035
16 3.85276e-02 2.64791 454.427 0.035
17 3.45595e-02 2.56528 384.483 0.036
18 5.03448e-02 2.81293 145.754 0.037
19 5.05349e-02 2.91518 295.312 0.038
20 4.73339e-02 4.25757 264.180 0.038
30 1.43500e-05 7.70008 0.312537 0.042
40 4.87850e-02 3.08133 0.186466 0.046
50 7.91838e-02 3.01181 5.58258 0.053
60 9.46068e-02 3.01330 0.999917 0.057
70 0.130093 3.00508 0.254586 0.063
80 0.303251 3.59867 2.73594 0.070
90 0.286306 4.14472 0.362659 0.075
100 0.284594 4.16112 0.348380 0.079
110 0.284595 4.16112 0.348383 0.083
120 1.43662e-02 5.36924 77.7058 0.091
130 5.71694e-02 3.26001 6.16358e-02 0.095
140 4.36659e-02 3.01068 0.235622 0.100
150 4.23991e-02 3.00078 0.243205 0.105
160 0.107554 3.00725 0.383172 0.111
170 0.124136 4.22160 0.389363 0.115
180 0.118185 4.21837 0.243145 0.120
190 0.118185 4.21837 0.243145 0.125
200 7.22958e-02 389.569 8.35045 0.137
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.193
Knitro deduced that the problem is non-convex.
2 2 0.103791 FCRD -inf 0.417
34 7 0.103009 FCRD -inf 1.276
189 13 8.47207e-02 FCRD -inf 4.110
475 13 8.47207e-02 -inf 8.683
494 14 8.38519e-02 FCRD -inf 8.942
500 14 8.38519e-02 -inf 9.092
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 8.38518582257527512e-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 (38.006s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 61834 (4.240s)
# of gradient evaluations = 54440 (3.106s)
# of hessian evaluations = 51653 (5.902s)
# of hessian-vector evaluations = 0
# of subproblems processed = 509 (38.101s)
Total program time (secs) = 9.09240 (38.987 CPU time)
Time spent in evaluations (secs) = 13.24779
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.014s
Rounding heuristic = 4 / 4 / 0.029s
MPEC heuristic = 1 / 0 / 0.017s
Local search heuristic = 9 / 0 / 0.009s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.08385185822575275; optimality
gap Infinity; 500 nodes; 509 subproblem solvesShow the full outputHide the full output
report_resolution(circle, circle_angle)Circle air traffic, deviation 0.0839
aircraft speed scale turn speed cost turn cost
-----------------------------------------------------------
1 180kt 0.900 +6.60° 0.0100 0.0133
2 198kt 0.988 -5.45° 0.0001 0.0091
3 202kt 1.010 -5.45° 0.0001 0.0091
4 202kt 1.008 -8.63° 0.0001 0.0227
5 216kt 1.078 +0.80° 0.0060 0.0002
6 193kt 0.963 -0.58° 0.0013 0.0001
7 220kt 1.100 -6.43° 0.0100 0.0126
8 220kt 1.100 -7.43° 0.0100 0.0168
Speed changes cost 0.0377 and turns cost 0.0839animate(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")Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 152 | 0 364 152
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: [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.64582e-02 3.55345 9.56076e-02 0.043
1 2.89968e-02 3.36309 2.42425e-02 0.044
2 0.111032 3.10960 2.46794e-02 0.045
3 0.156531 3.02600 4.85893e-02 0.046
4 0.157587 3.02063 2.46794e-02 0.047
5 0.183146 2.97526 2.46794e-02 0.048
6 0.183163 2.97293 2.46794e-02 0.048
7 0.183763 2.97202 2.46794e-02 0.049
8 0.183862 2.97183 2.46794e-02 0.049
9 0.184084 2.97143 2.46794e-02 0.050
10 0.184119 2.97138 2.46794e-02 0.050
11 0.434919 2.65131 0.189376 0.051
12 0.435258 2.65089 0.189349 0.052
13 0.435895 2.65026 0.189315 0.052
14 0.436067 2.65001 0.315052 0.053
15 0.436173 2.64985 1.64425 0.053
16 0.436285 2.64973 8.97207 0.054
17 0.436287 2.64973 53.8639 0.054
18 0.437769 2.64972 546.907 0.055
19 5.35226e-02 26.3312 5.91634e+06 0.059
20 5.35233e-02 3.20815 1.84563e+06 0.060
30 7.01979e-02 3.26362 0.605303 0.069
40 7.84917e-02 3.06139 0.142079 0.074
50 7.96963e-02 3.05912 4.61416e-02 0.084
60 7.99392e-02 3.02893 0.351342 0.090
70 7.99391e-02 4.92743 0.197595 0.095
80 7.99995e-02 5.00416 6.15365e-06 0.099
90 8.00000e-02 5.00452 2.22570e-07 0.102
100 8.00000e-02 5.00454 1.62253e-09 0.106
110 0.187109 4.56521 68.2682 0.116
120 3.47000e-02 5.47383 0.103518 0.122
130 7.89399e-02 4.89453 0.142372 0.127
140 7.97787e-02 4.92272 4.57307e-02 0.133
150 7.99228e-02 4.63955 1.16231 0.139
160 7.99223e-02 4.60992 1.00064 0.144
170 7.99115e-02 4.16155 0.377435 0.148
180 7.90585e-02 3.58765 1.85870 0.156
190 7.31498e-02 37.3257 0.614367 0.162
200 7.17068e-02 3.32718 1.43154 0.172
300 7.64279e-02 0.699863 0.185441 0.232
400 6.98108e-02 0.667224 0.159954 0.293
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.310
Knitro deduced that the problem is non-convex.
2 2 9.19528e-03 FCRD -inf 0.478
188 3 8.92649e-03 FCRD -inf 10.328
194 3 8.92649e-03 -inf 10.595
494 3 8.92649e-03 -inf 22.003
500 3 8.92649e-03 -inf 22.182
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 8.92648763319670024e-03
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.104s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 82147 (2.613s)
# of gradient evaluations = 74106 (1.547s)
# of hessian evaluations = 71178 (5.856s)
# of hessian-vector evaluations = 0
# of subproblems processed = 506 (47.255s)
Total program time (secs) = 22.18191 (47.414 CPU time)
Time spent in evaluations (secs) = 10.01600
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.017s
Rounding heuristic = 2 / 2 / 0.013s
MPEC heuristic = 1 / 0 / 0.028s
Local search heuristic = 7 / 0 / 0.011s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.0089264876331967; optimality
gap Infinity; 500 nodes; 506 subproblem solvesShow the full outputHide the full output
lines_angle = minimize_deviation(lines, "angle")Artelys Knitro 16.0.0: mip_multistart=1
mip_maxnodes=500
=======================================
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
datacheck 0
feastol 1e-06
feastol_abs 1e-06
findiff_numthreads 1
hessian_no_f 1
hessopt 1
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 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 152 | 0 364 152
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: [5e+00, 5e+00] | [5e+00, 5e+00]
Root node relaxation
--------------------
Iter Objective Feasibility Optimality Time
error error (secs)
---- --------- ----------- ---------- ------
0 6.34515e-02 3.55246 0.259314 0.044
1 5.63365e-02 3.36176 5.96429e-02 0.045
2 2.54442e-04 2.60680 8.98484e-02 0.045
3 4.97685e-04 2.56890 1.95182e-02 0.046
4 4.42707e-04 2.52782 3.27851e-02 0.047
5 4.39296e-04 2.52607 3.28619e-02 0.047
6 4.04701e-04 2.50925 3.33907e-02 0.048
7 3.98236e-04 2.50543 4.06809e-02 0.048
8 1.27064e-02 2.49682 0.191105 0.049
9 4.61649e-02 2.42097 0.306618 0.050
10 5.17559e-02 2.41760 0.298946 0.050
11 6.10161e-02 2.41635 0.314344 0.051
12 7.45036e-02 2.41451 0.331197 0.052
13 8.70863e-02 2.40552 0.348908 0.052
14 0.108999 2.40384 0.383796 0.053
15 0.114356 2.40344 0.390236 0.053
16 0.115248 2.40297 0.390286 0.054
17 9.83628e-02 3.00670 2663.68 0.057
18 8.58805e-02 3.00690 1082.18 0.058
19 6.95255e-02 2.94126 301.053 0.059
20 7.67778e-02 2.97949 389.753 0.060
30 2.19915e-02 5.48556 0.817077 0.066
40 0.207916 4.06515 0.157956 0.070
50 0.298666 4.02009 1.50688 0.076
60 0.267089 3.73357 0.313025 0.081
70 0.269655 3.73182 1.25382 0.086
80 0.275534 3.73171 0.208255 0.090
90 0.276191 3.73137 0.204574 0.094
100 0.276196 3.73137 0.204587 0.098
110 0.276197 3.73137 0.204587 0.102
120 0.367977 3.36706 8.31452 0.111
130 0.272091 3.86631 0.287831 0.116
140 0.286896 3.83006 0.546388 0.120
150 0.299938 3.77893 0.347178 0.126
160 0.306953 3.75441 0.384711 0.130
170 0.308321 3.74976 0.425284 0.136
180 0.308383 3.74981 0.392362 0.141
190 0.308384 3.74981 0.392482 0.146
200 0.159841 388.731 1387.09 0.157
300 2.39211e-02 20.3247 2181.85 0.219
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.255
Knitro deduced that the problem is non-convex.
17 2 8.28742e-02 FCRD -inf 3.611
19 3 8.28300e-02 FCRD -inf 3.769
211 3 8.28300e-02 -inf 11.120
502 3 8.28300e-02 -inf 21.993
EXIT: Node limit reached. Integer feasible point found.
Final Statistics for MIP
------------------------
Final objective value = 8.28299585717039710e-02
Final bound value = -inf
Final optimality gap (abs / rel) = inf / inf
# of root cutting plane rounds = 1
# of restarts = 0
# of nodes processed = 502 (48.713s)
# of strong branching evaluations = 0 (0.000s)
# of function evaluations = 89123 (3.145s)
# of gradient evaluations = 80419 (1.887s)
# of hessian evaluations = 76688 (6.689s)
# of hessian-vector evaluations = 0
# of subproblems processed = 508 (49.291s)
Total program time (secs) = 21.99298 (49.411 CPU time)
Time spent in evaluations (secs) = 11.72083
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.135s
Rounding heuristic = 2 / 2 / 0.131s
MPEC heuristic = 1 / 1 / 0.155s
Local search heuristic = 7 / 0 / 0.009s
===========================================================================
WARNING: Loading a SolverResults object with a warning status into
model.name="unknown";
- termination condition: maxIterations
- message from solver: Knitro 16.0.0\x3a MIP\x3a Node limit reached.
Integer feasible point found.; objective 0.08282995857170397; optimality
gap Infinity; 502 nodes; 508 subproblem solvesShow 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)