A Boeing 787-9 Dreamliner flying nonstop from Newark Liberty International Airport (EWR) to Leonardo da Vinci-Fiumicino Airport (FCO) could need $68K in jet fuel over the 8.5-hour flight. Adjusting the flight path for wind conditions could reduce fuel consumption and possibly save a few thousand dollars.

Firms like Jeppesen have offerings in this space, but Scikit-decide, together with a narrow- and wide-body fuel consumption model built by a professor at the Delft University of Technology and wind data from NOAA, offer an open source solution.

Scikit-decide has been in development for six years. It's a framework for reinforcement learning, automated planning and scheduling. The project can optimise flight paths, re-organise airline workforce schedules and calculate drone swarm paths.

OpenAP is an aircraft performance model and toolkit developed by Dr. Junzi Sun. Dr. Sun has a PhD in air traffic management and, among many other things, teaches a course on the subject as a tenured assistant professor at TU Delft in the Netherlands.

Scikit-decide's optimal flight path solver can be configured to use different fuel consumption models. In this post, I'll compare two flight paths flown using the Airbus A320 and OpenAP's fuel consumption model.

My Workstation

I'm using a 5.7 GHz AMD Ryzen 9 9950X CPU. It has 16 cores and 32 threads and 1.2 MB of L1, 16 MB of L2 and 64 MB of L3 cache. It has a liquid cooler attached and is housed in a spacious, full-sized Cooler Master HAF 700 computer case.

The system has 96 GB of DDR5 RAM clocked at 4,800 MT/s and a 5th-generation, Crucial T700 4 TB NVMe M.2 SSD which can read at speeds up to 12,400 MB/s. There is a heatsink on the SSD to help keep its temperature down. This is my system's C drive.

The system is powered by a 1,200-watt, fully modular Corsair Power Supply and is sat on an ASRock X870E Nova 90 Motherboard.

I'm running Ubuntu 24 LTS via Microsoft's Ubuntu for Windows on Windows 11 Pro. In case you're wondering why I don't run a Linux-based desktop as my primary work environment, I'm still using an Nvidia GTX 1080 GPU which has better driver support on Windows and ArcGIS Pro only supports Windows natively.

Installing Prerequisites

I'll use Python 3.12 along with jq in this post.

sudo add-apt-repository ppa:deadsnakes/ppa

$ sudo apt update

$ sudo apt install \

jq \

python3-pip \

python3.12-venv

I'll set up a Python Virtual Environment and install scikit-decide, along with the OpenAP open aircraft performance model and OpenTop, a flight trajectory toolkit that was also developed by Dr. Sun.

python3 -m venv ~/.flight_planning

$ source ~/.flight_planning/bin/activate

$ pip install \

'scikit-decide[all]' \

'openap[all]' \

opentop

The above will need at least 8 GB of storage capacity. These are the packages that were installed.

pip install pipdeptree

$ pipdeptree -d0

lz4==4.4.5

openevolve==0.3.2

opentop==2.6.0

pip==24.0

pipdeptree==4.2.5

plado==0.1.6

pygeodesy==26.9.9

pygrib==2.1.8

pyRDDLGym-gurobi==0.2

pyRDDLGym-jax==3.1

pyRDDLGym-rl==0.2

pytz==2026.3.post1

ray==2.37.0

rddlrepository==2.2

sb3_contrib==2.3.0

scikit-decide==1.1.1

scikit-image==0.26.0

tensorboardX==2.6.5

torch-geometric==2.8.0.post1

typer==0.27.2

unified-planning==1.2.0

up-enhsp==0.0.27

up_fast_downward==0.5.2

up-pyperplan==1.1.0

z3-solver==5.1.0.0

I'll use DuckDB, along with its H3, JSON, Lindel, Parquet and Spatial extensions in this post.

cd ~

$ wget -c https://github.com/duckdb/duckdb/releases/download/v1.5.4/duckdb_cli-linux-amd64.zip

$ unzip -j duckdb_cli-linux-amd64.zip

$ chmod +x duckdb

$ ~/duckdb

INSTALL h3 FROM community;

INSTALL lindel FROM community;

INSTALL json;

INSTALL parquet;

INSTALL spatial;

I'll set up DuckDB to load every installed extension each time it launches.

vi ~/.duckdbrc

The maps in this post were rendered with QGIS version 4.2.1. QGIS is a desktop application that runs on Windows, macOS and Linux. The application has grown in popularity in recent years and has ~22M application launches from users all around the world each month.

The boundaries and place names were sourced from Natural Earth. Maritime Boundaries were sourced from Marine Regions.

OpenAP's Aircraft Types

I'll first clone the OpenAP repository.

git clone https://github.com/junzis/openap

Excluding unit tests and utility scripts, there are 3,369 lines of Python in this package.

OpenAP's model relies on a large number of datasets that are packaged with its codebase. These cover a wide variety of aircraft. Below are the aircraft manufacturer counts.

grep -ho 'aircraft: .*[a-z] ' \

openap/data/aircraft/*.yml \

| cut -d' ' -f2 \

| sort \

| uniq -c \

| sort -rn

These are the properties for the Airbus A380-800.

cat openap/data/aircraft/a388.yml

aircraft: Airbus A380-800

mtow: 560000

mlw: 386000

oew: 277000

mfc: 320000

vmo: 340

mmo: 0.89

ceiling: 13100

pax:

max: 853

low: 410

high: 620

fuselage:

length: 72.72

height: 8.41

width: 7.14

wing:

area: 845

span: 79.75

mac: null

sweep: 33.5

t/c: 0.08

flaps:

type: single-slotted

area: null

bf/b: null

lambda_f: 0.900

cf/c: 0.150

Sf/S: 0.150

cruise:

height: 12800

mach: 0.85

range: 14800

engine:

type: turbofan

mount: wing

number: 4

default: GP7270

options:

A380-841: Trent 970-84

A380-842: Trent 972-84

A380-861: GP7270

drag:

cd0: 0.016

k: 0.050

e: 0.855

gears: 0.012

These are its drag coefficients.

cat openap/data/dragpolar/a388.yml

aircraft: Airbus A380-800

clean:

cd0: 0.016

k: 0.050

e: 0.855

gears: 0.012

flaps:

lambda_f: 0.900

cf/c: 0.150

Sf/S: 0.150

These are some additional properties.

echo "import pandas as pd; print(

pd.read_fwf('openap/data/wrap/a388.txt')

.to_csv(index=False))" \

| python3 \

| ~/duckdb \

-c '.maxwidth 150' \

-c "SELECT * EXCLUDE(parameters),

parameters: SPLIT(parameters, '|')

FROM READ_CSV('/dev/stdin')"

These are the aircraft type synonyms list.

~/duckdb -c "FROM READ_CSV('/dev/stdin')" \

< openap/data/aircraft/_synonym.csv

Aircraft Engines

Aircraft often have the option of at least two different engines to choose from. There are 427 engines listed in this package's dataset.

wc -l openap/data/engine/engines.csv # 427

These are the details for the Trent 970-84.

echo "FROM 'openap/data/engine/engines.csv'

WHERE name = 'Trent 970-84'

LIMIT 1" \

| ~/duckdb -json \

| jq -S .

[

{

"bpr": 8.45,

"cruise_alt": null,

"cruise_mach": null,

"cruise_sfc": null,

"cruise_thrust": null,

"ei_co_app": 1.16,

"ei_co_co": 0.31,

"ei_co_idl": 13.38,

"ei_co_to": 0.32,

"ei_hc_app": 0.08,

"ei_hc_co": 0.12,

"ei_hc_idl": 0.04,

"ei_hc_to": 0.02,

"ei_nox_app": 12.09,

"ei_nox_co": 29.42,

"ei_nox_idl": 5.44,

"ei_nox_to": 38.29,

"ff_app": 0.72,

"ff_co": 2.157,

"ff_idl": 0.255,

"ff_to": 2.605,

"fuel_lto": 965.0,

"manufacturer": "Rolls-Royce plc",

"max_thrust": 338700.0,

"name": "Trent 970-84",

"pr": 38.0,

"type": "TF",

"uid": "18RR081"

}

]

These are the engine manufacturer counts.

~/duckdb

CREATE OR REPLACE TABLE a AS

FROM 'openap/data/engine/engines.csv';

SELECT COUNT(*),

manufacturer

FROM a

GROUP BY 2

ORDER BY 1 DESC;

These are the engine-type counts for Turbofan (TF), Mixed-flow Turbofan (MTF), Turboprop (TP) and Piston (PS) engines in this dataset.

SELECT COUNT(*),

type

FROM a

GROUP BY 2

ORDER BY 1 DESC;

This is the engine list ranked by their maximum thrust.

SELECT manufacturer,

name,

type,

max_thrust

FROM a

ORDER BY 4 DESC

LIMIT 25;

These are the fuel model defaults and overrides.

~/duckdb -c "FROM READ_CSV('/dev/stdin')" \

< openap/data/fuel/fuel_models.csv

Toulouse to Berlin

Below, I'll find an optimal flight path from Toulouse-Blagnac Airport (LFBO / TLS) to Berlin Brandenburg Airport (EDDB / BER).

python3

import numpy as np

from openap.aero import cas2mach, ft, kts

from openap.extra.nav import airport

from pygeodesy.ellipsoidalVincenty import LatLon

from skdecide.hub.domain\

.flight_planning\

.aircraft_performance\

.bean.aircraft_state \

import AircraftState

from skdecide.hub.domain\

.flight_planning\

.aircraft_performance\

.performance.performance_model_enum \

import PerformanceModelEnum

from skdecide.hub.domain\

.flight_planning\

.aircraft_performance\

.performance.phase_enum \

import PhaseEnum

from skdecide.hub.domain\

.flight_planning\

.aircraft_performance\

.performance.rating_enum \

import RatingEnum

from skdecide.hub.domain\

.flight_planning\

.domain \

import FlightPlanningDomain, \

WeatherDate

from skdecide.hub.domain\

.flight_planning\

.flightplanning_utils \

import plot_network_adapted

from skdecide.hub.solver.astar import Astar

The heuristic parameter can be either "time", "distance", "lazy_fuel", "lazy_time", or None. If nothing is passed, A* will use a Dijkstra-like search algorithm.

origin = "LFPG"

destination = "LFBO"

aircraft = "A320"

weather_date = WeatherDate(day=1, month=5, year=2026)

heuristic = "lazy_fuel"

cost_function = "fuel"

acState = AircraftState(

model_type="A320",

performance_model_type=PerformanceModelEnum.OPENAP,

gw_kg=80_000,

zp_ft=10_000,

mach=cas2mach(250 * kts, h=10_000 * ft),

phase=PhaseEnum.CLIMB,

rating_level=RatingEnum.MCL,

cg=0.3)

domain_factory = lambda: FlightPlanningDomain(

aircraft_state=acState,

mach_cruise=0.78,

mach_climb=0.7,

mach_descent=0.65,

nb_forward_points=20,

nb_lateral_points=10,

nb_climb_descent_steps=5,

flight_levels_ft=list(np.arange(30_000, 38_000 + 2_000, 2_000)),

graph_width="medium",

origin=LatLon(43.629444, 1.363056),

destination="EDDB",

objective=cost_function,

heuristic_name=heuristic,

weather_date=weather_date)

domain = domain_factory()

When the above runs, if weather data hasn't been fetched from NOAA and if the date of the flight is within the past six months, GRB2 files will be downloaded.

du -hs ~/skdecide_data/weather/grib/nowcast/*/*.grb2

Each file has data covering the entire planet. These are the contents of gfs_4_20260501_1800_000.grb2 rendered on a globe in QGIS.

This is the solver's altitude and geographical search space.

plot_network_adapted(

graph=domain.network,

p0=LatLon(43.629444, 1.363056),

p1=LatLon(

airport("EDDB")["lat"],

airport("EDDB")["lon"],

airport("EDDB")["alt"] * ft))

This is the optimal flight path according to the solver.

solver = Astar(

domain_factory=domain_factory,

heuristic=lambda d, s: d.heuristic(s),

parallel=False)

solver.solve()

domain.custom_rollout(solver=solver, make_img=True)

({'time': 7666.281474928903, 'fuel': 5855.093906205222}, None)

I'll format each of the flight plan's steps so they're easier to read.

domain.observation.trajectory.to_csv('TLS-BER.csv', index=None)

~/duckdb

SELECT phase: UPPER(phase),

time_: ts::INT,

alt: alt::INT,

mass: mass::INT,

mach: ROUND(mach, 2),

cas: cas::INT,

fuel: fuel::INT,

geom: ST_POINT(lon, lat)

FROM 'TLS-BER.csv'

ORDER BY ts;

I'll export the flight plan to Parquet and render it on top of the ground-level wind data in QGIS.

COPY (

SELECT * EXCLUDE(lon, lat),

geometry: ST_POINT(lon, lat)

FROM 'TLS-BER.csv'

ORDER BY ts

) TO 'TLS-BER.parquet' (

FORMAT 'PARQUET',

CODEC 'ZSTD',

COMPRESSION_LEVEL 22,

ROW_GROUP_SIZE 15000);

Toulouse to Warsaw

Below, I'll find an optimal flight path from Toulouse-Blagnac Airport (LFBO / TLS) to Warsaw Chopin Airport (EPWA / WAW).

The initial target altitude will be much higher than in the previous example. The result is a flight that is able to take a much more direct route.

acState = AircraftState(

model_type="A320",

performance_model_type=PerformanceModelEnum.OPENAP,

gw_kg=80_000,

zp_ft=18000.0,

mach=cas2mach(250 * kts, h=10_000 * ft),

phase=PhaseEnum.CLIMB,

rating_level=RatingEnum.MCL,

cg=0.3,

x_graph=5,

y_graph=5,

z_graph=10)

domain_factory = lambda: FlightPlanningDomain(

aircraft_state=acState,

mach_cruise=0.78,

mach_climb=0.7,

mach_descent=0.65,

nb_forward_points=20,

nb_lateral_points=10,

nb_climb_descent_steps=5,

flight_levels_ft=list(np.arange(30_000, 38_000 + 2_000, 2_000)),

graph_width="medium",

origin=LatLon(43.629444, 1.363056),

destination="EPWA",

objective=cost_function,

heuristic_name=heuristic,

weather_date=weather_date)

domain = domain_factory()

solver = Astar(

domain_factory=domain_factory,

heuristic=lambda d, s: d.heuristic(s),

parallel=False)

solver.solve()

domain.custom_rollout(solver=solver, make_img=True)

({'time': 6153.660431613251, 'fuel': 5600.171145693044}, None)

Warsaw is 500 KM further away from Toulouse than Berlin. But the faster climb to cruising altitude under the given wind conditions meant the aircraft could take a more direct route. It made it to Warsaw almost 45 minutes faster and only needed 76% of the fuel that the Berlin flight needed.

These are the steps in the above flight plan.

domain.observation.trajectory.to_csv('TLS-WAW.csv', index=None)

~/duckdb

SELECT phase: UPPER(phase),

time_: ts::INT,

alt: alt::INT,

mass: mass::INT,

mach: ROUND(mach, 2),

cas: cas::INT,

fuel: fuel::INT,

geom: ST_POINT(lon, lat)

FROM 'TLS-WAW.csv'

ORDER BY ts;

Airbus A320 vs Boeing 737

OpenTop can be paired with OpenAP and used to figure out flight trajectories between two airports.

Its optimiser requires a grid cost file. I'll first download an example 142 MB NetCDF file provided by the project.

wget https://opendap.4tu.nl/thredds/fileServer/data2/djht/bea8a3fe-e34c-4598-9f94-c5a5c63348e5/1/contrail_original.nc

The cost file can be either in Casadi or Parquet format. I worked from an example in its documentation, which produced a 246 KB Casadi file.

import openap

import pandas as pd

from scipy.ndimage import gaussian_filter

from opentop.tools import cached_interpolant_from_dataframe

import xarray as xr

ds = xr.open_dataset('contrail_original.nc')\

.sel(time='2015-12-18')

level_pressure = [

0.0000,

10.0000,

30.0000,

50.0000,

70.0000,

90.0787,

110.6606,

132.3968,

155.7909,

181.1544,

208.6494,

238.3258,

270.1530,

304.0465,

339.8891,

377.5467,

416.8789,

457.7442,

500.0000,

543.4970,

588.0685,

633.5144,

679.5799,

725.9285,

772.1102,

817.5241,

861.3757,

902.6287,

939.9520,

971.6610,

995.6532,

1009.3396]

df = (

ds.to_dataframe()

.reset_index()

.assign(lev=lambda x: x.lev.astype(int))

.merge(

pd.DataFrame(level_pressure, columns=["hPa"]).reset_index(names="lev"),

on="lev",

)

.assign(height=lambda x: openap.aero.h_isa(x.hPa * 100).round(-2))

.assign(longitude=lambda x: ((x.lon + 180) % 360 - 180))

.query("height<15000"))

df_cost_world = df.rename(

columns={

"lat": "latitude",

"atr20_contrail": "cost",

}

)[["time",

"latitude",

"longitude",

"hPa",

"height",

"cost"]]

df_cost = df_cost_world.query(

"-20<longitude<40 and 30<latitude<70 and time.dt.hour==12"

).sort_values(["height", "latitude", "longitude"])

cost = df_cost.cost.values.reshape(

df_cost.height.nunique(),

df_cost.latitude.nunique(),

df_cost.longitude.nunique())

cost_ = gaussian_filter(cost, sigma=1, mode="nearest")

df_cost = df_cost.assign(cost=cost_.flatten())

interpolant = cached_interpolant_from_dataframe(

df_cost,

"contrail.casadi",

shape="bspline")

These are the first and last few bytes of its contents.

hexdump -C contrail.casadi | head

hexdump -C contrail.casadi | tail

I noticed the contents are repetitive and compress well.

gzip -9 < contrail.casadi | wc -c

I'll get the metrics of an optimal flight between Amsterdam's Schiphol (EHAM / AMS) and Frankfurt (EDDF / FRA) on an Airbus A320.

opentop optimize \

EHAM EDDF \

-a A320 \

--phase all \

--obj "0.3*fuel+0.7*grid" \

--grid contrail.casadi

I'll then do the same using a Boeing 737.

opentop optimize \

EHAM EDDF \

-a B737 \

--phase all \

--obj "0.3*fuel+0.7*grid" \

--grid contrail.casadi