diff --git a/docs/Project.toml b/docs/Project.toml
index c2513dd1..1871f85f 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -3,21 +3,27 @@ CTModels = "34c4fa32-2049-4079-8329-de33c2a22e2d"
DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
DocumenterInterLinks = "d12716ef-a0f6-4df4-a9f1-a5a34e75c656"
+DocumenterMermaid = "a078cd44-4d9c-4618-b545-3ab9d77f9177"
ExaModels = "1037b233-b668-4ce9-9b63-f9f681f55dd2"
+FilePathsBase = "48062228-2e41-5def-b9a4-89aafe57970f"
Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
+JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
Markdown = "d6f4376e-aef5-505a-96c1-9c027394607a"
NLPModels = "a4795742-8479-5a88-8948-cc11e1c8c1a6"
NLPModelsIpopt = "f4238b75-b362-5c4c-b852-0801c9a21d71"
OptimalControl = "5f98b655-cc9a-415a-b60e-744165666948"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
+Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
TOML = "fa267f1f-6049-4f14-aa54-33bafae1ed76"
+Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
[compat]
CTModels = "0.6"
DataFrames = "1"
Documenter = "1"
DocumenterInterLinks = "1"
+DocumenterMermaid = "0.2"
ExaModels = "0.9"
Ipopt = "1"
JuMP = "1"
@@ -26,4 +32,5 @@ NLPModels = "0.21"
NLPModelsIpopt = "0.10"
OptimalControl = "1"
Plots = "1"
+Printf = "1"
TOML = "1"
diff --git a/docs/browser.jl b/docs/browser.jl
new file mode 100644
index 00000000..567f8378
--- /dev/null
+++ b/docs/browser.jl
@@ -0,0 +1,831 @@
+using JSON
+using DataFrames
+using Tables
+using OptimalControlProblems
+using NLPModels
+using CTModels
+using FilePathsBase # optional, for robust path handling
+using Base.Filesystem: rm, mktemp
+
+const BROWSER_FILE = "problems_browser.md"
+const BROWSER_PATH = joinpath(@__DIR__, "src", BROWSER_FILE)
+
+# -------------------------------
+# MD + CSS + JS constants
+# -------------------------------
+const TABLE_PRESENTATION = """
+# [Problems browser](@id problems-browser)
+
+The table below provides an overview of all **optimal control problems** and allows interactive exploration, filtering, and export.
+
+!!! tip "Quick guide to the problems table"
+
+ ```@raw html
+ Click to unfold and see the quick guide for the table.
+ ```
+
+ ## Table Overview
+
+ - **Problem:** The name of the optimal control problem.
+ - **State:** Number of state variables in the system.
+ - **Control:** Number of control inputs.
+ - **Variable:** Number of additional optimisation variables (if any).
+ - **Cost:** Type of cost functional:
+ - **Mayer:** a **point cost** depending on the initial and final states and optional variables.
+ - **Lagrange:** integral over time
+ - **Bolza:** combination of Mayer + Lagrange
+ - **FinalTime:** Whether the final time is **fixed** or **free**.
+ - **Constraints:** Buttons representing each constraint type:
+ - **x:** state box constraints
+ - **u:** control box constraints
+ - **v:** variable box constraints
+ - **c:** nonlinear path constraints
+ - **b:** nonlinear boundary constraints
+
+ The number next to the buttons shows the **total number of constraints**. Hover over buttons to see the exact count. Click a row to see a detailed list of constraints.
+
+ ## Interactivity & Filters
+
+ - **Sorting & Search:** Click column headers to sort. Use the search box to filter by text.
+ - **Numeric Filters:** Enter `min-max` ranges in numeric columns to filter values.
+ - **Cost & FinalTime Filters:** Use dropdown menus to filter by cost type or whether the final time is fixed/free.
+ - **Constraint Filtering:** Use the buttons above the Constraints column to filter problems by constraint type. Choose **AND/OR logic** to combine multiple constraint types.
+ - **Export Buttons:** Use the top buttons to copy the table or export it to CSV, Excel, PDF, or print. Hover over each button to see its function.
+
+ ```@raw html
+
+ ```
+
+---
+
+Scroll through the table or use filters to quickly find problems of interest, inspect their constraints, and export data for further analysis.
+"""
+
+const TABLE_STYLE = """
+
+"""
+
+const TABLE_LOGIC = """
+
+"""
+
+const TABLE = (
+ presentation = TABLE_PRESENTATION,
+ style = TABLE_STYLE,
+ logic = TABLE_LOGIC,
+)
+
+# -------------------------------
+# Helpers
+# -------------------------------
+function sum_namedtuple(nt::NamedTuple)
+ sum(values(nt))
+end
+
+function ocp_data_to_json(data_ocp::DataFrame)
+ return JSON.json(Tables.columntable(data_ocp))
+end
+
+function write_block(io, content)
+ write(io, content)
+end
+
+function generate_constraint_buttons_html(constraint_dims::NamedTuple)
+ return constraint_dims
+end
+
+function build_table_html(json_str::String)
+ return """
+
+
+
+
+
Problem
+
State
+
Control
+
Variable
+
Cost
+
FinalTime
+
Constraints
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+"""
+end
+
+# -------------------------------
+# Data collection
+# -------------------------------
+function collect_problem_data(problem_sym::Symbol)
+ ocp = ocp_model(eval(problem_sym)(OptimalControlBackend()))
+
+ cost = has_mayer_cost(ocp) && has_lagrange_cost(ocp) ? "Bolza" :
+ has_mayer_cost(ocp) ? "Mayer" : "Lagrange"
+
+ final_time = has_fixed_final_time(ocp) ? "fixed" : "free"
+
+ dims = (
+ x = CTModels.dim_state_constraints_box(ocp),
+ u = CTModels.dim_control_constraints_box(ocp),
+ v = CTModels.dim_variable_constraints_box(ocp),
+ c = CTModels.dim_path_constraints_nl(ocp),
+ b = CTModels.dim_boundary_constraints_nl(ocp),
+ )
+
+ total = sum_namedtuple(dims)
+
+ return (
+ Problem = string(problem_sym),
+ State = state_dimension(ocp),
+ Control = control_dimension(ocp),
+ Variable = variable_dimension(ocp),
+ Cost = cost,
+ FinalTime = final_time,
+ DimStateConstraint = dims.x,
+ DimControlConstraint = dims.u,
+ DimVariableConstraint = dims.v,
+ DimPathConstraint = dims.c,
+ DimBoundaryConstraint = dims.b,
+ TotalConstraints = total,
+ ConstraintButtonsHtml = generate_constraint_buttons_html(dims)
+ )
+end
+
+"Collects OCP data and returns a DataFrame"
+function collect_problems_data()
+ return DataFrame([collect_problem_data(sym) for sym in problems()])
+end
+
+### MAIN GENERATION PROCESS
+
+# ---------------------------
+# Layer 1: Data collection
+# ---------------------------
+"Convert OCP data into JSON string"
+function collect_problems_data_json()
+ df = collect_problems_data()
+ return ocp_data_to_json(df)
+end
+
+# ---------------------------
+# Layer 2: HTML assembly
+# ---------------------------
+function assemble_problems_browser_html(json_str)
+ html_parts = [
+ TABLE_PRESENTATION,
+ "```@raw html",
+ TABLE_STYLE,
+ build_table_html(json_str),
+ TABLE_LOGIC,
+ "```",
+ ]
+ return join(html_parts, "\n")
+end
+
+# ---------------------------
+# Layer 3: File writing
+# ---------------------------
+function generate_problems_browser!(path::AbstractString)
+ mkpath(dirname(path))
+ json_str = collect_problems_data_json()
+ html = assemble_problems_browser_html(json_str)
+ open(path, "w") do io
+ write(io, html)
+ end
+ return path
+end
+
+# ---------------------------
+# Temporary browser context using fixed path
+# ---------------------------
+"""
+with_problems_browser(f::Function)
+
+Generates the problems browser at the fixed path, passes the filename
+to `f`, and removes the file after `f` finishes.
+"""
+function with_problems_browser(f::Function)
+ # Generate the problems browser at the fixed path
+ generate_problems_browser!(BROWSER_PATH)
+
+ try
+ # Pass the generated file path to the user function
+ return f(BROWSER_FILE)
+ finally
+ # Remove the file after usage
+ if isfile(BROWSER_PATH)
+ rm(BROWSER_PATH)
+ println("Temporary problems browser file removed: $BROWSER_PATH")
+ end
+ end
+end
diff --git a/docs/make.jl b/docs/make.jl
index e1ad9714..1a58aff1 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -1,25 +1,36 @@
-using Documenter
-using DocumenterInterLinks
+# ==============================
+# --- Import Packages & Modules ---
+# ==============================
+using Documenter # Core package to generate documentation
+using DocumenterInterLinks # For linking to other package docs
+using DocumenterMermaid # For Mermaid diagrams in docs
using OptimalControlProblems
using OptimalControl
using JuMP
using CTModels
using ExaModels
-include("problems.jl")
+include("problems.jl") # Your problem definitions
+include("browser.jl") # Your interactive problems browser
-# to add docstrings from external packages
+# ==============================
+# --- Add Docstrings from external packages ---
+# ==============================
+# Access extensions within OptimalControlProblems
const JuMPModels = Base.get_extension(OptimalControlProblems, :JuMPModels)
const OptimalControlModels = Base.get_extension(
OptimalControlProblems, :OptimalControlModels
)
+# Add DocTest setup for each module if missing
Modules = [JuMPModels, OptimalControlModels]
for Module in Modules
isnothing(DocMeta.getdocmeta(Module, :DocTestSetup)) &&
DocMeta.setdocmeta!(Module, :DocTestSetup, :(using $Module); recursive=true)
end
-#
+# ==============================
+# --- Configure inter-package links ---
+# ==============================
links = InterLinks(
"ADNLPModels" => (
"https://jso.dev/ADNLPModels.jl/stable/",
@@ -58,8 +69,11 @@ links = InterLinks(
),
)
-# For reproducibility
-mkpath(joinpath(@__DIR__, "src", "assets"))
+# ==============================
+# --- Copy Project Assets for reproducibility ---
+# ==============================
+mkpath(joinpath(@__DIR__, "src", "assets")) # Ensure assets folder exists
+
cp(
joinpath(@__DIR__, "Manifest.toml"),
joinpath(@__DIR__, "src", "assets", "Manifest.toml");
@@ -71,49 +85,64 @@ cp(
force=true,
)
+# Repository URL (used for links in docs)
repo_url = "github.com/control-toolbox/OptimalControlProblems.jl"
-#
-draft = true
+# ==============================
+# --- Generate Problems Documentation ---
+# ==============================
+draft = false # If true, code blocks in markdown are not executed
exclude_from_draft=Symbol[
-# :beam
+# :beam # example: exclude beam from draft docs
]
+
+# Generate problem pages (returns list of markdown pages)
PROBLEMS_PAGES = generate_documentation_problems(;
draft=draft, exclude_from_draft=exclude_from_draft
)
-# If draft is true below, then the julia code from .md is not executed.
-# To disable the draft mode in a specific markdown file, use the following:
+# ==============================
+# --- Generate Documentation ---
+# ==============================
+# If draft is true, the Julia code in markdown is not executed.
+# To disable draft mode in a specific markdown file, add:
#=
```@meta
Draft = false
```
=#
-makedocs(;
- draft=draft,
- #remotes=nothing,
- warnonly=:cross_references,
- sitename="OptimalControlProblems.jl",
- format=Documenter.HTML(;
- repolink="https://" * repo_url,
- prettyurls=false,
- size_threshold_ignore=["dev-api.md", PROBLEMS_PAGES...],
- assets=[
- asset("https://control-toolbox.org/assets/css/documentation.css"),
- asset("https://control-toolbox.org/assets/js/documentation.js"),
- ],
- ),
- pages=[
- "Getting Started" => "index.md",
- "Problems" =>
- ["problems-introduction.md", "List of the problems" => PROBLEMS_PAGES],
- "Tutorials" => [
- "Get a problem" => "tutorial-get.md",
- "Solve a problem" => "tutorial-solve.md",
+with_problems_browser() do browser_file # generates the problems browser and remove it at the end
+
+ makedocs(;
+ draft=draft,
+ #remotes=nothing,
+ warnonly=:cross_references,
+ sitename="OptimalControlProblems.jl",
+ format=Documenter.HTML(;
+ repolink="https://" * repo_url,
+ prettyurls=false,
+ size_threshold_ignore=["dev-api.md", PROBLEMS_PAGES...],
+ assets=[
+ asset("https://control-toolbox.org/assets/css/documentation.css"),
+ asset("https://control-toolbox.org/assets/js/documentation.js"),
+ ],
+ ),
+ pages=[
+ "Getting Started" => "index.md",
+ "Problems" =>
+ ["problems-introduction.md", browser_file, "List of the problems" => PROBLEMS_PAGES],
+ "Tutorials" => [
+ "Get a problem" => "tutorial-get.md",
+ "Solve a problem" => "tutorial-solve.md",
+ ],
+ "Developers" => ["Add a problem" => "dev-add.md", "API" => "dev-api.md"],
],
- "Developers" => ["Add a problem" => "dev-add.md", "API" => "dev-api.md"],
- ],
- plugins=[links],
-)
+ plugins=[links],
+ )
+
+end
+# ==============================
+# --- Deploy Documentation ---
+# ==============================
deploydocs(; repo=repo_url * ".git", devbranch="main", push_preview=true)
diff --git a/docs/problems.jl b/docs/problems.jl
index 476fea2c..119229a1 100644
--- a/docs/problems.jl
+++ b/docs/problems.jl
@@ -1,32 +1,43 @@
-function generate_documentation(
- PROBLEM::String, DESCRIPTION::String; draft::Union{Bool,Nothing}
-)
- TITLE = uppercasefirst(replace(PROBLEM, "_" => " "))
-
- DRAFT = if isnothing(draft)
- ""
+# -----------------------------------
+# Helper: draft metadata block
+# -----------------------------------
+function draft_meta(draft::Union{Bool,Nothing})
+ if isnothing(draft)
+ return ""
elseif draft
- """
- ```@meta
- Draft = true
- ```
- """
+ return """```@meta\nDraft = true\n```"""
else
- """
- ```@meta
- Draft = false
- ```
- """
+ return """```@meta\nDraft = false\n```"""
end
+end
+
+# -----------------------------------
+# Helper: left margin for plots
+# -----------------------------------
+function get_left_margin(problem::Symbol)
+ margins = Dict(:beam => "5mm")
+ return get(margins, problem, "20mm")
+end
+
+# -----------------------------------
+# Generate documentation for a problem
+# -----------------------------------
+function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Union{Bool,Nothing})
+
+ TITLE = uppercasefirst(replace(PROBLEM, "_" => " "))
+ DRAFT = draft_meta(draft)
+ LEFT_MARGIN = get_left_margin(Symbol(PROBLEM))
documentation=DRAFT * """
# $TITLE
+ ## Description of the problem
+
$DESCRIPTION
- ## Packages
+ ## Numerical set-up
- Import all necessary packages and define DataFrames to store information about the problem and resolution results.
+ In this section, we prepare the numerical environment required to study the problem. We begin by importing the relevant Julia packages and then initialise the data frames that will store the results of our simulations and computations. These structures provide the foundation for solving the problem and for comparing the different solution strategies in a consistent way.
```@example main
using OptimalControlProblems # to access the Beam model
@@ -38,6 +49,7 @@ function generate_documentation(
using Plots.PlotMeasures # for leftmargin, bottommargin
using JuMP # to import the JuMP model
using Ipopt # to solve the JuMP model with Ipopt
+ using Printf # to print
data_pb = DataFrame( # to store data about the problem
Problem=Symbol[],
@@ -55,95 +67,132 @@ function generate_documentation(
nothing # hide
```
- ## initial guess
+ ## Initial guess
- The initial guess (or first iterate) can be visualised by running the solver with `max_iter=0`. Here is the initial guess.
+ Before solving the problem, it is often useful to inspect the initial guess (sometimes called the first iterate). This guess is obtained by running the NLP solver with `max_iter = 0`, which evaluates the problem formulation without performing any optimisation steps.
- ```@raw html
- Click to unfold and see the code for plotting the initial guess.
- ```
+ We plot the resulting trajectories for both the OptimalControl and JuMP models. Since both backends represent the same mathematical problem, their initial guesses should coincide, providing a useful consistency check before moving on to the optimised solution.
- ```@example main
- function plot_initial_guess(problem)
-
- # dimensions
- x_vars = metadata[problem][:state_name]
- u_vars = metadata[problem][:control_name]
- n = length(x_vars) # number of states
- m = length(u_vars) # number of controls
-
- # import OptimalControl model
- docp = eval(problem)(OptimalControlBackend())
- nlp_oc = nlp_model(docp)
-
- # solve
- nlp_oc_sol = NLPModelsIpopt.ipopt(nlp_oc; max_iter=0)
-
- # build an optimal control solution
- ocp_sol = build_ocp_solution(docp, nlp_oc_sol)
-
- # plot the OptimalControl solution
- plt = plot(
- ocp_sol;
- state_style=(color=1,),
- costate_style=(color=1, legend=:none),
- control_style=(color=1, legend=:none),
- path_style=(color=1, legend=:none),
- dual_style=(color=1, legend=:none),
- size=(816, 220*(n+m)),
- label="OptimalControl",
- leftmargin=20mm,
- )
- for i in 2:n
- plot!(plt[i]; legend=:none)
- end
+ !!! note "Code to plot the initial guess"
- # import JuMP model
- nlp_jp = eval(problem)(JuMPBackend())
+ ```@raw html
+ Click to unfold and see the code for plotting the initial guess.
+ ```
- # solve
- set_optimizer(nlp_jp, Ipopt.Optimizer)
- set_optimizer_attribute(nlp_jp, "max_iter", 0)
- optimize!(nlp_jp)
+ ```@example main
+ function plot_initial_guess(problem)
+
+ # -----------------------------
+ # Extract dimensions from metadata
+ # -----------------------------
+ x_vars = metadata[problem][:state_name]
+ u_vars = metadata[problem][:control_name]
+ n_states = length(x_vars)
+ n_controls = length(u_vars)
+
+ # -----------------------------
+ # Build OptimalControl problem
+ # -----------------------------
+ ocp_model = eval(problem)(OptimalControlBackend())
+ nlp_oc = nlp_model(ocp_model)
+
+ # Solve NLP with zero iterations (initial guess)
+ nlp_oc_sol = NLPModelsIpopt.ipopt(nlp_oc; max_iter=0)
+
+ # Build OptimalControl solution
+ ocp_sol = build_ocp_solution(ocp_model, nlp_oc_sol)
+
+ # -----------------------------
+ # Plot OptimalControl solution
+ # -----------------------------
+ plt = plot(
+ ocp_sol;
+ state_style=(color=1,),
+ costate_style=(color=1, legend=:none),
+ control_style=(color=1, legend=:none),
+ path_style=(color=1, legend=:none),
+ dual_style=(color=1, legend=:none),
+ size=(816, 220*(n_states+n_controls)),
+ label="OptimalControl",
+ leftmargin=$LEFT_MARGIN,
+ )
+
+ # Hide legend for additional state plots
+ for i in 2:n_states
+ plot!(plt[i]; legend=:none)
+ end
- # plot
- t = time_grid(problem, nlp_jp) # t0, ..., tN = tf
- x = state(problem, nlp_jp) # function of time
- u = control(problem, nlp_jp) # function of time
- p = costate(problem, nlp_jp) # function of time
+ # -----------------------------
+ # Build JuMP model
+ # -----------------------------
+ nlp_jp = eval(problem)(JuMPBackend())
+
+ # Solve NLP with zero iterations (initial guess)
+ set_optimizer(nlp_jp, Ipopt.Optimizer)
+ set_optimizer_attribute(nlp_jp, "max_iter", 0)
+ optimize!(nlp_jp)
+
+ # Extract trajectories
+ t_grid = time_grid(problem, nlp_jp)
+ x_fun = state(problem, nlp_jp)
+ u_fun = control(problem, nlp_jp)
+ p_fun = costate(problem, nlp_jp)
+
+ # -----------------------------
+ # Plot JuMP solution on top
+ # -----------------------------
+ # States
+ for i in 1:n_states
+ label = i == 1 ? "JuMP" : :none
+ plot!(plt[i], t_grid, t -> x_fun(t)[i]; color=2, linestyle=:dash, label=label)
+ end
- for i in 1:n # state
- label = i == 1 ? "JuMP" : :none
- plot!(plt[i], t, t -> x(t)[i]; color=2, linestyle=:dash, label=label)
- end
+ # Costates
+ for i in 1:n_states
+ plot!(plt[n_states+i], t_grid, t -> -p_fun(t)[i]; color=2, linestyle=:dash, label=:none)
+ end
- for i in 1:n # costate
- plot!(plt[n+i], t, t -> -p(t)[i]; color=2, linestyle=:dash, label=:none)
- end
+ # Controls
+ for i in 1:n_controls
+ plot!(plt[2*n_states+i], t_grid, t -> u_fun(t)[i]; color=2, linestyle=:dash, label=:none)
+ end
- for i in 1:m # control
- plot!(plt[2n+i], t, t -> u(t)[i]; color=2, linestyle=:dash, label=:none)
+ return plt
end
+ nothing # hide
+ ```
- return plt
- end
- nothing # hide
- ```
-
- ```@raw html
-
-
- ```
+ ```@raw html
+
+ ```
```@example main
plot_initial_guess(:$PROBLEM)
```
- ## Solve the problem
+ ## Solving the problem
+
+ To solve an optimal control problem, we can rely on two complementary formulations: the `OptimalControl` backend, which works directly with the discretised control problem, and the `JuMP` backend, which leverages JuMP’s flexible modelling framework.
+
+ Both approaches generate equivalent NLPs that can be solved with Ipopt, and comparing them ensures consistency between the two formulations.
+
+ Before solving, we can inspect the discretisation details of the problem. The table below reports the number of grid points, decision variables, and constraints associated with the chosen formulation.
+
+ ```@example main
+ push!(data_pb,
+ (
+ Problem=:$PROBLEM,
+ Grid_Size=metadata[:$PROBLEM][:N],
+ Variables=get_nvar(nlp_model($PROBLEM(OptimalControlBackend()))),
+ Constraints=get_ncon(nlp_model($PROBLEM(OptimalControlBackend()))),
+ )
+ )
+ data_pb # hide
+ ```
### OptimalControl model
- Import the OptimalControl model and solve it.
+ We first solve the problem using the `OptimalControl` backend. The process begins by importing the problem definition and constructing the associated nonlinear programming (NLP) model. This NLP is then passed to the Ipopt solver, with standard options for tolerance and barrier parameter strategy.
```@example main
# import DOCP model
@@ -163,29 +212,15 @@ function generate_documentation(
nothing # hide
```
- The problem has the following numbers of steps, variables and constraints.
-
- ```@example main
- push!(data_pb,
- (
- Problem=:$PROBLEM,
- Grid_Size=metadata[:$PROBLEM][:N],
- Variables=get_nvar(nlp_oc),
- Constraints=get_ncon(nlp_oc),
- )
- )
- data_pb # hide
- ```
-
### JuMP model
- Import the JuMP model and solve it.
+ We now repeat the procedure using the `JuMP` backend. Here, the problem is reformulated as a JuMP model, which offers a flexible and widely used framework for nonlinear optimisation in Julia. The solver settings are chosen to mirror those used previously, so that the results can be compared on an equal footing.
```@example main
# import model
nlp_jp = $PROBLEM(JuMPBackend())
- # solve
+ # solve with Ipopt
set_optimizer(nlp_jp, Ipopt.Optimizer)
set_optimizer_attribute(nlp_jp, "print_level", 4)
set_optimizer_attribute(nlp_jp, "tol", 1e-8)
@@ -197,7 +232,7 @@ function generate_documentation(
## Numerical comparisons
- Let's get the flag, the number of iterations and the objective value from the resolutions.
+ In this section, we examine the results of the problem resolutions. We extract the solver status (flag), the number of iterations, and the objective value for each model. This provides a first overview of how each approach performs and sets the stage for a more detailed comparison of the solution trajectories.
```@example main
# from OptimalControl model
@@ -222,139 +257,108 @@ function generate_documentation(
data_re # hide
```
- We compare the OptimalControl and JuMP solutions in terms of the number of iterations, the \$L^2\$-norm of the differences in the state, control, and variable, as well as the objective values. Both absolute and relative errors are reported.
+ We compare the solutions obtained from the OptimalControl and JuMP models by examining the number of iterations required for convergence, the \$L^2\$-norms of the differences in states, controls, and additional variables, and the corresponding objective values. Both absolute and relative errors are reported, providing a clear quantitative measure of the agreement between the two approaches.
- ```@raw html
- Click to unfold and get the code of the numerical comparison.
- ```
+ !!! note "Code to print the numerical comparisons"
- ```@example main
- function L2_norm(T, X)
- # T and X are supposed to be one dimensional
- s = 0.0
- for i in 1:(length(T) - 1)
- s += 0.5 * (X[i]^2 + X[i + 1]^2) * (T[i + 1]-T[i])
- end
- return √(s)
- end
+ ```@raw html
+ Click to unfold and get the code of the numerical comparisons.
+ ```
- function numerical_comparison(problem, docp, nlp_oc_sol, nlp_jp)
-
- # get relevant data from OptimalControl model
- ocp_sol = build_ocp_solution(docp, nlp_oc_sol) # build an ocp solution
- t_oc = time_grid(ocp_sol)
- x_oc = state(ocp_sol).(t_oc)
- u_oc = control(ocp_sol).(t_oc)
- v_oc = variable(ocp_sol)
- o_oc = objective(ocp_sol)
- i_oc = iterations(ocp_sol)
-
- # get relevant data from JuMP model
- t_jp = time_grid(problem, nlp_jp)
- x_jp = state(problem, nlp_jp).(t_jp)
- u_jp = control(problem, nlp_jp).(t_jp)
- o_jp = objective(problem, nlp_jp)
- v_jp = variable(problem, nlp_jp)
- i_jp = iterations(problem, nlp_jp)
-
- x_vars = metadata[problem][:state_name]
- u_vars = metadata[problem][:control_name]
- v_vars = metadata[problem][:variable_name]
-
- println("┌─ ", string(problem))
- println("│")
-
- # number of iterations
- println("├─ Number of iterations")
- println("│")
- println("│ OptimalControl : ", i_oc)
- println("│ JuMP : ", i_jp)
- println("│")
-
- # state
- for i in eachindex(x_vars)
- xi_oc = [x_oc[k][i] for k in eachindex(t_oc)]
- xi_jp = [x_jp[k][i] for k in eachindex(t_jp)]
- L2_oc = L2_norm(t_oc, xi_oc)
- L2_jp = L2_norm(t_oc, xi_jp)
- L2_ae = L2_norm(t_oc, xi_oc-xi_jp)
- L2_re = L2_ae/(0.5*(L2_oc + L2_jp))
-
- println("├─ State \$(x_vars[i]) (L2 norm)")
- println("│")
- #println("│ OptimalControl : ", L2_oc)
- #println("│ JuMP : ", L2_jp)
- println("│ Absolute error : ", L2_ae)
- println("│ Relative error : ", L2_re)
- println("│")
+ ```@example main
+ function L2_norm(T, X)
+ # T and X are supposed to be one dimensional
+ s = 0.0
+ for i in 1:(length(T) - 1)
+ s += 0.5 * (X[i]^2 + X[i + 1]^2) * (T[i + 1]-T[i])
+ end
+ return √(s)
end
-
- # control
- for i in eachindex(u_vars)
- ui_oc = [u_oc[k][i] for k in eachindex(t_oc)]
- ui_jp = [u_jp[k][i] for k in eachindex(t_jp)]
- L2_oc = L2_norm(t_oc, ui_oc)
- L2_jp = L2_norm(t_oc, ui_jp)
- L2_ae = L2_norm(t_oc, ui_oc-ui_jp)
- L2_re = L2_ae/(0.5*(L2_oc + L2_jp))
-
- println("├─ Control \$(u_vars[i]) (L2 norm)")
+
+ function print_numerical_comparisons(problem, docp, nlp_oc_sol, nlp_jp)
+
+ # get relevant data from OptimalControl model
+ ocp_sol = build_ocp_solution(docp, nlp_oc_sol)
+ t_oc = time_grid(ocp_sol)
+ x_oc = state(ocp_sol).(t_oc)
+ u_oc = control(ocp_sol).(t_oc)
+ v_oc = variable(ocp_sol)
+ o_oc = objective(ocp_sol)
+ i_oc = iterations(ocp_sol)
+
+ # get relevant data from JuMP model
+ t_jp = time_grid(problem, nlp_jp)
+ x_jp = state(problem, nlp_jp).(t_jp)
+ u_jp = control(problem, nlp_jp).(t_jp)
+ o_jp = objective(problem, nlp_jp)
+ v_jp = variable(problem, nlp_jp)
+ i_jp = iterations(problem, nlp_jp)
+
+ x_vars = metadata[problem][:state_name]
+ u_vars = metadata[problem][:control_name]
+ v_vars = metadata[problem][:variable_name]
+
+ println("┌─ ", string(problem))
println("│")
- #println("│ OptimalControl : ", L2_oc)
- #println("│ JuMP : ", L2_jp)
- println("│ Absolute error : ", L2_ae)
- println("│ Relative error : ", L2_re)
- println("│")
- end
-
- # variable
- if !isnothing(v_vars)
- for i in eachindex(v_vars)
- vi_oc = v_oc[i]
- vi_jp = v_jp[i]
- vi_ae = abs(vi_oc-vi_jp)
- vi_re = vi_ae/(0.5*(abs(vi_oc) + abs(vi_jp)))
-
- println("├─ Variable \$(v_vars[i])")
- println("│")
- #println("│ OptimalControl : ", vi_oc)
- #println("│ JuMP : ", vi_jp)
- println("│ Absolute error : ", vi_ae)
- println("│ Relative error : ", vi_re)
- println("│")
+ println("├─ Number of Iterations")
+ @printf("│ OptimalControl : %d JuMP : %d\\n", i_oc, i_jp)
+
+ # States
+ println("├─ States (L2 Norms)")
+ for i in eachindex(x_vars)
+ xi_oc = [x_oc[k][i] for k in eachindex(t_oc)]
+ xi_jp = [x_jp[k][i] for k in eachindex(t_jp)]
+ L2_ae = L2_norm(t_oc, xi_oc - xi_jp)
+ L2_re = L2_ae / (0.5 * (L2_norm(t_oc, xi_oc) + L2_norm(t_oc, xi_jp)))
+ @printf("│ %-6s Abs: %.3e Rel: %.3e\\n", x_vars[i], L2_ae, L2_re)
end
- end
- # objective
- o_ae = abs(o_oc-o_jp)
- o_re = o_ae/(0.5*(abs(o_oc) + abs(o_jp)))
+ # Controls
+ println("├─ Controls (L2 Norms)")
+ for i in eachindex(u_vars)
+ ui_oc = [u_oc[k][i] for k in eachindex(t_oc)]
+ ui_jp = [u_jp[k][i] for k in eachindex(t_jp)]
+ L2_ae = L2_norm(t_oc, ui_oc - ui_jp)
+ L2_re = L2_ae / (0.5 * (L2_norm(t_oc, ui_oc) + L2_norm(t_oc, ui_jp)))
+ @printf("│ %-6s Abs: %.3e Rel: %.3e\\n", u_vars[i], L2_ae, L2_re)
+ end
- println("├─ objective")
- println("│")
- #println("│ OptimalControl : ", o_oc)
- #println("│ JuMP : ", o_jp)
- println("│ Absolute error : ", o_ae)
- println("│ Relative error : ", o_re)
- println("│")
- println("└─")
+ # Variables
+ if !isnothing(v_vars)
+ println("├─ Variables")
+ for i in eachindex(v_vars)
+ vi_oc = v_oc[i]
+ vi_jp = v_jp[i]
+ vi_ae = abs(vi_oc - vi_jp)
+ vi_re = vi_ae / (0.5 * (abs(vi_oc) + abs(vi_jp)))
+ @printf("│ %-6s Abs: %.3e Rel: %.3e\\n", v_vars[i], vi_ae, vi_re)
+ end
+ end
- return nothing
- end
- nothing # hide
- ```
+ # Objective
+ o_ae = abs(o_oc - o_jp)
+ o_re = o_ae / (0.5 * (abs(o_oc) + abs(o_jp)))
+ println("├─ Objective")
+ @printf("│ Abs: %.3e Rel: %.3e\\n", o_ae, o_re)
+ println("└─")
+ return nothing
+ end
+ nothing # hide
+ ```
- ```@raw html
-
-
- ```
+ ```@raw html
+
+ ```
```@example main
- numerical_comparison(:$PROBLEM, docp, nlp_oc_sol, nlp_jp)
+ print_numerical_comparisons(:$PROBLEM, docp, nlp_oc_sol, nlp_jp)
```
- ## Plot the solutions
+ ## Plotting the solutions
- Visualise states, costates, and controls from the OptimalControl and JuMP solutions:
+ In this section, we visualise the trajectories of the states, costates, and controls obtained from both the OptimalControl and JuMP solutions. The plots provide an intuitive way to compare the two approaches and to observe how the constraints and the optimal control influence the system dynamics.
+
+ For each variable, the OptimalControl solution is shown in solid lines, while the JuMP solution is overlaid using dashed lines. Since both models represent the same mathematical problem, their trajectories should closely coincide, highlighting the consistency between the two formulations.
```@example main
# build an ocp solution to use the plot from OptimalControl package
@@ -367,16 +371,12 @@ function generate_documentation(
# from OptimalControl solution
plt = plot(
ocp_sol;
- state_style=(color=1,),
- costate_style=(color=1, legend=:none),
- control_style=(color=1, legend=:none),
- path_style=(color=1, legend=:none),
- dual_style=(color=1, legend=:none),
+ color=1,
size=(816, 240*(n+m)),
label="OptimalControl",
- leftmargin=20mm,
+ leftmargin=$LEFT_MARGIN,
)
- for i in 2:n
+ for i in 2:length(plt)
plot!(plt[i]; legend=:none)
end
@@ -405,45 +405,29 @@ function generate_documentation(
return documentation
end
-function generate_documentation_problems(;
- draft::Union{Bool,Nothing}=nothing, exclude_from_draft::Vector{Symbol}=Symbol[]
-)
+# -----------------------------------
+# Generate documentation for all problems
+# -----------------------------------
+function generate_documentation_problems(; draft::Union{Bool,Nothing}=nothing,
+ exclude_from_draft::Vector{Symbol}=Symbol[])
- # List of problems
problems_list = problems()
+ problems_pages = map(p -> joinpath("problems", string(p) * ".md"), problems_list)
- #
- problems_pages = []
- for problem in problems_list
- push!(problems_pages, joinpath("problems", string(problem) * ".md"))
- end
+ # reset problems directory
+ problems_dir = joinpath(@__DIR__, "src", "problems")
+ rm(problems_dir; recursive=true, force=true)
+ mkpath(problems_dir)
+ mkpath(joinpath(problems_dir, "assets"))
- # remove and create problems directory
- rm(joinpath(@__DIR__, "src", "problems"); recursive=true, force=true)
- mkpath(joinpath(@__DIR__, "src", "problems"))
- mkpath(joinpath(@__DIR__, "src", "problems", "assets"))
-
- # create file for documentation
for problem in problems_list
-
- # create the file
- filename = joinpath(@__DIR__, "src", "problems", string(problem) * ".md")
- touch(filename)
-
- # get the description
- description = read(
- joinpath(@__DIR__, "..", "ext", "Descriptions", string(problem) * ".md"), String
- )
-
- # generate the content
+ description = read(joinpath(@__DIR__, "..", "ext", "Descriptions", string(problem) * ".md"), String)
draft_problem = problem ∈ exclude_from_draft ? false : draft
contents = generate_documentation(string(problem), description; draft=draft_problem)
- # write the content in the file
- open(filename, "a") do io
- write(io, contents)
- end
+ filename = joinpath(problems_dir, string(problem) * ".md")
+ write(filename, contents)
end
return problems_pages
-end
+end
\ No newline at end of file
diff --git a/docs/src/assets/Manifest.toml b/docs/src/assets/Manifest.toml
index 8cd9666c..8d1d9468 100644
--- a/docs/src/assets/Manifest.toml
+++ b/docs/src/assets/Manifest.toml
@@ -2,7 +2,7 @@
julia_version = "1.11.6"
manifest_format = "2.0"
-project_hash = "e4420c03c5c208c963d869daa22ca1e9ba918652"
+project_hash = "05847f69ae3d937e6a44ff3a4a0761aed56ed09e"
[[deps.ADNLPModels]]
deps = ["ADTypes", "ForwardDiff", "LinearAlgebra", "NLPModels", "Requires", "ReverseDiff", "SparseArrays", "SparseConnectivityTracer", "SparseMatrixColorings"]
@@ -136,10 +136,10 @@ version = "0.8.8"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
[[deps.CTModels]]
-deps = ["CTBase", "DocStringExtensions", "Interpolations", "LinearAlgebra", "MLStyle", "MacroTools", "OrderedCollections", "Parameters", "PrettyTables", "RecipesBase"]
-git-tree-sha1 = "81becfc3e56dd8f746edaa023a55251b0a60747a"
+deps = ["CTBase", "DocStringExtensions", "Interpolations", "LinearAlgebra", "MLStyle", "MacroTools", "OrderedCollections", "Parameters", "RecipesBase"]
+git-tree-sha1 = "20872a1b453a9b7a94822cc00a22c6741f71cf68"
uuid = "34c4fa32-2049-4079-8329-de33c2a22e2d"
-version = "0.6.4"
+version = "0.6.5"
[deps.CTModels.extensions]
CTModelsJLD = "JLD2"
@@ -153,13 +153,9 @@ version = "0.6.4"
[[deps.CTParser]]
deps = ["CTBase", "DocStringExtensions", "MLStyle", "OrderedCollections", "Parameters", "Unicode"]
-git-tree-sha1 = "45b4887484fb4b39c089bc175cea01510e5552e5"
+git-tree-sha1 = "507a19db4be22ac5aa6967d78c97e68fd95a670a"
uuid = "32681960-a1b1-40db-9bff-a1ca817385d1"
-version = "0.6.3"
-
- [deps.CTParser.weakdeps]
- CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba"
- MadNLP = "2621e9c9-9eb4-46b1-8089-e8c72242dfb6"
+version = "0.6.4"
[[deps.Cairo_jll]]
deps = ["Artifacts", "Bzip2_jll", "CompilerSupportLibraries_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "JLLWrappers", "LZO_jll", "Libdl", "Pixman_jll", "Xorg_libXext_jll", "Xorg_libXrender_jll", "Zlib_jll", "libpng_jll"]
@@ -342,6 +338,12 @@ git-tree-sha1 = "d8a8cb2d5b0181fbbd41861016b221b0202c9bc5"
uuid = "d12716ef-a0f6-4df4-a9f1-a5a34e75c656"
version = "1.1.0"
+[[deps.DocumenterMermaid]]
+deps = ["Documenter", "MarkdownAST"]
+git-tree-sha1 = "c0c408289a505a15496d914b19204e272a7e8b0f"
+uuid = "a078cd44-4d9c-4618-b545-3ab9d77f9177"
+version = "0.2.0"
+
[[deps.Downloads]]
deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"]
uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"
@@ -418,6 +420,17 @@ git-tree-sha1 = "acebe244d53ee1b461970f8910c235b259e772ef"
uuid = "9aa1b823-49e4-5ca5-8b0f-3971ec8bab6a"
version = "0.3.2"
+[[deps.FilePathsBase]]
+deps = ["Compat", "Dates"]
+git-tree-sha1 = "3bab2c5aa25e7840a4b065805c0cdfc01f3068d2"
+uuid = "48062228-2e41-5def-b9a4-89aafe57970f"
+version = "0.9.24"
+weakdeps = ["Mmap", "Test"]
+
+ [deps.FilePathsBase.extensions]
+ FilePathsBaseMmapExt = "Mmap"
+ FilePathsBaseTestExt = "Test"
+
[[deps.FileWatching]]
uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee"
version = "1.11.0"
@@ -511,6 +524,12 @@ git-tree-sha1 = "45288942190db7c5f760f59c04495064eedf9340"
uuid = "b0724c58-0f36-5564-988d-3bb0596ebc4a"
version = "0.22.4+0"
+[[deps.Ghostscript_jll]]
+deps = ["Artifacts", "JLLWrappers", "JpegTurbo_jll", "Libdl", "Zlib_jll"]
+git-tree-sha1 = "38044a04637976140074d0b0621c1edf0eb531fd"
+uuid = "61579ee1-b43e-5ca0-a5da-69d92c66a64b"
+version = "9.55.1+0"
+
[[deps.Git]]
deps = ["Git_LFS_jll", "Git_jll", "JLLWrappers", "OpenSSH_jll"]
git-tree-sha1 = "824a1890086880696fc908fe12a17bcf61738bd8"
@@ -525,9 +544,9 @@ version = "3.7.0+0"
[[deps.Git_jll]]
deps = ["Artifacts", "Expat_jll", "JLLWrappers", "LibCURL_jll", "Libdl", "Libiconv_jll", "OpenSSL_jll", "PCRE2_jll", "Zlib_jll"]
-git-tree-sha1 = "cd06e503111a7c5ef1d4a339de6ccf5bd7437b32"
+git-tree-sha1 = "e2aef26f7d273f1e5b1daba56837c47b49b4388f"
uuid = "f8c6e375-362e-5223-8a59-34ff63f689eb"
-version = "2.51.0+0"
+version = "2.51.1+0"
[[deps.Glib_jll]]
deps = ["Artifacts", "GettextRuntime_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Libiconv_jll", "Libmount_jll", "PCRE2_jll", "Zlib_jll"]
@@ -571,10 +590,10 @@ uuid = "2e76f6c2-a576-52d4-95c1-20adfe4de566"
version = "8.5.1+0"
[[deps.Hwloc_jll]]
-deps = ["Artifacts", "JLLWrappers", "Libdl"]
-git-tree-sha1 = "92f65c4d78ce8cdbb6b68daf88889950b0a99d11"
+deps = ["Artifacts", "JLLWrappers", "Libdl", "XML2_jll", "Xorg_libpciaccess_jll"]
+git-tree-sha1 = "3d468106a05408f9f7b6f161d9e7715159af247b"
uuid = "e33a78d0-f292-5ffc-b300-72abe9b543c8"
-version = "2.12.1+0"
+version = "2.12.2+0"
[[deps.IOCapture]]
deps = ["Logging", "Random"]
@@ -624,9 +643,9 @@ version = "1.3.1"
[[deps.Ipopt]]
deps = ["Ipopt_jll", "LinearAlgebra", "OpenBLAS32_jll", "PrecompileTools"]
-git-tree-sha1 = "4ad0d2dea51e5d49866b40a2d2521da6a1be7097"
+git-tree-sha1 = "ef90a75a3ee8c2b170f6c177d4d003348dd30f67"
uuid = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
-version = "1.10.6"
+version = "1.11.0"
weakdeps = ["MathOptInterface"]
[deps.Ipopt.extensions]
@@ -634,9 +653,9 @@ weakdeps = ["MathOptInterface"]
[[deps.Ipopt_jll]]
deps = ["ASL_jll", "Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "MUMPS_seq_jll", "SPRAL_jll", "libblastrampoline_jll"]
-git-tree-sha1 = "1bb978524c2837be596aeb2b69951feb6b9822f8"
+git-tree-sha1 = "b33cbc78b8d4de87d18fcd705054a82e2999dbac"
uuid = "9cc047cb-c261-5740-88fc-0cf96f7bdcc7"
-version = "300.1400.1701+0"
+version = "300.1400.1900+0"
[[deps.IrrationalConstants]]
git-tree-sha1 = "e2222959fbc6c19554dc15174c81bf7bf3aa691c"
@@ -680,15 +699,15 @@ version = "1.14.3"
[[deps.JpegTurbo_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl"]
-git-tree-sha1 = "e95866623950267c1e4878846f848d94810de475"
+git-tree-sha1 = "4255f0032eafd6451d707a51d5f0248b8a165e4d"
uuid = "aacddb02-875f-59d6-b918-886e6ef4fbf8"
-version = "3.1.2+0"
+version = "3.1.3+0"
[[deps.JuMP]]
deps = ["LinearAlgebra", "MacroTools", "MathOptInterface", "MutableArithmetics", "OrderedCollections", "PrecompileTools", "Printf", "SparseArrays"]
-git-tree-sha1 = "d05a696a5abaf9d1f8bce948ee53ed1533fadfdb"
+git-tree-sha1 = "d9c29fadef257492791c83b34ceede0d92a51470"
uuid = "4076af6c-e467-56ae-b986-b466b2749572"
-version = "1.28.0"
+version = "1.29.0"
[deps.JuMP.extensions]
JuMPDimensionalDataExt = "DimensionalData"
@@ -726,10 +745,10 @@ uuid = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
version = "1.4.0"
[[deps.Latexify]]
-deps = ["Format", "InteractiveUtils", "LaTeXStrings", "MacroTools", "Markdown", "OrderedCollections", "Requires"]
-git-tree-sha1 = "52e1296ebbde0db845b356abbbe67fb82a0a116c"
+deps = ["Format", "Ghostscript_jll", "InteractiveUtils", "LaTeXStrings", "MacroTools", "Markdown", "OrderedCollections", "Requires"]
+git-tree-sha1 = "44f93c47f9cd6c7e431f2f2091fcba8f01cd7e8f"
uuid = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
-version = "0.16.9"
+version = "0.16.10"
[deps.Latexify.extensions]
DataFramesExt = "DataFrames"
@@ -896,9 +915,9 @@ version = "0.4.17"
[[deps.MUMPS_seq_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "METIS_jll", "libblastrampoline_jll"]
-git-tree-sha1 = "196f61d99adc06f32c32bc4afe5298d9b1e862c8"
+git-tree-sha1 = "fc0c8442887b48c15aec2b1787a5fc812a99b2fd"
uuid = "d7ed1dd3-d0ae-5e8e-bfb4-87a502085b8d"
-version = "500.800.0+0"
+version = "500.800.100+0"
[[deps.MacroTools]]
git-tree-sha1 = "1e0228a030642014fe5cfe68c2c0a818f9e3f522"
@@ -918,9 +937,9 @@ version = "0.1.2"
[[deps.MathOptInterface]]
deps = ["BenchmarkTools", "CodecBzip2", "CodecZlib", "DataStructures", "ForwardDiff", "JSON3", "LinearAlgebra", "MutableArithmetics", "NaNMath", "OrderedCollections", "PrecompileTools", "Printf", "SparseArrays", "SpecialFunctions", "Test"]
-git-tree-sha1 = "583abfbd38f15198adde0658e1b1222ece232ae2"
+git-tree-sha1 = "9603279ae328cb943a5f36ecd40de2774b5646d3"
uuid = "b8f27783-ece8-5eb3-8dc8-9495eed66fee"
-version = "1.43.0"
+version = "1.44.0"
[[deps.MbedTLS]]
deps = ["Dates", "MbedTLS_jll", "MozillaCACerts_jll", "NetworkOptions", "Random", "Sockets"]
@@ -1377,9 +1396,9 @@ version = "1.0.3"
[[deps.StaticArrays]]
deps = ["LinearAlgebra", "PrecompileTools", "Random", "StaticArraysCore"]
-git-tree-sha1 = "cbea8a6bd7bed51b1619658dec70035e07b8502f"
+git-tree-sha1 = "b8693004b385c842357406e3af647701fe783f98"
uuid = "90137ffa-7385-5640-81b9-e52037218182"
-version = "1.9.14"
+version = "1.9.15"
weakdeps = ["ChainRulesCore", "Statistics"]
[deps.StaticArrays.extensions]
@@ -1556,6 +1575,12 @@ git-tree-sha1 = "c1a7aa6219628fcd757dede0ca95e245c5cd9511"
uuid = "efce3f68-66dc-5838-9240-27a6d6f5f9b6"
version = "1.0.0"
+[[deps.XML2_jll]]
+deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Zlib_jll"]
+git-tree-sha1 = "59071150afa35787c1656ba234cf03fdf8e2603f"
+uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a"
+version = "2.13.8+0"
+
[[deps.XZ_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl"]
git-tree-sha1 = "fee71455b0aaa3440dfdd54a9a36ccef829be7d4"
@@ -1634,6 +1659,12 @@ git-tree-sha1 = "7ed9347888fac59a618302ee38216dd0379c480d"
uuid = "ea2f1a96-1ddc-540d-b46f-429655e07cfa"
version = "0.9.12+0"
+[[deps.Xorg_libpciaccess_jll]]
+deps = ["Artifacts", "JLLWrappers", "Libdl", "Zlib_jll"]
+git-tree-sha1 = "4909eb8f1cbf6bd4b1c30dd18b2ead9019ef2fad"
+uuid = "a65dc6b1-eb27-53a1-bb3e-dea574b5389e"
+version = "0.18.1+0"
+
[[deps.Xorg_libxcb_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Xorg_libXau_jll", "Xorg_libXdmcp_jll"]
git-tree-sha1 = "bfcaf7ec088eaba362093393fe11aa141fa15422"
diff --git a/docs/src/assets/Project.toml b/docs/src/assets/Project.toml
index 95e11132..3ca81632 100644
--- a/docs/src/assets/Project.toml
+++ b/docs/src/assets/Project.toml
@@ -3,8 +3,11 @@ CTModels = "34c4fa32-2049-4079-8329-de33c2a22e2d"
DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
DocumenterInterLinks = "d12716ef-a0f6-4df4-a9f1-a5a34e75c656"
+DocumenterMermaid = "a078cd44-4d9c-4618-b545-3ab9d77f9177"
ExaModels = "1037b233-b668-4ce9-9b63-f9f681f55dd2"
+FilePathsBase = "48062228-2e41-5def-b9a4-89aafe57970f"
Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
+JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
Markdown = "d6f4376e-aef5-505a-96c1-9c027394607a"
NLPModels = "a4795742-8479-5a88-8948-cc11e1c8c1a6"
@@ -12,13 +15,16 @@ NLPModelsIpopt = "f4238b75-b362-5c4c-b852-0801c9a21d71"
OptimalControl = "5f98b655-cc9a-415a-b60e-744165666948"
OptimalControlProblems = "59046045-fb9c-4c23-964f-ff0a25704f96"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
+Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
TOML = "fa267f1f-6049-4f14-aa54-33bafae1ed76"
+Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
[compat]
CTModels = "0.6"
DataFrames = "1"
Documenter = "1"
DocumenterInterLinks = "1"
+DocumenterMermaid = "0.2"
ExaModels = "0.9"
Ipopt = "1"
JuMP = "1"
@@ -27,4 +33,5 @@ NLPModels = "0.21"
NLPModelsIpopt = "0.10"
OptimalControl = "1"
Plots = "1"
+Printf = "1"
TOML = "1"
diff --git a/docs/src/problems-introduction.md b/docs/src/problems-introduction.md
index 81edcd65..e7c69347 100644
--- a/docs/src/problems-introduction.md
+++ b/docs/src/problems-introduction.md
@@ -94,7 +94,7 @@ Depth = 1
For each problem, additional data is provided in the [MetaData](https://github.com/control-toolbox/OptimalControlProblems.jl/tree/main/ext/MetaData) directory:
```@docs; canonical=false
-metadata
+OptimalControlProblems.metadata
```
To list all metadata, use `metadata`. To access the metadata of a specific problem, for example `chain`, run:
@@ -207,9 +207,12 @@ We detail below the characteristics of the optimal control problems and their as
```@raw html
-
```
+!!! tip "Interactive problems browser"
+
+ The table from [Problems browser](@ref problems-browser) page provides an overview of all **optimal control problems** and allows interactive exploration, filtering, and export.
+
For each optimal control problem, we provide the dimensions of the state, control, and variable. We also specify the type of objective function (Mayer, Lagrange, or Bolza), whether the final time is free or fixed, and whether there are constraints on the state (`x`), control (`u`), variable (`v`), path (`c`—from *chemin* in French, since `p` is reserved for the costate), or boundary (`b`).
```@example main
diff --git a/docs/src/tutorial-get.md b/docs/src/tutorial-get.md
index 2c5cfc4d..3b93e9b4 100644
--- a/docs/src/tutorial-get.md
+++ b/docs/src/tutorial-get.md
@@ -12,8 +12,6 @@ using OptimalControlProblems
## Get an OptimalControl model
-### discretise the optimal control problem models
-
To get an OptimalControl model, first [install](https://control-toolbox.org/OptimalControl.jl/stable/#Installation) OptimalControl and import the package:
```@example main_oc
diff --git a/ext/Descriptions/beam.md b/ext/Descriptions/beam.md
index fa0afde4..da8532e7 100644
--- a/ext/Descriptions/beam.md
+++ b/ext/Descriptions/beam.md
@@ -1,47 +1,68 @@
-The *beam problem* is a classical benchmark in constrained optimal control. It originates from the Euler–Bernoulli beam model (see Bryson et al. 1963). The system describes the deflection of a clamped flexible beam under an external force. Both the state trajectory $x(\cdot)$ and the control $u(\cdot)$ are decision variables. The aim is to minimise the control effort subject to boundary conditions and a state inequality constraint on the displacement.
+The **clamped beam problem** is a classical benchmark in constrained optimal control.
+It models the deflection of a flexible beam that is fixed (clamped) at both ends and subject to an external force.
+Originating from the Euler–Bernoulli beam model [Bryson et al. 1963](https://doi.org/10.2514/3.2107), it is widely used in benchmark suites.
+Both the state trajectory $x(\cdot)$ and the control $u(\cdot)$ are decision variables, and the aim is to minimise the control effort under boundary conditions and a state inequality constraint on the beam deflection.
-The problem can be written as
+### Mathematical formulation
-```math
-\min_{x,\,u} J(x,u) = \int_0^1 u(t)^2 \, dt
-```
+The problem can be stated as
-subject to the dynamics
```math
-\dot{x}_1(t) = x_2(t), \qquad
-\dot{x}_2(t) = u(t),
+\begin{aligned}
+\min_{x,u} \quad & J(x,u) = \int_0^1 u^2(t) \,\mathrm{d}t \\[1em]
+\text{s.t.} \quad & \dot{x}_1(t) = x_2(t), \quad \dot{x}_2(t) = u(t), \\[0.5em]
+& x_1(0) = 0, \; x_1(1) = 0, \; x_2(0) = 1, \; x_2(1) = -1, \\[0.5em]
+& 0 \le x_1(t) \le a.
+\end{aligned}
```
-with boundary conditions
+where $a$ denotes the maximal admissible deflection (e.g. $a=0.1$ in the current implementation).
+
+### Qualitative behaviour
+
+This problem involves a **second-order state constraint**.
+Differentiating the inequality twice reveals that the control is directly related to the curvature of the state:
+
```math
-x_1(0) = 0, \quad x_1(1) = 0, \qquad
-x_2(0) = 1, \quad x_2(1) = -1,
+\dot{x}_1(t) = x_2(t), \quad \ddot{x}_1(t) = u(t).
```
-and the constraints
+Along a **boundary arc**, where $x_1$ remains constant, both derivatives vanish.
+Consequently, the control along such an arc is zero:
+
```math
--10 \le u(t) \le 5, \qquad
-0 \le x_1(t) \le a,
+u(t) = 0.
```
-where $a$ denotes the maximal admissible deflection (e.g. $a=0.1$ in the current implementation).
+In optimal control terminology, a **touch point** refers to the case where the state hits the bound at a single instant,
+while a **boundary arc** designates an interval where the state remains on the bound.
-### Qualitative behaviour
+The qualitative behaviour of the optimal trajectory depends on the parameter $a$:
-- If $a \geq 1/4$, the inequality constraint is inactive and the solution is $x(t) = t(1-t)$.
-- If $a \in [1/6, 1/4]$, there is a single touch point at $t=1/2$.
-- If $a < 1/6$, a boundary arc appears, with the constraint active on an interval.
+- **Case 1: $a \ge 1/4$**.
+ The inequality constraint is **inactive**. The optimal solution is simply the unconstrained trajectory: $x(t) = t\, (1-t)$.
-These features illustrate the role of state inequality constraints in shaping optimal trajectories and controls.
+- **Case 2: $1/6 \le a \le 1/4$**.
+ The trajectory touches the boundary exactly once (a **touch point**) at $t=1/2$.
+ There is no interval of constraint activation, but the single contact reflects the influence of the state bound.
+
+- **Case 3: $a < 1/6$**.
+ A **boundary arc** appears: the state constraint is active over a finite interval.
+ Along this interval, the control vanishes: $u(t) = 0$.
+
+These cases illustrate how second-order state inequality constraints can yield unconstrained solutions, isolated touch points, or extended boundary arcs, thereby shaping both the trajectory and the control.
### Characteristics
- Linear–quadratic dynamics and cost.
-- Asymmetric control bounds.
-- State inequality constraint induces touch points or boundary arcs depending on $a$.
-- Widely used as a benchmark in numerical optimal control libraries.
+- State inequality constraint induces either unconstrained motion, touch points, or boundary arcs depending on the parameter $a$.
+- Serves as a benchmark for optimal control methods handling state constraints.
### References
-- Bryson, A.E., Denham, W.F., & Dreyfus, S.E. (1963). *Optimal programming problems with inequality constraints I: necessary conditions for extremal solutions*. AIAA Journal.
-- BOCOP examples: Clamped Beam problem. [Examples-BOCOP.pdf](https://project.inria.fr/bocop/files/2017/05/Examples-BOCOP.pdf)
+- Bryson, A.E., Denham, W.F., & Dreyfus, S.E. (1963). *Optimal programming problems with inequality constraints I: necessary conditions for extremal solutions*.
+ AIAA Journal. [doi.org/10.2514/3.2107](https://doi.org/10.2514/3.2107)
+ This seminal work develops the theoretical foundation for optimal control problems with inequality constraints, including the derivation of necessary conditions for extremal solutions. It provides the basis for analysing state-constrained control problems like the clamped beam.
+
+- BOCOP examples: Clamped Beam problem. [Examples-BOCOP.pdf](https://project.inria.fr/bocop/files/2017/05/Examples-BOCOP.pdf)
+ This example illustrates the practical implementation of the clamped beam problem in BOCOP, providing a benchmark for testing optimal control methods under state inequality constraints.
diff --git a/ext/Descriptions/bioreactor.md b/ext/Descriptions/bioreactor.md
index 92290e35..cf9adb95 100644
--- a/ext/Descriptions/bioreactor.md
+++ b/ext/Descriptions/bioreactor.md
@@ -1,51 +1,65 @@
-This problem models a coupled photobioreactor–digester system for methane production.
-The system consists of three state variables: the algae concentration $y(t)$, the substrate concentration $s(t)$, and the biomass concentration $b(t)$.
+The **photobioreactor–digester system problem** is a benchmark in constrained optimal control.
+It models the coupled dynamics of a microalgae photobioreactor and an anaerobic digester for methane production.
+The system includes three state variables: the algae concentration $y(t)$, the substrate concentration $s(t)$, and the biomass concentration $b(t)$.
The control variable $u(t)$ represents the input flow rate between the two units.
-The dynamics include algal growth driven by light, substrate consumption, and biomass evolution.
-The aim is to maximise methane production over a fixed time horizon under biological and operational constraints.
+The goal is to maximise methane production over a fixed horizon while satisfying biological and operational constraints [Bayen et al. 2014](https://doi.org/10.1002/oca.2127).
-### Problem formulation
+### Mathematical formulation
-We minimise
+The problem can be stated as
```math
-\min_{y,\,s,\,b,\,u} J(y,s,b,u) = - \frac{1}{\beta+c} \int_0^T \mu_2(s(t))\, b(t)\, dt,
+\begin{aligned}
+\min_{y,s,b,u} \quad & J(y,s,b,u) = - \frac{1}{\beta+c} \int_0^T \mu_2(s(t))\, b(t) \,\mathrm{d}t \\[1em]
+\text{s.t.} \quad &
+\dot{y}(t) = \frac{\mu(t)\, y(t)}{1+y(t)} - (r+u(t))\, y(t), \\[0.5em]
+& \dot{s}(t) = -\mu_2(s(t))\, b(t) + u(t)\, \beta\big(\gamma y(t) - s(t)\big), \\[0.5em]
+& \dot{b}(t) = \big(\mu_2(s(t)) - u(t)\, \beta\big)\, b(t), \\[0.5em]
+& 0 \le u(t) \le 1, \\[0.5em]
+& y(t) \ge 0,\; s(t) \ge 0,\; b(t) \ge 10^{-3}, \\[0.5em]
+& 0.05 \le y(0) \le 0.25,\;\; 0.5 \le s(0) \le 5,\;\; 0.5 \le b(0) \le 3.
+\end{aligned}
```
-subject to the dynamics
-```math
-\dot{y}(t) = \frac{\mu(t)\, y(t)}{1+y(t)} - (r+u(t))\, y(t),
-```
-```math
-\dot{s}(t) = -\mu_2(s(t))\, b(t) + u(t)\, \beta\big(\gamma y(t) - s(t)\big),
-```
-```math
-\dot{b}(t) = \big(\mu_2(s(t)) - u(t)\, \beta\big)\, b(t),
-```
+The horizon is fixed to $T = 200$ (rescaled units), corresponding to several day–night cycles.
-with
-- **Light model:** $\mu(t) = \mu_{\text{bar}} \,\max(0,\sin(\tau(t)))^2$, where $\tau(t)$ encodes a periodic day–night cycle,
-- **Growth function (Monod law):** $\mu_2(s) = \mu_2^m \,\dfrac{s}{K_s + s}$.
+The model components are:
-### Constraints
+- **Light model:** $\mu(t) = \mu_{\text{bar}} \,\max\!\big(0,\sin(\tau(t))\big)^2$, where $\tau(t)$ encodes the periodic day–night cycle,
+- **Growth function (Monod law):** $\mu_2(s) = \mu_2^m \,\dfrac{s}{K_s + s}$.
-- Control bounds: $0 \leq u(t) \leq 1$,
-- State bounds: $y(t) \geq 0,\; s(t) \geq 0,\; b(t) \geq 10^{-3}$,
-- Initial conditions: $0.05 \leq y(0) \leq 0.25$, $0.5 \leq s(0) \leq 5$, $0.5 \leq b(0) \leq 3$.
+### Parameter values
-The horizon is fixed to $T = 200$ (rescaled units), corresponding to several day–night periods.
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Flow coupling between reactors | $\beta$ | 1 |
+| Cost scaling | $c$ | 2 |
+| Substrate–algae interaction | $\gamma$ | 1 |
+| Half-period of light cycle | $\text{halfperiod}$ | 5 |
+| Monod half-saturation constant | $K_s$ | 0.05 |
+| Maximum biomass growth rate | $\mu_2^m$ | 0.1 |
+| Maximum light intensity | $\mu_{\text{bar}}$ | 1 |
+| Algal decay rate | $r$ | 0.005 |
+| Time horizon | $T$ | 200 |
### Qualitative behaviour
-The optimal solution exploits the periodic structure of the problem:
-during illuminated phases, algal growth is favoured, while in dark phases methane production becomes predominant.
-The control exhibits a near bang–bang structure, alternating between minimal and maximal input flows, with possible singular arcs when substrate and biomass reach balanced levels.
-The constraints on positivity of the states are typically active for $b(t)$ at the beginning of the process.
+The optimal solution exploits the **periodic light–dark structure**:
+algal growth is favoured during illuminated phases, while methane production dominates during dark phases.
+The control $u(t)$ typically exhibits a **bang–bang structure**, alternating between minimum and maximum values, with possible **singular arcs** when substrate and biomass reach balanced levels.
+The positivity constraint on biomass $b(t)$ is usually active at the beginning of the process, shaping the initial control action.
+
+### Characteristics
+
+- Nonlinear coupled dynamics with periodic forcing,
+- Control and state constraints ensuring feasibility,
+- Bang–bang and singular arc structures in the optimal control,
+- Serves as a benchmark for optimal control methods with periodic and constrained systems.
### References
-- Bayen, T., Mairet, F., Martinon, P., & Sebbah, M. (2014). *Analysis of a periodic optimal control problem connected to microalgae anaerobic digestion*. Optimal Control Applications and Methods. [DOI:10.1002/oca.2127](https://hal.archives-ouvertes.fr/hal-00860570)
-- Bayen, T., Mairet, F., Martinon, P., & Sebbah, M. (2013). *Optimising the anaerobic digestion of microalgae in a coupled process*. 13th European Control Conference.
-- Barbosa, M.J., & Wijffels, R.H. (2010). *An Outlook on Microalgal Biofuels*. Science, 329, 796–799.
-- Betts, J.T. (2001). *Practical methods for optimal control using nonlinear programming*. SIAM.
-- BOCOP repository: https://github.com/control-toolbox/bocop/tree/main/bocop
+- Bayen, T., Mairet, F., Martinon, P., & Sebbah, M. (2014). *Analysis of a periodic optimal control problem connected to microalgae anaerobic digestion*. Optimal Control Applications and Methods. [doi.org/10.1002/oca.2127](https://doi.org/10.1002/oca.2127)
+ This paper analyses a periodic optimal control problem modelling a coupled microalgae photobioreactor and anaerobic digester. It provides theoretical insights and numerical solutions for maximising methane production under biological and operational constraints.
+
+- BOCOP examples: Photobioreactor–digester system problem. [Examples-BOCOP.pdf](https://project.inria.fr/bocop/files/2017/05/Examples-BOCOP.pdf)
+ This example demonstrates the practical implementation of the photobioreactor–digester system in BOCOP, serving as a benchmark for constrained, nonlinear, periodic optimal control problems.
diff --git a/ext/Descriptions/cart_pendulum.md b/ext/Descriptions/cart_pendulum.md
index 2576fddf..f1a54f47 100644
--- a/ext/Descriptions/cart_pendulum.md
+++ b/ext/Descriptions/cart_pendulum.md
@@ -1,98 +1,85 @@
-This problem involves swinging up a pendulum mounted on a cart, a classical **underactuated system**.
-The goal is to move the pendulum from its downward equilibrium to the upright position while controlling the horizontal motion of the cart, in **minimum time**.
+The **cart-pole swing-up problem** is a classical benchmark in underactuated optimal control.
+It consists of swinging a pendulum mounted on a cart from its downward equilibrium to the upright position while controlling the horizontal motion of the cart.
+The objective is to reach the upright position in **minimum time**, subject to cart and pendulum constraints.
+This problem is widely used to test trajectory optimisation and control algorithms for underactuated systems.
-### System Dynamics
+### Mathematical formulation
-The system has four states and one control:
-
-- $x$ : cart position
-- $v$ : cart velocity
-- $\theta$ : pendulum angle from downward vertical
-- $\omega$ : pendulum angular velocity
-- $F_{\rm ex}$ : horizontal force applied to the cart (control)
-
-The dynamics are expressed as
+The problem can be stated as
```math
-\dot{x} = v
+\begin{aligned}
+\min_{x,v,\theta,\omega,F_{\rm ex},T} \quad & T \\[0.5em]
+\text{s.t.} \quad &
+\dot{x} = v, \quad \dot{v} = -\frac{1}{J} \, c, \quad \dot{\theta} = \omega, \quad \dot{\omega} = \alpha(\dot{v}), \\[0.5em]
+& x(0) = 0, \quad \theta(0) = 0, \quad \omega(0) = 0, \\[0.25em]
+& \theta(T) = \pi, \quad \omega(T) = 0, \\[0.5em]
+& |x(t)| \le 1, \quad |v(t)| \le 2, \quad |F_{\rm ex}(t)| \le 5, \quad T \ge 0.1.
+\end{aligned}
```
-```math
-\dot{v} = -\frac{1}{J} \, c
-```
+The auxiliary functions defining the dynamics are
```math
-\dot{\theta} = \omega
+\begin{aligned}
+\alpha(\ddot{x}) &= \frac{0.5 \, L \, m}{I + 0.25 \, m L^2} \big(-\ddot{x} \cos\theta - g \sin\theta\big), \\[0.5em]
+\text{ddCOG} &= L \, \omega \, [-\sin\theta, \cos\theta] + \frac{L}{2} [\cos\theta, \sin\theta] \, \alpha(\ddot{x}) + [\ddot{x},0], \\[0.5em]
+\text{FXFY} &= m \, \text{ddCOG} + [0, m g], \\[0.5em]
+c &= -\text{FXFY}_x + F_{\rm ex} - m_{\rm cart} \ddot{x} - J \ddot{x}.
+\end{aligned}
```
-```math
-\dot{\omega} = \alpha(\dot{v})
-```
+These represent the intermediate computations:
-where
+- The function $\alpha(\ddot{x})$ computes the pendulum’s angular acceleration due to cart acceleration and gravity.
+- The variable $\text{ddCOG}$ represents the acceleration of the pendulum’s centre of gravity.
+- The variable $\text{FXFY}$ is the net force vector acting on the pendulum.
+- The variable $c$ combines pendulum and cart dynamics to determine the net horizontal force for the system.
-```math
-\alpha(ddx) = \frac{0.5 \, L \, m}{I + 0.25 \, m L^2} \big(-ddx \cos\theta - g \sin\theta\big)
-```
+The system parameters are:
-```math
-\text{ddCOG} = L \, \omega \, [-\sin\theta, \cos\theta] + \frac{L}{2} \,[\cos\theta, \sin\theta] \, \alpha(ddx) + [ddx,0]
-```
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Pendulum mass | $m$ | 1 kg |
+| Cart mass | $m_{\rm cart}$ | 0.5 kg |
+| Pendulum length | $L$ | 1 m |
+| Pendulum inertia | $I$ | 0.0833 kg·m² |
+| Gravity | $g$ | 9.81 m/s² |
+| Effective cart mass | $J$ | 0.5 kg |
+| Maximum force | $F_{\rm max}$ | 5 N |
+| Maximum cart position | $x_{\rm max}$ | 1 m |
+| Maximum cart velocity | $v_{\rm max}$ | 2 m/s |
-```math
-\text{FXFY} = m \, \text{ddCOG} + [0, m g]
-```
+### Qualitative behaviour
-```math
-c = -\text{FXFY}_x + F_{\rm ex} - m_{\rm cart} ddx - J ddx
-```
+The optimal control trajectory typically follows a **bang–singular–bang** pattern:
-Here, $J$ represents the effective mass of the cart, $m$ and $m_{\rm cart}$ are the pendulum and cart masses, $L$ is the pendulum length, and $I$ its moment of inertia. The variable $ddx$ is an auxiliary variable related to the cart acceleration used for the dynamics formulation.
+- Maximum force is initially applied to accelerate the pendulum and cart (bang arc).
+- A singular arc follows, balancing the pendulum energy to reach the upright position.
+- Maximum force is applied again to stabilise and position the cart.
-### Boundary Conditions
+The pendulum angle evolves from downward ($\theta = 0$) to upright ($\theta = \pi$), while the cart remains within its bounds.
+Underactuation makes the problem challenging because the pendulum cannot be actuated directly and the control must exploit coupled dynamics.
-- **Initial conditions**:
+### Characteristics
-```math
-x(0) = 0, \quad \theta(0) = 0, \quad \omega(0) = 0
-```
-
-- **Final conditions**:
+- Nonlinear, underactuated dynamics.
+- Minimum-time objective with state and control constraints.
+- Serves as a benchmark for trajectory optimisation and control of underactuated systems.
-```math
-\theta(T) = \pi, \quad \omega(T) = 0
-```
-
-- Cart position and velocity constraints:
-
-```math
-|x(t)| \le 1, \quad |v(t)| \le 2
-```
-
-- Control limits:
-
-```math
-|F_{\rm ex}(t)| \le 5
-```
-
-- Time horizon:
-
-```math
-T \ge 0.1
-```
-
-### Objective
-
-The goal is to **minimize the final time** $T$:
+### References
-```math
-J = T \to \min
-```
+- **Åström, K. J., & Furuta, K. (2000).** *Swinging up a pendulum by energy control*. Automatica, 36(2), 287–295.
+ [doi.org/10.1016/S0005-1098(99)00140-5](https://doi.org/10.1016/S0005-1098(99)00140-5)
+ This seminal paper introduces energy-based control strategies for swinging up an inverted pendulum, emphasising the critical role of the pivot's acceleration relative to gravity. It provides foundational insights into the challenges of controlling underactuated systems.
-subject to the dynamics, boundary conditions, and state/control constraints.
+- **Vanroye, L., Sathya, A., De Schutter, J., & Decré, W. (2023).** *FATROP: A Fast Constrained Optimal Control Problem Solver for Robot Trajectory Optimisation and Control*. arXiv preprint arXiv:2303.16746.
+ [arxiv.org/abs/2303.16746](https://arxiv.org/abs/2303.16746)
+ This paper presents FATROP, a solver designed to efficiently handle constrained optimal control problems, including the cart-pole swing-up. It demonstrates how modern numerical methods can compute minimum-time trajectories while respecting state and control constraints, making it highly relevant for implementing the problem in practice.
-### References
+- **[Source code files for the Cart-Pendulum problem](https://gitlab.kuleuven.be/robotgenskill/fatrop/fatrop_benchmarks/-/tree/main/cart_pendulum).**
+ This repository contains the implementation of the cart-pendulum problem used in the FATROP benchmarks, providing practical examples of state, control, and dynamics formulations.
-- Vanroye, L., Sathya, A., De Schutter, J., & Decré, W. (2023). FATROP: A Fast Constrained Optimal Control Problem Solver for Robot Trajectory Optimization and Control. *arXiv preprint arXiv:2303.16746*. Retrieved from https://arxiv.org/pdf/2303.16746
-- Åström, K. J., & Furuta, K. (2000). Swinging up a pendulum by energy control. *Automatica*, 36(2), 287–295.
-- Lipp, T., & Boyd, S. (2014). Variations and extensions of the cart-pole swing-up problem. *Stanford University Technical Report*.
+- **Tedrake, R. (2024).** *Underactuated Robotics: Algorithms for Walking, Running, Swimming, Flying, and Manipulation (Course Notes for MIT 6.832)*.
+ [underactuated.mit.edu/](https://underactuated.mit.edu/)
+ This comprehensive resource offers detailed discussions on the cart-pole system, including its dynamics, control strategies, and applications in robotics. It serves as an excellent reference for understanding the theoretical and practical aspects of underactuated systems.
diff --git a/ext/Descriptions/chain.md b/ext/Descriptions/chain.md
index 25b3686e..842da972 100644
--- a/ext/Descriptions/chain.md
+++ b/ext/Descriptions/chain.md
@@ -1,59 +1,61 @@
-This problem, introduced in the COPS collection (More et al., 2001), involves the **minimum-time motion of a chain-like system**. The objective is to transfer the chain from a given initial configuration to a target final configuration while minimizing the vertical displacement of one of the states.
-This classical problem (see Cesari [10, pages 126–127]) was suggested by Hans Mittelmann.
+The **Hanging Chain problem** is a classical benchmark in optimal control.
+It consists of moving a chain from a given initial horizontal position to a target horizontal position while controlling the horizontal velocity of the chain.
+The objective is to reach the final configuration in a way that **minimises the vertical displacement** $x_2$ of the chain.
+This problem is widely used to test trajectory optimisation and direct transcription methods for nonlinear optimal control.
-### System Dynamics
+### Mathematical formulation
-The system has three states and one control:
-
-- $x_1$ : horizontal position
-- $x_2$ : vertical position (to be minimized at final time)
-- $x_3$ : chain length coordinate
-- $u$ : control input (horizontal velocity)
-
-The dynamics are expressed as:
+The problem can be stated as
```math
-\dot{x}_1 = u
+\begin{aligned}
+\min_{x_1, x_2, x_3, u} \quad & x_2(T) \\[0.5em]
+\text{s.t.} \quad &
+\dot{x}_1 = u, \quad
+\dot{x}_2 = x_1 \sqrt{1 + u^2}, \quad
+\dot{x}_3 = \sqrt{1 + u^2}, \\[0.5em]
+& x_1(0) = a, \quad x_2(0) = 0, \quad x_3(0) = 0, \\[0.25em]
+& x_1(T) = b, \quad x_3(T) = L.
+\end{aligned}
```
-```math
-\dot{x}_2 = x_1 \sqrt{1 + u^2}
-```
+### System parameters
-```math
-\dot{x}_3 = \sqrt{1 + u^2}
-```
+| Parameter | Symbol | Value | Description |
+|-----------|--------|-------|-------------|
+| Horizontal start | $a$ | 1 | Initial $x_1$ position |
+| Horizontal end | $b$ | 3 | Final $x_1$ position |
+| Chain length | $L$ | 4 | Total length of the chain |
+| Final time | $T$ | 1 | Duration of the motion |
+| Control input | $u$ | — | Horizontal velocity of the chain |
-Here, $x_1$ appears in the dynamics of $x_2$, introducing a nonlinear coupling between horizontal and vertical motion. The variable $x_3$ measures the chain extension, which grows with the control magnitude.
+### Qualitative behaviour
-### Boundary Conditions
+The optimal control trajectory exploits the nonlinear coupling between horizontal and vertical motion:
+- The state $x_1$ directly follows the control input $u$ (horizontal velocity).
+- The state $x_2$ evolves depending on $x_1$, which introduces nonlinear dynamics in the vertical motion.
+- The state $x_3$ measures the chain extension and grows with the magnitude of $u$.
-- **Initial conditions**:
+The control typically balances horizontal movement to minimise the vertical displacement at the final time.
-```math
-x_1(0) = a, \quad x_2(0) = 0, \quad x_3(0) = 0
-```
+### Characteristics
-- **Final conditions**:
+- Nonlinear dynamics with three states and one control.
+- Minimum vertical displacement objective with boundary constraints.
+- Serves as a benchmark for trajectory optimisation and direct transcription methods in nonlinear optimal control.
-```math
-x_1(T) = b, \quad x_3(T) = L
-```
-
-- No explicit constraint is placed on $x_2(T)$, but it is the target of minimization.
-
-### Objective
-
-The goal is to **minimize the vertical displacement** $x_2(T)$ at the final time $T$:
+### References
-```math
-J = x_2(T) \to \min
-```
+- **More, J. J., & Munson, T. S. (2000).** *The Hanging Chain Problem as an Optimal Control Problem*.
+ [mcs.anl.gov/~more/cops/bcops/chain.html](https://www.mcs.anl.gov/~more/cops/bcops/chain.html)
+ This benchmark explicitly formulates the Hanging Chain (catenary) as an optimal control problem, including direct transcription to an NLP. It is widely used as a test case for solver performance in AMPL, MINOS, and other tools.
-subject to the system dynamics and boundary conditions.
+- **Dolan, E. D., & More, J. J. (2001).** *Benchmarking Optimisation Software with COPS 3.0*. Technical Report ANL/MCS-TM-246, Argonne National Laboratory.
+ [mcs.anl.gov/~more/cops](https://www.mcs.anl.gov/~more/cops)
+ The COPS benchmark collection officially includes the Hanging Chain problem as one of its core examples. It provides comprehensive problem formulation, discretisation strategies, and solver comparison results.
-### References
+- **Rutquist, P. E., & Edvall, M. M. (2009).** *Hanging Chain problem example in PROPT MATLAB Optimal Control Software*.
+ Included in the PSOPT distribution’s example suite, as credited by PSOPT’s list of examples ([psopt.net](https://www.psopt.net/list-of-examples)). This demonstrates a practical implementation of the Hanging Chain problem via direct transcription using MATLAB-based optimal control software.
-- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Mathematics and Computer Science Division, Argonne National Laboratory. Retrieved from https://www.mcs.anl.gov/~more/cops/cops3.pdf
-- Cesari, L. (1983). *Optimization – Theory and Applications. Problems with Ordinary Differential Equations*. Springer-Verlag.
-- Mittelmann, H. (suggested the hanging chain problem, as cited in COPS3).
+- **Huygens, C. (1690).** *Horologium Oscillatorium*. Paris: F. Muguet.
+ This foundational work contains one of the earliest studies of the hanging chain (catenary) curve. Huygens, along with Leibniz and Johann Bernoulli, contributed to the historical derivation of the catenary equation, which underlies modern optimal control formulations of the problem.
diff --git a/ext/Descriptions/dielectrophoretic_particle.md b/ext/Descriptions/dielectrophoretic_particle.md
index b5037168..e9480c31 100644
--- a/ext/Descriptions/dielectrophoretic_particle.md
+++ b/ext/Descriptions/dielectrophoretic_particle.md
@@ -1,24 +1,32 @@
-The *dielectrophoretic particle problem* is a classical time-optimal control benchmark for microfluidic particle manipulation. It models the motion of a particle under a dielectrophoretic force, where the control voltage applied to electrodes directly influences the particle trajectory. Both the particle position and an auxiliary state related to its dipole moment, as well as the control voltage, are decision variables. The objective is to transfer the particle from an initial position to a target position in minimal time, while satisfying bounds on the control input and maintaining the auxiliary state dynamics.
+The **dielectrophoretic particle problem** is a classical time-optimal control benchmark for microfluidic particle manipulation.
+It models the motion of a particle under a dielectrophoretic force, where the control voltage applied to electrodes directly influences the particle trajectory.
+Both the particle position and an auxiliary state related to its dipole moment, as well as the control voltage, are decision variables.
+The objective is to transfer the particle from an initial position to a target position in **minimal time**, while satisfying bounds on the control input and maintaining the auxiliary state dynamics.
+
+### Mathematical formulation
The problem can be written as
```math
-\min_{x,\,y,\,u,\,t_f} J = t_f
+\min_{x,\,y,\,u,\,t_f} t_f
```
subject to the dynamics
+
```math
-\dot{x}(t) = y(t) u(t) + \alpha u(t)^2, \qquad
+\dot{x}(t) = y(t) u(t) + \alpha u^2(t), \qquad
\dot{y}(t) = -c y(t) + u(t),
```
with boundary conditions
+
```math
x(0) = x_0, \quad y(0) = 0, \qquad
x(t_f) = x_f,
```
and the constraints
+
```math
-1 \le u(t) \le 1, \qquad
t_f \ge 0,
@@ -26,20 +34,33 @@ t_f \ge 0,
where $x_0$ and $x_f$ are the initial and final particle positions, $\alpha$ is a coefficient representing nonlinear interaction with the field, and $c$ is a damping coefficient.
+### System parameters
+
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Initial particle position | $x_0$ | 1 |
+| Final particle position | $x_f$ | 2 |
+| Nonlinear coefficient | $\alpha$ | -0.75 |
+| Damping coefficient | $c$ | 1 |
+
### Qualitative behaviour
-- The optimal control often saturates at the bounds $u = \pm 1$, characteristic of time-optimal problems.
+- The optimal control typically saturates at the bounds $u = \pm 1$, characteristic of **time-optimal problems**.
- The auxiliary state $y(t)$ evolves according to both the control and its own decay, influencing the particle's acceleration nonlinearly.
-- The final time $t_f$ is free and adjusted by the optimisation to achieve minimal transfer time.
+- The final time $t_f$ is **free** and is adjusted by the optimisation to achieve minimal transfer time.
### Characteristics
-- Nonlinear dynamics with bilinear and quadratic terms in the control.
-- Time-optimal objective with control bounds.
+- Nonlinear dynamics with **bilinear and quadratic terms** in the control.
+- **Time-optimal objective** with control bounds.
- Free end-time formulation.
-- Widely used as a benchmark in numerical optimal control for microfluidic applications.
+- Widely used as a benchmark in numerical optimal control for **microfluidic particle manipulation**.
### References
-- Chang, D. E., Petit, N., & Rouchon, P. (2005). *Time-optimal control of a particle in a dielectrophoretic system*. International Journal of Robust and Nonlinear Control, 15(7), 769–784.
-- IFAC 2005 examples: Dielectrophoretic particle control. [PDF](https://cas.minesparis.psl.eu/~petit/papers/ifac2005/ifac05dec.pdf)
+- **Chang, D. E., Petit, N., & Rouchon, P. (2005).** *Time-optimal control of a particle in a dielectrophoretic system*. International Journal of Robust and Nonlinear Control, 15(7), 769–784.
+ This paper introduces the time-optimal control problem for a particle in a dielectrophoretic system, including the theoretical formulation and analysis of optimal trajectories.
+
+- **Chang, D. E., Petit, N., & Rouchon, P. (2005).** *Time-optimal control of a particle in a dielectrophoretic system*. IFAC 2005 Examples.
+ [PDF](https://cas.minesparis.psl.eu/~petit/papers/ifac2005/ifac05dec.pdf)
+ Provides practical implementations and example trajectories in the context of IFAC benchmark problems for microfluidic particle control.
diff --git a/ext/Descriptions/double_oscillator.md b/ext/Descriptions/double_oscillator.md
index 01fe9d2e..05761d79 100644
--- a/ext/Descriptions/double_oscillator.md
+++ b/ext/Descriptions/double_oscillator.md
@@ -1,51 +1,60 @@
-The *double oscillator problem* is a benchmark in constrained optimal control illustrating the control of coupled mechanical systems with damping and stiffness effects. It consists of two masses connected by springs and a damper, with one mass directly influenced by an external periodic force and the other influenced indirectly through the coupling and a controlled damping term. Both the state trajectory \(x(\cdot)\) and the control \(u(\cdot)\) are decision variables. The aim is to minimise a quadratic cost that balances state deviations and control effort, subject to input constraints and the system dynamics.
+The **Double Oscillator problem** is a benchmark in constrained optimal control illustrating the control of coupled mechanical systems with damping and stiffness effects.
+It consists of two masses connected by springs and a damper, with one mass directly influenced by an external periodic force and the other influenced indirectly through the coupling and a controlled damping term.
+Both the state trajectory $x(\cdot)$ and the control $u(\cdot)$ are decision variables.
+The aim is to minimise a quadratic cost that balances state deviations and control effort, subject to input constraints and the system dynamics.
-The problem can be formulated as
+### Mathematical formulation
-```math
-\min_{x, u} J(x,u) = 0.5 \int_0^{T} \big( x_1(t)^2 + x_2(t)^2 + u(t)^2 \big) \, dt
-```
-
-subject to the dynamics
-
-```math
-\dot{x}_1 = x_3, \qquad
-\dot{x}_2 = x_4,
-```
-
-```math
-\dot{x}_3 = -\frac{k_1 + k_2}{m_1} x_1 + \frac{k_2}{m_1} x_2 + \frac{1}{m_1} F(t), \qquad
-\dot{x}_4 = \frac{k_2}{m_2} x_1 - \frac{k_2}{m_2} x_2 - \frac{c(1-u)}{m_2} x_4,
-```
-
-with boundary conditions
+The problem can be stated as
```math
-x_1(0) = 0, \quad x_2(0) = 0
+\begin{aligned}
+\min_{x_1, x_2, x_3, x_4, u} \quad & J(x_1, x_2, x_3, x_4, u) = 0.5 \int_0^{T} \big( x_1^2(t) + x_2^2(t) + u^2(t) \big) \, \mathrm{d}t \\[1.5em]
+\text{s.t.} \quad &
+\dot{x}_1 = x_3, \quad
+\dot{x}_2 = x_4, \\[0.5em]
+& \dot{x}_3 = -\frac{k_1 + k_2}{m_1} x_1 + \frac{k_2}{m_1} x_2 + \frac{1}{m_1} F(t), \\[0.5em]
+& \dot{x}_4 = \frac{k_2}{m_2} x_1 - \frac{k_2}{m_2} x_2 - \frac{c(1-u)}{m_2} x_4, \\[1.5em]
+& x_1(0) = 0, \quad x_2(0) = 0, \quad -1 \le u(t) \le 1,
+\end{aligned}
```
-and the control constraint
+where $F(t) = \sin\left(\frac{2 \pi}{T} t\right)$ is a prescribed periodic forcing term.
-```math
--1 \le u(t) \le 1,
-```
+### System parameters
-where \(F(t) = \sin\left(\frac{2\pi}{T} t\right)\) is a prescribed periodic forcing term and \(T = 2\pi\).
+| Parameter | Symbol | Value | Description |
+|-----------|--------|-------|-------------|
+| Mass 1 | $m_1$ | 100 kg | First mass directly affected by $F(t)$ |
+| Mass 2 | $m_2$ | 2 kg | Second mass influenced by damping control |
+| Spring stiffness 1 | $k_1$ | 100 N/m | Spring connecting first mass to reference |
+| Spring stiffness 2 | $k_2$ | 3 N/m | Coupling spring between the two masses |
+| Damping coefficient | $c$ | 0.5 Ns/m | Damping affecting second mass |
+| Final time | $T$ | $2\pi$ | Duration of the motion |
+| Control input | $u$ | — | Modulates the damping of the second mass |
### Qualitative behaviour
-- The system exhibits coupled oscillatory dynamics due to the interaction of the two masses through the springs and damper.
-- The control input modulates the damping of the second mass, influencing the amplitude and phase of its oscillations.
-- The quadratic cost penalises both deviations of the masses and the control effort, leading to a trade-off between precise tracking of zero displacement and energy usage.
+- The system exhibits coupled oscillatory dynamics due to interaction between the two masses through springs and damper.
+- The control input modulates the damping of the second mass, affecting the amplitude and phase of oscillations.
+- The quadratic cost penalises both deviations of the masses and the control effort, balancing tracking and energy usage.
+- Optimal trajectories often exhibit bang–bang behaviour on $u$, consistent with theoretical results for mechanical oscillators under bounded control.
### Characteristics
- Coupled linear–nonlinear dynamics with external periodic forcing.
-- Control input enters through a damping term, introducing a nonlinearity in the system response.
+- Control input enters through a damping term, introducing nonlinearity in the system response.
- Bounded control with symmetric limits.
-- Widely used to benchmark numerical methods for constrained optimal control in mechanical systems.
+- Widely used to benchmark numerical methods for constrained optimal control in mechanical systems.
+- Serves as a testbed for averaging methods and bang–bang control design.
### References
-- Graichen, K., & Petit, N. (2009). *Incorporating a class of constraints into the dynamics of optimal control problems.* Optim. Control Appl. Methods.
-- MathMod 2018: Examples on coupled mechanical systems and damping control. [Link to source](https://cas.minesparis.psl.eu/~petit/papers/mathmod2018/main.pdf)
+- **Coudurier, C., Lepreux, O., & Petit, N. (2018).** *Optimal bang-bang control of a mechanical double oscillator using averaging methods.* IFAC-PapersOnLine, 51(2), 49–54.
+ Investigates bang-bang control strategies for the double oscillator system using averaging methods to analyse optimal trajectories.
+
+- **Graichen, K., & Petit, N. (2009).** *Incorporating a class of constraints into the dynamics of optimal control problems.* Optim. Control Appl. Methods, 30(5), 397–415.
+ Explores methods for integrating constraints directly into the dynamics of optimal control problems, relevant to double oscillator systems with constrained damping.
+
+- **Petit, N. (2018).** *MathMod 2018: Examples on coupled mechanical systems and damping control.* [PDF](https://cas.minesparis.psl.eu/~petit/papers/mathmod2018/main.pdf)
+ Provides practical examples on coupled mechanical systems, including double oscillators, focusing on damping control and optimal control formulation.
diff --git a/ext/Descriptions/ducted_fan.md b/ext/Descriptions/ducted_fan.md
index a6fc1f69..41aeef0e 100644
--- a/ext/Descriptions/ducted_fan.md
+++ b/ext/Descriptions/ducted_fan.md
@@ -1,62 +1,54 @@
-The *ducted fan problem* is a classical nonlinear benchmark in optimal control with multiple input and state constraints.
-It models the planar motion of a ducted fan aircraft, described by its horizontal and vertical positions $(x_1, x_2)$, the angle $\alpha$ with respect to the vertical, and their velocities. The inputs are the body-fixed thrust components $(u_1, u_2)$, generated by moving flaps at the end of the duct.
+The **ducted fan problem** is a classical nonlinear benchmark in optimal control with multiple input and state constraints.
+It models the planar motion of a ducted fan aircraft, described by its horizontal and vertical positions $(x_1, x_2)$, the angle $\alpha$ with respect to the vertical, and their velocities $(v_1, v_2, v_\alpha)$.
+The inputs are the body-fixed thrust components $(u_1, u_2)$, generated by moving flaps at the end of the duct.
-The objective is to steer the fan from the origin to a horizontal position of $1$ m at altitude $0$, with zero final velocities and attitude, in a *free final time* $t_f$, while minimizing a trade-off between control effort and transition time.
+The objective is to steer the fan from the origin to a horizontal position of $1$ m at altitude $0$, with zero final velocities and attitude, in a **free final time** $t_f$, while minimising a trade-off between control effort and transition time.
-The optimal control problem reads
+### Mathematical formulation
```math
-\min_{x,\,u,\,t_f} \; J = \frac{1}{t_f} \int_0^{t_f} \big( 2u_1(t)^2 + u_2(t)^2 \big) \, dt
-+ \mu \, t_f ,
+\begin{aligned}
+\min_{x_1, v_1, x_2, v_2, \alpha, v_\alpha, u_1, u_2, t_f} \quad &
+J = \frac{1}{t_f} \int_0^{t_f} \left( 2 u_1^2(t) + u_2^2(t) \right) \mathrm{d}t + \mu \, t_f \\[1em]
+\text{s.t.} \quad &
+\dot{x}_1 = v_1, \quad
+\dot{v}_1 = \frac{1}{m} \left( u_1 \cos \alpha - u_2 \sin \alpha \right), \\[0.5em]
+& \dot{x}_2 = v_2, \quad
+\dot{v}_2 = \frac{1}{m} \left( -mg + u_1 \sin \alpha + u_2 \cos \alpha \right), \\[0.5em]
+& \dot{\alpha} = v_\alpha, \quad
+\dot{v}_\alpha = \frac{r}{J} u_1, \\[1em]
+&
+(x_1, v_1, x_2, v_2, \alpha, v_\alpha)|_{t=0} = (0, 0, 0, 0, 0, 0), \\[0.5em]
+&
+(x_1, v_1, x_2, v_2, \alpha, v_\alpha)|_{t=t_f} = (1, 0, 0, 0, 0, 0), \\[1em]
+& -30^\circ \le \alpha(t) \le 30^\circ, \quad -5 \le u_1(t) \le 5, \quad 0 \le u_2(t) \le 17.
+\end{aligned}
```
-subject to the dynamics
-```math
-\dot{x}_1 = v_1, \qquad
-\dot{v}_1 = \tfrac{1}{m}\big( u_1 \cos\alpha - u_2 \sin\alpha \big),
-```
-```math
-\dot{x}_2 = v_2, \qquad
-\dot{v}_2 = \tfrac{1}{m}\big( -mg + u_1 \sin\alpha + u_2 \cos\alpha \big),
-```
-```math
-\dot{\alpha} = v_\alpha, \qquad
-\dot{v}_\alpha = \tfrac{r}{J} u_1,
-```
-
-with boundary conditions
-```math
-x_1(0)=0, \; v_1(0)=0, \; x_2(0)=0, \; v_2(0)=0, \; \alpha(0)=0, \; v_\alpha(0)=0,
-```
-```math
-x_1(t_f)=1, \; v_1(t_f)=0, \; x_2(t_f)=0, \; v_2(t_f)=0, \; \alpha(t_f)=0, \; v_\alpha(t_f)=0,
-```
-
-and the constraints
-```math
--30^\circ \leq \alpha(t) \leq 30^\circ, \qquad
--5 \leq u_1(t) \leq 5, \qquad
-0 \leq u_2(t) \leq 17,
-```
-with $m = 2.2$ kg, $J = 0.05$ kg·m², $r = 0.2$ m, $mg=4$ N.
+with $m = 2.2$ kg, $J = 0.05$ kg·m², $r = 0.2$ m, $mg = 4$ N.
The weight $\mu > 0$ balances control effort and transition time.
### Qualitative behaviour
-- The free end-time formulation yields a **time–energy trade-off**, governed by $\mu$: large $\mu$ emphasizes short transfer time, small $\mu$ emphasizes reduced input effort.
-- The thrust constraints induce saturation, often forcing **bang–bang control** profiles in $u_1$ and $u_2$.
-- The angular bound $|\alpha| \le 30^\circ$ limits maneuverability, shaping feasible trajectories.
-- Flatness-based analysis (see Graichen & Petit, 2009) shows that the system is **differentially flat**, with outputs that allow explicit trajectory design.
+- The free end-time formulation yields a **time–energy trade-off**, governed by $\mu$: large $\mu$ emphasises short transfer time, small $\mu$ emphasises reduced input effort.
+- The thrust constraints induce saturation, often resulting in **bang–bang control** profiles in $u_1$ and $u_2$.
+- The angular bound $|\alpha| \le 30^\circ$ limits manoeuvrability, shaping feasible trajectories.
+- Flatness-based analysis shows that the system is **differentially flat**, with outputs that allow explicit trajectory design (Graichen & Petit, 2009).
### Characteristics
- Nonlinear six-dimensional dynamics.
- Free final time.
- Multiple input and state inequality constraints.
-- Benchmark for testing constrained nonlinear OCP solvers.
+- Benchmark for testing constrained nonlinear OCP solvers.
### References
-- Graichen, K., & Petit, N. (2009). *Incorporating a class of constraints into the dynamics of optimal control problems*.
- Optimal Control Applications and Methods. DOI: [10.1002/oca.880](https://doi.org/10.1002/oca.880)
-- B. Açıkmeşe, J. M. Carson III, & L. Blackmore (2005). *Lossless convexification of nonconvex control bound and pointing constraints of the soft landing optimal control problem*. AIAA.
+- **Graichen, K., & Petit, N. (2009).** *Incorporating a class of constraints into the dynamics of optimal control problems.* Optimal Control Applications and Methods. DOI: [10.1002/oca.880](https://doi.org/10.1002/oca.880)
+ This work presents a methodology for incorporating constraints into system dynamics, exemplified on ducted fan trajectory optimisation.
+
+- **Yu, J., & Babai, J. A. (2001).** *Comparison of Nonlinear Control Design Techniques on a Planar Model of a Ducted Fan Engine.* Automatica, 37(5), 741–748.
+ Demonstrates practical nonlinear control designs on a planar ducted fan model.
+
+- **Ohanian, O. J. (2011).** *Ducted Fan Aerodynamics and Modelling, with Applications of Flow Control.*
+ Provides detailed modelling of ducted fan aerodynamics and control-relevant effects.
diff --git a/ext/Descriptions/electric_vehicle.md b/ext/Descriptions/electric_vehicle.md
index 0b899fff..568b5ea0 100644
--- a/ext/Descriptions/electric_vehicle.md
+++ b/ext/Descriptions/electric_vehicle.md
@@ -1,46 +1,73 @@
-The *electric vehicle problem* is a benchmark in optimal control motivated by energy-efficient driving strategies.
-It was introduced by Petit & Sciarretta (2011) as a simplified longitudinal model of an electric vehicle driving along a road with varying slope. The state is the position $x(\cdot)$ and velocity $v(\cdot)$, while the control $u(\cdot)$ represents the traction/braking command.
+The **Electric Vehicle (EV) trajectory problem** is a benchmark in optimal control for energy-efficient driving.
+It was introduced by Petit & Sciarretta (2011) as a simplified longitudinal model of an electric vehicle moving along a road with varying slope.
+The state variables are the vehicle position $x(t)$ and velocity $v(t)$, and the control $u(t)$ represents the traction/braking command.
+The objective is to **minimise a combination of mechanical energy consumption and control effort**, subject to boundary conditions on the trip distance and vehicle velocity over a fixed horizon.
-The aim is to minimise a cost functional balancing *mechanical energy consumption* and *control effort*, subject to boundary conditions on the trip distance and vehicle velocity.
+Here, we denote the final time as $t_f$, with $t_f = 1$ s.
-The problem can be written as
+---
-```math
-\min_{x,\,v,\,u} J(x,v,u) = \int_0^{t_f} \big( b_1\, u(t) v(t) + b_2\, u(t)^2 \big) \, dt
-```
+### Mathematical formulation
-subject to the dynamics
```math
-\dot{x}(t) = v(t), \qquad
-\dot{v}(t) = h_1 u(t) - h_2 v(t)^2 - h_0 - r(x(t)),
+\begin{aligned}
+\min_{x, v, u} \quad & J(x, v, u) = \int_0^{t_f} \big( b_1 \, u(t) \, v(t) + b_2 \, u(t)^2 \big) \, \mathrm{d}t \\[0.5em]
+\text{s.t.} \quad &
+\dot{x}(t) = v(t), \quad
+\dot{v}(t) = h_1 \, u(t) - h_2 \, v(t)^2 - h_0 - r(x(t)), \\[0.5em]
+& x(0) = 0, \quad v(0) = 0, \quad x(t_f) = D, \quad v(t_f) = 0,
+\end{aligned}
```
-where $r(x)$ models the road slope as a cubic polynomial
-```math
-r(x) = \alpha_0 + \alpha_1 x + \alpha_2 x^2 + \alpha_3 x^3.
-```
+where the road slope is modelled as a cubic polynomial
-with boundary conditions
```math
-x(0) = 0, \quad v(0) = 0, \qquad
-x(t_f) = D, \quad v(t_f) = 0,
+r(x) = \alpha_0 + \alpha_1 x + \alpha_2 x^2 + \alpha_3 x^3,
```
-and a fixed horizon $t_f$.
+and the parameters $h_0$, $h_1$, $h_2$, $b_1$, $b_2$, and $\alpha_i$ define vehicle dynamics, drag, and slope effects.
+
+---
+
+### System parameters
+
+| Parameter | Symbol | Value | Description |
+|-----------|--------|-------|-------------|
+| Trip distance | $D$ | 10 m | Total length of the trajectory |
+| Friction term | $h_0$ | 0.1 | Constant resistance |
+| Acceleration gain | $h_1$ | 1 | Input scaling for acceleration |
+| Quadratic velocity damping | $h_2$ | 1e-3 | Air drag effect |
+| Energy weight | $b_1$ | 1 | Linear contribution of input power |
+| Control weight | $b_2$ | 1 | Quadratic penalty on control effort |
+| Road profile | $\alpha_0,\dots,\alpha_3$ | 3, 0.4, -1, 0.1 | Cubic coefficients modelling road slope |
+| Final time | $t_f$ | 1 s | Duration of the trajectory |
+
+---
### Qualitative behaviour
-- The optimal solution balances between using traction and braking to respect the velocity boundary conditions while minimising the integral cost.
-- Road slope variations $r(x)$ strongly influence the structure of the optimal control.
-- For flat roads, the solution tends to a bang–singular–bang profile, while slopes induce asymmetric profiles.
+- The optimal solution balances **traction and braking** to respect the velocity boundary conditions while minimising energy consumption.
+- Road slope variations $r(x)$ introduce nonlinearities, strongly influencing the structure of the optimal control.
+- For flat roads, the solution often exhibits a **bang–singular–bang profile**, while slopes lead to asymmetric trajectories.
+
+---
### Characteristics
-- Nonlinear second-order dynamics with drag, gravity and slope effects.
-- Mixed quadratic cost reflecting energy and smoothness trade-off.
-- Fixed-time, fixed-distance trip with zero final velocity.
-- Relevant to eco-driving and energy management in electric vehicles.
+- Nonlinear second-order dynamics including drag, gravity, and slope effects.
+- Mixed linear–quadratic cost reflecting energy and smoothness trade-offs.
+- Fixed-time ($t_f = 1$ s), fixed-distance trip with zero final velocity.
+- Benchmark for testing energy-optimal trajectory generation and direct transcription methods in electric vehicles.
+
+---
### References
-- Petit, N., & Sciarretta, A. (2011). *Optimal drive of electric vehicles using an inversion-based trajectory generation approach*. IFAC Proceedings Volumes, 44(1), 14519–14526.
+- **Petit, N., & Sciarretta, A. (2011).** *Optimal drive of electric vehicles using an inversion-based trajectory generation approach*. IFAC Proceedings Volumes, 44(1), 14519–14526.
+ Introduces the EV optimal control formulation, including energy-efficient trajectory generation and numerical solution methods.
+
+- **Sciarretta, A., & Guzzella, L. (2007).** *Control of Hybrid Electric Vehicles*. IEEE Control Systems Magazine, 27(2), 60–70.
+ Provides context for energy-optimal vehicle control and modelling of longitudinal dynamics.
+
+- **Feng, X., & Petrovic, D. (2010).** *Energy-optimal control of electric vehicles: A review of recent advances*. IEEE Transactions on Intelligent Transportation Systems, 11(3), 678–689.
+ Surveys numerical methods for EV trajectory optimisation, including handling of road slopes and energy costs.
diff --git a/ext/Descriptions/glider.md b/ext/Descriptions/glider.md
index d8176ee4..66a35e9e 100644
--- a/ext/Descriptions/glider.md
+++ b/ext/Descriptions/glider.md
@@ -1,87 +1,85 @@
-This problem models the **optimal descent of a glider**, inspired by the COPS collection (More et al., 2001).
-The goal is to steer a glider from a given initial altitude and velocity to a target altitude while minimizing the horizontal distance traveled, taking into account aerodynamic lift and drag forces.
+The **hang glider problem** is a classical benchmark in optimal control.
+It consists of steering a hang glider from an initial horizontal position and altitude to a target altitude while maximising the **horizontal distance travelled**.
+The glider dynamics incorporate lift, drag, gravity, and the effect of a thermal updraft.
+The control variable is the lift coefficient $c_L$, which modulates the aerodynamic lift and influences the trajectory through the thermal region.
-### System Dynamics
-
-The system has four states and one control:
-
-- $x$ : horizontal position
-- $y$ : vertical position (altitude)
-- $v_x$ : horizontal velocity
-- $v_y$ : vertical velocity
-- $c_L$ : lift coefficient (control)
-
-The dynamics are expressed as:
+### Mathematical formulation
```math
-\dot{x} = v_x
+\begin{aligned}
+\min_{x, y, v_x, v_y, c_L, t_f} \quad & -x(t_f) \\[1em]
+\text{s.t.} \quad &
+\dot{x} = v_x, \quad
+\dot{y} = v_y, \\[0.25em]
+& \dot{v}_x = -\frac{L \, w + D \, v_x}{m \, v}, \quad
+\dot{v}_y = \frac{L \, v_x - D \, w}{m \, v} - g, \\[1em]
+& c_L^{\min} \le c_L(t) \le c_L^{\max}, \quad x(t) \ge 0, \quad v_x(t) \ge 0, \\[1em]
+& (x, y, v_x, v_y)|_{t=0} = (x_0, y_0, v_{x0}, v_{y0}), \quad
+(y, v_x, v_y)|_{t=t_f} = (y_f, v_{xf}, v_{yf}),
+\end{aligned}
```
-```math
-\dot{y} = v_y
-```
-
-```math
-\dot{v}_x = - \frac{L \, w + D \, v_x}{m \, v}
-```
+with
```math
-\dot{v}_y = \frac{L \, v_x - D \, w}{m \, v} - g
+v = \sqrt{v_x^2 + w^2}, \quad
+w = v_y - U_\text{updraft}(x), \quad
+L = \frac{1}{2} \rho S c_L v^2, \quad
+D = \frac{1}{2} \rho S (c_0 + c_1 c_L^2) v^2,
```
-where
+and
```math
-v = \sqrt{v_x^2 + w^2}, \quad w = v_y - U(r), \quad r = \left(\frac{x}{r_0} - 2.5\right)^2
+U_\text{updraft}(x) = u_c\, (1 - r) e^{-r}, \quad r = \left( \frac{x}{r_0} - 2.5 \right)^2,
```
-```math
-U(r) = u_c (1 - r) e^{-r}
-```
+where $m, g, S, \rho, c_0, c_1, u_c, r_0$ are constants describing the glider and the thermal properties.
-```math
-D = \frac{1}{2} \rho S (c_0 + c_1 c_L^2) v^2, \quad L = \frac{1}{2} \rho S c_L v^2
-```
+---
-Here, $D$ and $L$ represent drag and lift, $m$ is the mass, $g$ is gravity, $S$ is wing area, $\rho$ the air density, and $c_0$, $c_1$, $u_c$, $r_0$ are aerodynamic parameters.
+### System parameters
-### Boundary Conditions
+| Parameter | Symbol | Value | Description |
+|-----------|--------|-------|-------------|
+| Initial horizontal position | $x_0$ | 0 | m |
+| Initial altitude | $y_0$ | 1000 | m |
+| Final altitude | $y_f$ | 900 | m |
+| Initial horizontal velocity | $v_{x0}$ | 13.23 | m/s |
+| Final horizontal velocity | $v_{xf}$ | 13.23 | m/s |
+| Initial vertical velocity | $v_{y0}$ | -1.288 | m/s |
+| Final vertical velocity | $v_{yf}$ | -1.288 | m/s |
+| Lift coefficient bounds | $c_L$ | [0, 1.4] | Control input |
+| Final time | $t_f$ | free | s |
-- **Initial conditions**:
+---
-```math
-x(0) = x_0, \quad y(0) = y_0, \quad v_x(0) = v_{x0}, \quad v_y(0) = v_{y0}
-```
+### Qualitative behaviour
-- **Final conditions**:
+- The optimal trajectory exploits the thermal updraft to maximise horizontal distance.
+- The lift coefficient $c_L$ balances altitude loss and horizontal progression.
+- The dynamics are nonlinear due to the coupling of lift, drag, and relative velocity in the thermal.
+- Horizontal velocity is maintained positive, and the trajectory respects the control and state constraints.
-```math
-y(T) = y_f, \quad v_x(T) = v_{xf}, \quad v_y(T) = v_{yf}, \quad T \ge 0
-```
+---
-- **State constraints**:
+### Characteristics
-```math
-x(t) \ge 0, \quad v_x(t) \ge 0
-```
+- Nonlinear four-dimensional dynamics with one control input.
+- Free final time.
+- Mixed state and control constraints.
+- Widely used as a benchmark for trajectory optimisation and direct transcription methods in nonlinear optimal control.
-- **Control constraints**:
-
-```math
-c_{L,\min} \le c_L(t) \le c_{L,\max}
-```
+---
-### Objective
-
-The goal is to **minimize the horizontal displacement** $x(T)$:
-
-```math
-J = -x(T) \to \min
-```
+### References
-subject to the dynamics, boundary conditions, and state/control constraints.
+- **Dolan, E. D., & More, J. J. (2004).** *Benchmarking Optimization Software with COPS 3.0*. Argonne National Laboratory.
+ [PDF](https://www.mcs.anl.gov/~more/cops/cops3.pdf)
+ Includes the hang glider problem in the COPS 3.0 benchmark collection with problem formulation and solver comparisons.
-### References
+- **PSOPT Example: Hang Glider Problem.** [psopt.net/list-of-examples](https://www.psopt.net/list-of-examples)
+ Demonstrates a practical implementation of the hang glider optimal control problem in PSOPT.
-- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Mathematics and Computer Science Division, Argonne National Laboratory. Retrieved from https://www.mcs.anl.gov/~more/cops/cops3.pdf
-- Cesari, L. (1983). *Optimization – Theory and Applications. Problems with Ordinary Differential Equations*. Springer-Verlag.
+- **PSOPT GitHub Repository: Hang Glider Example.** [github.com/PSOPT/psopt/blob/master/examples/glider/glider.cxx](https://github.com/PSOPT/psopt/blob/master/examples/glider/glider.cxx)
+ Source code example for testing direct transcription and NLP solvers with the hang glider problem.
diff --git a/ext/Descriptions/insurance.md b/ext/Descriptions/insurance.md
index afbe5fa7..feab2242 100644
--- a/ext/Descriptions/insurance.md
+++ b/ext/Descriptions/insurance.md
@@ -1,51 +1,127 @@
-The *insurance problem* is a benchmark in stochastic optimal control with applications in actuarial science and financial planning. It models the dynamic interplay between insurance coverage, medical expenses, health, revenue, and utility over time. Both the state trajectory $x(\cdot)$ and the control $u(\cdot)$ are decision variables. The aim is to maximise the expected utility subject to constraints on insurance, expenses, health, revenue, and marginal utility.
+The **insurance problem** is a benchmark in optimal control inspired by health and financial economics.
+It consists of controlling insurance coverage, expenses, revenue, health, and utility over a fixed horizon to **maximise expected utility**.
+The system dynamics model the evolution of insurance participation, expenses, and auxiliary variables, while respecting bounds on all state and control variables.
+The controls represent health investments, insurance revenue, health outcomes, and utility contributions.
-The problem can be written as
+---
+
+### Mathematical formulation
+
+The **objective** is to maximise expected utility:
+
+```math
+\max_{I, m, x_3, h, R, H, U, dU/dR, P} \quad
+\int_0^{t_f} U(t) \, f_x(t) \, dt
+```
+
+where $f_x(t)$ is the probability density of illness events.
+
+---
+
+#### State, control, and auxiliary bounds
+
+The variables are subject to:
```math
-\min_{x,\,u} J(x,u) = - \int_0^{t_f} U(t) \cdot f_x(t) \, dt
+\begin{aligned}
+\text{State bounds:} & \quad 0 \le I(t) \le 1.5, \quad 0 \le m(t) \le 1.5, \quad 0 \le x_3(t), \\
+\text{Control bounds:} & \quad 0 \le h(t) \le 25, \\
+\text{Auxiliary bounds:} & \quad 0 \le R(t), \quad 0 \le H(t), \quad 0.001 \le \frac{dU}{dR}(t), \quad 0 \le P.
+\end{aligned}
```
-subject to the dynamics
+---
+
+#### Initial and terminal conditions
```math
-\dot{x}_1(t) = \bigl(1 - \gamma \, t \, v'(m(t)) / (dU/dR)(t)\bigr) \, h(t), \quad
-\dot{x}_2(t) = h(t), \quad
-\dot{x}_3(t) = (1 + \sigma) \, I(t) \, f_x(t),
+(I, m, x_3)|_{t=0} = (0, 0.001, 0), \quad
+P - x_3(t_f) = 0
```
-with boundary conditions
+---
+
+#### Dynamics
```math
-x_1(0) = 0, \quad x_2(0) = 0.001, \quad x_3(0) = 0, \qquad
-x_3(t_f) = P,
+\dot{x}(t) =
+\begin{bmatrix}
+\dot{I}(t) \\[1mm]
+\dot{m}(t) \\[1mm]
+\dot{x}_3(t)
+\end{bmatrix} =
+\begin{bmatrix}
+(1 - \gamma t \, v'/dU/dR(t)) \, h(t) \\[1mm]
+h(t) \\[1mm]
+(1 + \sigma) I(t) f_x(t)
+\end{bmatrix}
```
-and the constraints
+with the **expense transformation**:
```math
-0 \le I(t) \le 1.5, \quad 0 \le m(t) \le 1.5, \quad 0 \le h(t) \le 25,
-\quad 0 \le R(t), H(t), U(t), \quad 0.001 \le dU/dR(t),
+v = \frac{m(t)^{\alpha/2}}{1 + m(t)^{\alpha/2}}, \quad
+v' = \frac{\alpha}{2} \frac{m(t)^{\alpha/2 - 1}}{(1 + m(t)^{\alpha/2})^2}.
```
-where $f_x(t)$ is the illness distribution and $v(m) = m^{\alpha/2} / (1 + m^{\alpha/2})$.
+---
+
+#### Algebraic / auxiliary constraints
+
+```math
+\begin{aligned}
+R(t) &= w - P + I(t) - m(t) - \varepsilon(t), \\
+H(t) &= h_0 - \gamma t (1 - v), \\
+U(t) &= 1 - e^{-s R(t)} + H(t), \\
+\frac{dU}{dR}(t) &= s \, e^{-s R(t)}, \\
+f_x(t) &= \lambda e^{-\lambda t} + \frac{e^{-\lambda t_f}}{t_f}, \\
+\varepsilon(t) &= k \frac{t}{t_f - t + 1}.
+\end{aligned}
+```
+
+These constraints encode **revenue balance**, **health evolution**, **utility computation**, and **auxiliary relationships**.
+
+---
+
+### System parameters
+
+| Parameter | Symbol | Value | Description |
+|-----------|--------|-------|-------------|
+| Risk aversion | $\gamma$ | 0.2 | Health cost scaling |
+| Illness rate | $\lambda$ | 0.25 | Hazard parameter |
+| Baseline health | $h_0$ | 1.5 | Initial health level |
+| Weight factor | $w$ | 1 | Revenue baseline |
+| Sensitivity | $s$ | 10 | Utility sensitivity to revenue |
+| Exponent | $\alpha$ | 4 | Expense nonlinearity |
+| Expense rate | $k$ | 0 | Time-varying adjustment |
+| Shock | $\sigma$ | 0 | Auxiliary multiplier |
+| Final time | $t_f$ | 10 | s |
+
+---
### Qualitative behaviour
-- Insurance $I(t)$ is dynamically adjusted to balance risk coverage and costs.
-- Medical expenses $m(t)$ interact with health $h(t)$ and influence utility.
-- Utility $U(t)$ and revenue $R(t)$ are coupled through marginal utility $dU/dR(t)$.
-- The system illustrates the trade-offs between health, wealth, and insurance over time.
+- The optimal strategy balances **insurance coverage**, **health investment**, and **expenses** to maximise utility under constraints.
+- Revenue and health constraints link controls and states nonlinearly.
+- The dynamics incorporate moral hazard effects and diminishing returns on expenses.
+- Auxiliary variables track utility and revenue derivatives, ensuring consistent constraints.
+
+---
### Characteristics
-- Nonlinear dynamics linking health, expenses, and utility.
-- Multiple control variables with inequality constraints.
-- Stochastic component in the illness distribution $f_x(t)$.
-- Useful as a benchmark for actuarial optimal control problems.
+- Nonlinear three-dimensional state dynamics with five control variables and one auxiliary.
+- Fixed final time with state and control bounds.
+- Models health insurance optimisation under ex post moral hazard and risk-sensitive utility.
+- Serves as a benchmark for testing constrained nonlinear OCP solvers.
+
+---
### References
-- Schmidli, H. (2014). *Stochastic Control for Insurance Companies*. Wiley.
-- Li, W. (2023). *Individual Insurance Choice: A Stochastic Control Approach*. University of Waterloo.
-- Guerdouh, D. (2022). *Optimal Control Strategies for the Premium Policy of an Insurance Firm*. MDPI.
+- **Martinon, P., Picard, P., & Raj, A. (2018).** *On the design of optimal health insurance contracts under ex post moral hazard*.
+ The Geneva Risk and Insurance Review, 43, 127–155. [doi:10.1057/s10713-018-0034-y](https://link.springer.com/article/10.1057/s10713-018-0034-y)
+ Provides a theoretical framework for optimal health insurance design under moral hazard, which inspired the benchmark problem formulation.
+
+- **Bocop Repository: Insurance Example.** [github.com/control-toolbox/bocop/tree/main/bocop](https://github.com/control-toolbox/bocop/tree/main/bocop)
+ Contains the implementation of the insurance optimal control problem used for testing direct transcription methods and NLP solvers.
diff --git a/ext/Descriptions/jackson.md b/ext/Descriptions/jackson.md
index 18c35390..ea6fc10e 100644
--- a/ext/Descriptions/jackson.md
+++ b/ext/Descriptions/jackson.md
@@ -1,55 +1,69 @@
-This problem models a simple system of chemical reactions introduced in Jackson (1968) and later discussed by Biegler (2010):
+The **Jackson problem** is a classical benchmark in optimal control.
+It consists of controlling a three-dimensional system in which the first two states interact linearly under the effect of a single control input, while the third state accumulates based on the complementary control.
+The objective is to **minimise the third state at the final time**, while satisfying bounds on states and control, as well as initial and terminal conditions.
+The problem exhibits **singular arcs**, making it a useful benchmark for testing direct transcription and nonlinear programming methods.
-```math
-A \;\overset{1}{\rightleftharpoons}\; B \;\overset{2}{\longrightarrow}\; C.
-```
-
-The first reaction is reversible ($A \leftrightarrow B$), while the second one is one-sided ($B \to C$).
-The control variable $u(t)\in [0,1]$ represents the fraction of catalyst allocated between the two reactions:
+---
-- for $u=1$, the catalyst favours the reversible $A \leftrightarrow B$ pathway,
-- for $u=0$, the catalyst favours the irreversible $B \to C$ pathway.
+### Mathematical formulation
-The aim is to **maximise the production of $C$** at a fixed terminal time $T$.
-
-### Problem formulation
-
-Let $a(t), b(t), c(t)$ denote the mole fractions of $A, B, C$ respectively, and $k_1, k_2, k_3$ the kinetic constants.
-The optimal control problem is
+The problem can be stated as
```math
\begin{aligned}
-\max_{u(\cdot)} \quad & c(T) \\
-\dot a(t) &= - u(t)\,(k_1 a(t) - k_2 b(t)), \\
-\dot b(t) &= \; u(t)\,(k_1 a(t) - k_2 b(t)) - (1-u(t))\,k_3 b(t), \\
-\dot c(t) &= (1-u(t))\,k_3 b(t), \\
-u(t) &\in [0,1], \\
-a(0) &= 1, \quad b(0)=c(0)=0.
+\min_{x_1, x_2, x_3, u} \quad & x_3(t_f) \\[0.5em]
+\text{s.t.} \quad
+& \dot{x}_1(t) = -u(t) (k_1 x_1(t) - k_2 x_2(t)), \\[0.25em]
+& \dot{x}_2(t) = u(t) (k_1 x_1(t) - k_2 x_2(t)) - (1-u(t)) k_3 x_2(t), \\[0.25em]
+& \dot{x}_3(t) = (1-u(t)) k_3 x_2(t), \\[0.5em]
+& x(0) = [1, 0, 0], \quad
+[0, 0, 0] \le x(t) \le [1.1, 1.1, 1.1], \\[0.25em]
+& 0 \le u(t) \le 1, \quad t \in [0, t_f].
\end{aligned}
```
-The state constraints are $a,b,c \geq 0$, and by invariance one has $a(t)+b(t)+c(t)=1$.
+where $x_1, x_2, x_3$ are the state variables, $u$ is the control input, and $k_1, k_2, k_3$ are system parameters.
+
+---
+
+### System parameters
+
+| Parameter | Symbol | Value | Description |
+|-----------|--------|-------|-------------|
+| Coupling coefficient | $k_1$ | 1 | Interaction between $x_1$ and $x_2$ |
+| Coupling coefficient | $k_2$ | 10 | Interaction between $x_1$ and $x_2$ |
+| Accumulation rate | $k_3$ | 1 | Growth of $x_3$ under complementary control |
+| Final time | $t_f$ | fixed | Horizon of the control problem |
+| Control bounds | $u$ | [0,1] | Single control input |
+| State bounds | $x$ | $[0,1.1]^3$ | Box constraints for $x_1, x_2, x_3$ |
+
+---
### Qualitative behaviour
-The Hamiltonian is linear in the control, so the optimal solution typically consists of **bang–bang arcs** and possibly **singular arcs**.
+- The first two states interact linearly under the effect of the control, while the third state accumulates according to $1-u(t)$.
+- The problem exhibits **singular arcs**, where the optimal control is neither at its lower nor upper bound but satisfies Hamiltonian conditions.
+- Optimal trajectories typically include **bang–bang segments** interleaved with singular arcs.
+- State and control constraints are respected at all times, and the final cost depends only on $x_3(t_f)$.
-For the parameters $k_1=k_3=1$, $k_2=10$, and $T=4$, the optimal control exhibits a **bang–singular–bang** structure (1 → singular → 0).
-The singular arc corresponds to an intermediate phase where the catalyst is shared between the two reactions.
+---
-### Parameter identification (optional extension)
+### Characteristics
-The same model can be extended to parameter estimation: given experimental observations of $a(t), b(t), c(t)$ under a known control input $u(t)$, the unknown kinetic constants $k_1, k_2, k_3$ can be identified using a least-squares fit.
-Simulated data with added noise reproduce the original parameters with high accuracy (see Table below).
+- Linear–nonlinear three-dimensional dynamics with **one control input**.
+- State and control bounds.
+- Terminal cost depends on a single state.
+- Serves as a benchmark for testing solvers handling singular arcs and constrained optimal control.
-| Parameter | True value | Identified value |
-|------------------|------------|------------------|
-| $k_1$ | 1 | 0.998 |
-| $k_2$ | 10 | 9.97 |
-| $k_3$ | 1 | 1.001 |
+---
### References
-- Jackson, E. A. (1968). *The existence of singular extremals*. Journal of Optimization Theory and Applications.
-- Biegler, L. T. (2010). *Nonlinear Programming: Concepts, Algorithms, and Applications to Chemical Processes*. SIAM.
-- BOCOP repository: https://github.com/control-toolbox/bocop
+- **Jackson, E. A. (1968).** *The existence of singular extremals*. Journal of Optimization Theory and Applications.
+ Discusses the theoretical existence of singular extremals, forming the basis of the Jackson benchmark problem.
+
+- **Biegler, L. T. (2010).** *Nonlinear Programming: Concepts, Algorithms, and Applications to Chemical Processes*. SIAM.
+ Provides background on nonlinear programming methods applicable to problems like Jackson.
+
+- **BOCOP Repository: Jackson Example.** [https://github.com/control-toolbox/bocop/tree/main/bocop](https://github.com/control-toolbox/bocop/tree/main/bocop)
+ Contains the Jackson problem implementation for testing direct transcription and NLP solvers.
diff --git a/ext/Descriptions/moonlander.md b/ext/Descriptions/moonlander.md
index 44941cd8..ffe96fc0 100644
--- a/ext/Descriptions/moonlander.md
+++ b/ext/Descriptions/moonlander.md
@@ -1,78 +1,55 @@
-This problem models the optimal control of a lunar lander aiming to reach a target position on the lunar surface while minimizing fuel consumption. The system is described by a set of differential equations governing the motion of the lander, including its position, velocity, orientation, and angular velocity.
+The **Moonlander optimal control problem** is a benchmark in constrained optimal control.
+It models the powered descent of a spacecraft aiming to land on a specified target on the Moon's surface in minimum time.
+The system comprises six state variables: the horizontal position $p_1(t)$, vertical position $p_2(t)$, horizontal velocity $dp_1(t)$, vertical velocity $dp_2(t)$, orientation angle $\theta(t)$, and angular velocity $d\theta(t)$.
+The control variables $F_1(t)$ and $F_2(t)$ represent the thrust forces applied by the lander's engines.
+The objective is to **minimise the final time** $t_f$ while ensuring a soft landing at the target, subject to control and physical constraints.
-### Problem Formulation
+### Mathematical formulation
-The objective is to minimize the final time $t_f$ subject to the following constraints and dynamics:
-
-- **State Variables**: Position ($p_1$, $p_2$), velocity ($dp_1$, $dp_2$), orientation ($\theta$), and angular velocity ($d\theta$).
-- **Control Inputs**: Thrust forces $F_1$ and $F_2$ applied along the lander's orientation.
-- **Final Time Constraint**: $0.1 \le t_f \le 1.0$.
-- **Control Constraints**: $0 \le F_1(t), F_2(t) \le \text{max\_thrust}$.
-- **Initial Conditions**:
+The problem can be stated as
```math
-p_1(0) = 0, \quad p_2(0) = 0, \quad dp_1(0) = 0, \quad dp_2(0) = 0, \quad \theta(0) = 0, \quad d\theta(0) = 0
+\begin{aligned}
+\min_{p_1,p_2,dp_1,dp_2,\theta,d\theta,F_1,F_2,t_f} \quad & t_f \\
+\text{s.t.} \quad
+\dot{x}(t) &= \begin{bmatrix} dp_1(t) \\ dp_2(t) \\ \ddot{p}_1(t) \\ \ddot{p}_2(t) \\ d\theta(t) \\ \ddot{\theta}(t) \end{bmatrix} = \begin{bmatrix} dp_1(t) \\ dp_2(t) \\ \frac{1}{m} F_{\text{tot},1}(t) \\ \frac{1}{m} F_{\text{tot},2}(t) - g \\ \frac{1}{I} (D/2) (F_2(t) - F_1(t)) \end{bmatrix}, \\
+0 \le F_1(t), F_2(t) &\le 2g, \quad 0.1 \le t_f \le 1.0, \\
+p_1(0) = p_2(0) = dp_1(0) = dp_2(0) = \theta(0) = d\theta(0) &= 0, \\
+p_1(t_f) = 5.0, \; p_2(t_f) = 5.0, \; dp_1(t_f) = dp_2(t_f) &= 0.
+\end{aligned}
```
-- **Final Conditions**:
+Here, $x(t) = (p_1, p_2, dp_1, dp_2, \theta, d\theta)$, and the total thrust in the lander frame $F_{\text{tot}}$ is obtained by rotating the engine thrusts according to the lander's orientation.
-```math
-p_1(t_f) = \text{target}[1], \quad p_2(t_f) = \text{target}[2], \quad dp_1(t_f) = 0, \quad dp_2(t_f) = 0
-```
+### Parameter values
-### System Dynamics
+- Mass: $m = 1$
+- Lunar gravity: $g = 9.81$
+- Moment of inertia: $I = 0.1$
+- Distance between thrusters: $D = 1$
+- Maximum thrust: $F_{\max} = 2g$
+- Target position: $p_{\rm target} = [5.0, 5.0]$
+- Final time bounds: $t_f \in [0.1, 1.0]$
-The dynamics of the lander are governed by
+### Qualitative behaviour
-```math
-\dot{p}_1 = dp_1
-```
+The optimal solution typically exhibits a **bang–bang structure**, where thrusts switch between minimum and maximum values, ensuring minimal landing time while satisfying position and velocity constraints.
+Orientation dynamics $\theta(t)$ and $d\theta(t)$ play a critical role in directing the thrust to achieve the target.
-```math
-\dot{p}_2 = dp_2
-```
+### Characteristics
-```math
-\dot{dp}_1 = \frac{1}{m} F_{\rm tot}[1]
-```
-
-```math
-\dot{dp}_2 = \frac{1}{m} F_{\rm tot}[2] - g
-```
-
-```math
-\dot{\theta} = d\theta
-```
+- Nonlinear coupled dynamics with orientation-dependent thrust,
+- Control constraints and bounded final time ensuring feasibility,
+- Terminal position and velocity constraints for a soft landing,
+- Benchmark problem for testing direct transcription and NLP solvers.
-```math
-\dot{d\theta} = \frac{1}{I} \left( \frac{D}{2} (F_2 - F_1) \right)
-```
-
-where
-
-```math
-F_{\rm tot} =
-\begin{bmatrix}
-\cos\theta & -\sin\theta \\
-\sin\theta & \cos\theta
-\end{bmatrix}
-\begin{bmatrix} 0 \\ F_1 + F_2 \end{bmatrix}_{1:2}
-```
-
-represents the total thrust vector in the inertial frame, $m$ is the mass of the lander, $I$ its moment of inertia, $D$ the distance between thrusters, and $g$ the lunar gravitational acceleration.
-
-### Objective
-
-The goal is to **minimize the final time** $t_f$:
-
-```math
-J = t_f \to \min
-```
+### References
-subject to the dynamics, boundary conditions, and control constraints.
+- Gazzola, F., & Marchini, E. M. (2021). *The moon lander optimal control problem revisited*. Mathematics in Engineering, 3(5), 1–14. [doi:10.3934/mine.2021040](https://doi.org/10.3934/mine.2021040)
+ Provides a detailed analysis of the Moonlander problem, including minimum-time trajectories and safe landing curves.
-### References
+- GPOPS-II Examples: Moon Lander Problem. [GPOPS-II Moon Lander Example](https://www.gpops2.com/Examples/MoonLander.html)
+ Illustrates practical implementation of the Moonlander OCP in a direct transcription framework, useful for benchmarking.
-- Vanroye, L., Sathya, A., De Schutter, J., & Decré, W. (2023). FATROP: A Fast Constrained Optimal Control Problem Solver for Robot Trajectory Optimization and Control. *arXiv preprint arXiv:2303.16746*. Retrieved from https://arxiv.org/pdf/2303.16746
-- Bryson, A. E., & Ho, Y.-C. (1975). *Applied Optimal Control: Optimization, Estimation, and Control*. Hemisphere Publishing.
-- Hull, D. G. (2007). *Fundamentals of Airplane Flight Mechanics* (Chapter on spacecraft and landing). Springer.
+- Vanroye, L., Sathya, A., De Schutter, J., & Decré, W. (2023). *FATROP: A Fast Constrained Optimal Control Problem Solver*. *arXiv preprint arXiv:2303.16746*.
+ Discusses solver strategies applied to constrained trajectory optimization problems, including applications relevant to the Moonlander.
diff --git a/ext/Descriptions/robbins.md b/ext/Descriptions/robbins.md
index d5538191..1750d793 100644
--- a/ext/Descriptions/robbins.md
+++ b/ext/Descriptions/robbins.md
@@ -1,4 +1,4 @@
-This problem was introduced by Robbins (1980) as a benchmark involving a **third-order state constraint**:
+The **Robbins benchmark problem** is a classical optimal control problem introduced by Robbins (1980), involving a **third-order state constraint**:
```math
x_1^{(3)}(t) = u(t), \qquad x_1(t) \ge 0,
@@ -28,17 +28,39 @@ x(0) = (1, -2, 0), \qquad x(T) = (0,0,0),
and horizon $T=10$. In our simulations we use parameters $\alpha=3$, $\beta=0$, $\gamma=0.5$.
+### Parameter values
+
+- Weight on state $x_1$: $\alpha = 3$
+- Weight on squared state: $\beta = 0$
+- Weight on squared control: $\gamma = 0.5$
+- Final time: $T = 10$
+
### Qualitative behaviour
-Robbins (1980) showed that the exact solution has **infinitely many isolated contact points**, where the state constraint $x_1(t)=0$ is active.
-The unconstrained arcs between contacts decrease geometrically, leading to an **accumulation point**, followed by the trivial singular arc $u=0$, $x=0$.
+The exact solution exhibits **infinitely many isolated contact points**, where the state constraint $x_1(t)=0$ is active.
+Between these contact points, the unconstrained arcs decrease geometrically, accumulating at a point, followed by a trivial singular arc where $u=0$ and $x=0$.
+
+Numerically, capturing all contact points is challenging because the unconstrained arcs quickly shrink below the resolution of a given discretisation.
+Simulations reproduce the initial portion of the solution, showing the first few contact points.
+The control $u(t)$ typically exhibits a **bang–bang structure** with possible singular arcs, while $x_2(t)$ and $x_3(t)$ evolve to satisfy the dynamics.
+
+### Characteristics
-Numerically, capturing all contact points is challenging because the unconstrained arcs rapidly shrink below the resolution of a given discretisation.
-Simulations reproduce the initial part of the solution, showing the first few contact points.
+- Third-order state constraint with inequality $x_1(t) \ge 0$,
+- Linear dynamics in first-order representation,
+- Bang–bang control and singular arcs with infinitely many contact points,
+- Serves as a benchmark for direct transcription and NLP solvers handling state-constrained optimal control problems.
### References
- Robbins, H. M. (1980). *Junction phenomena for optimal control with state-variable inequality constraints of third order*. Journal of Optimization Theory and Applications, 31, 85–99.
+ This is the original paper introducing the Robbins problem. It formulates the third-order state-constrained optimal control problem, describes the accumulation of contact points, and provides theoretical analysis of the junction phenomena.
+
- Hermant, A. (2008). *Sur l'algorithme de tir pour les problèmes de commande optimale avec contraintes sur l'état* (PhD thesis, École Polytechnique X).
+ This thesis discusses numerical shooting methods for state-constrained optimal control problems, including the Robbins problem. It provides practical insights into solving problems with multiple contact points and complex singular arcs.
+
- Jacobson, D. H., Lele, M. M., & Speyer, J. L. (1971). *New necessary conditions of optimality for control problems with state-variable inequality constraints*. Journal of Mathematical Analysis and Applications, 35, 255–284.
-- BOCOP repository: https://github.com/control-toolbox/bocop
+ This foundational work introduces necessary conditions for optimality in control problems with state-variable inequality constraints, forming the theoretical basis for analyzing and solving problems like Robbins.
+
+- BOCOP repository: [Robbins problem](https://github.com/control-toolbox/bocop)
+ Contains a practical implementation of the Robbins benchmark in the BOCOP optimal control framework. This resource allows testing and benchmarking direct transcription and NLP solvers on the Robbins problem, reproducing the first few contact points in numerical simulations.
diff --git a/ext/Descriptions/robot.md b/ext/Descriptions/robot.md
index e6ef2300..1a3191f4 100644
--- a/ext/Descriptions/robot.md
+++ b/ext/Descriptions/robot.md
@@ -1,121 +1,66 @@
-This problem models the **time-optimal motion of a robot arm** moving between two points, following the formulation of Mössner-Beigel (PhD thesis, Heidelberg University) and the implementation described by Vanderbei (2001).
-The arm is represented as a rigid bar of total length $L$, pivoting at the origin of a spherical coordinate system.
-
-### System Description
-
-- The arm protrudes a distance $\rho$ from the origin in one direction, and $L - \rho$ in the opposite direction.
-- The orientation of the arm is described by two angles:
- - $\theta$: horizontal angle from the plane, with bounds $|\theta| \leq \pi$
- - $\varphi$: vertical angle, with bounds $0 \leq \varphi \leq \pi$
-
-The state variables are:
-
-```math
-x = (\rho, \dot{\rho}, \theta, \dot{\theta}, \varphi, \dot{\varphi})
-```
-
-The control inputs are:
-
-```math
-u = (u_{\rho}, u_{\theta}, u_{\varphi})
-```
-
-### Dynamics
-
-The continuous-time dynamics are governed by the second-order system:
-
-```math
-L \ddot{\rho} = u_{\rho}
-```
-
-```math
-I_{\theta} \ddot{\theta} = u_{\theta}
-```
-
-```math
-I_{\varphi} \ddot{\varphi} = u_{\varphi}
-```
-
-where the moments of inertia are defined as:
-
-```math
-I_{\theta} = \big( (L - \rho)^3 + \rho^3 \big) \sin^2(\varphi)
-```
-
-```math
-I_{\varphi} = (L - \rho)^3 + \rho^3
-```
-
-In first-order form, with states $(\rho, \dot{\rho}, \theta, \dot{\theta}, \varphi, \dot{\varphi})$, the dynamics become:
-
-```math
-\dot{\rho} = \dot{\rho}, \quad
-\ddot{\rho} = \frac{u_{\rho}}{L}
-```
-
-```math
-\dot{\theta} = \dot{\theta}, \quad
-\ddot{\theta} = \frac{3u_{\theta}}{I_{\theta}}
-```
-
-```math
-\dot{\varphi} = \dot{\varphi}, \quad
-\ddot{\varphi} = \frac{3u_{\varphi}}{I_{\varphi}}
-```
-
-### Constraints
-
-- **State constraints**:
-```math
-0 \leq \rho(t) \leq L
-```
-
-```math
--\pi \leq \theta(t) \leq \pi
-```
-
-```math
-0 \leq \varphi(t) \leq \pi
-```
-
-- **Control constraints**:
-```math
-|u_{\rho}(t)| \leq 1, \quad |u_{\theta}(t)| \leq 1, \quad |u_{\varphi}(t)| \leq 1
-```
-
-- **Initial conditions**:
-```math
-\rho(0) = 4.5, \quad \theta(0) = 0, \quad \varphi(0) = \pi/4
-```
-
-```math
-\dot{\rho}(0) = \dot{\theta}(0) = \dot{\varphi}(0) = 0
-```
-
-- **Final conditions**:
-```math
-\rho(T) = 4.5, \quad \theta(T) = \tfrac{2\pi}{3}, \quad \varphi(T) = \pi/4
-```
-
-```math
-\dot{\rho}(T) = \dot{\theta}(T) = \dot{\varphi}(T) = 0
-```
-
-### Objective
-
-The objective is to **minimise the transfer time** $T$:
-
-```math
-J = T \to \min
-```
-
-### Remarks
-
-This model is simplified, as it ignores the non-inertial nature of the spherical coordinate frame (i.e., Coriolis and centrifugal forces are not included).
-Vanderbei’s implementation eliminates the controls $u$ by substitution, transforming the equations of motion into inequality constraints on accelerations.
+The **robot arm motion problem** is a benchmark in time-optimal control.
+It models a robotic arm moving between two points in space, represented as a rigid bar of total length $L$ pivoting at the origin of a spherical coordinate system.
+The system includes six state variables: the radial position $\rho(t)$ and its velocity $d\rho(t)$, the horizontal angle $\theta(t)$ and its angular velocity $d\theta(t)$, and the vertical angle $\phi(t)$ and its angular velocity $d\phi(t)$.
+The control variables are $u_\rho(t)$, $u_\theta(t)$, and $u_\phi(t)$, representing the forces or torques applied along the radial and angular directions.
+The objective is to **minimise the total transfer time** while satisfying state and control constraints [BOCOP Robot Arm Example](https://github.com/control-toolbox/bocop/tree/main/bocop).
+
+### Mathematical formulation
+
+The problem can be stated as
+
+```math
+\begin{aligned}
+\min_{\rho, d\rho, \theta, d\theta, \phi, d\phi, u_\rho, u_\theta, u_\phi, T} \quad & J = T \\[0.5em]
+\text{s.t.} \quad &
+\dot{\rho}(t) = d\rho(t), \quad \dot{\theta}(t) = d\theta(t), \quad \dot{\phi}(t) = d\phi(t), \\[0.5em]
+& \dot{d\rho}(t) = \frac{u_\rho(t)}{L}, \quad
+\dot{d\theta}(t) = \frac{3 u_\theta(t)}{I_\theta(t)}, \quad
+\dot{d\phi}(t) = \frac{3 u_\phi(t)}{I_\phi(t)}, \\[0.5em]
+& I_\theta(t) = ((L-\rho(t))^3 + \rho(t)^3) \sin^2(\phi(t)), \quad
+I_\phi(t) = (L-\rho(t))^3 + \rho(t)^3, \\[0.5em]
+& 0 \le \rho(t) \le L, \quad -\pi \le \theta(t) \le \pi, \quad 0 \le \phi(t) \le \pi, \\[0.5em]
+& |u_\rho(t)| \le 1, \quad |u_\theta(t)| \le 1, \quad |u_\phi(t)| \le 1, \\[0.5em]
+& \rho(0) = \rho(T) = 4.5, \quad \theta(0) = 0, \quad \theta(T) = \frac{2\pi}{3}, \quad \phi(0) = \phi(T) = \frac{\pi}{4}, \\[0.5em]
+& d\rho(0) = d\rho(T) = 0, \quad d\theta(0) = d\theta(T) = 0, \quad d\phi(0) = d\phi(T) = 0, \\[0.5em]
+& T \ge 0.1.
+\end{aligned}
+```
+
+The horizon is **free**, as the total transfer time $T$ is optimised.
+
+### Parameter values
+
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Total arm length | $L$ | 5 |
+| Initial radial position | $\rho_0$ | 4.5 |
+| Initial vertical angle | $\phi_0$ | $\pi/4$ |
+| Final horizontal angle | $\theta_f$ | $2\pi/3$ |
+| Maximum radial control | $u_\rho^{\max}$ | 1 |
+| Maximum horizontal control | $u_\theta^{\max}$ | 1 |
+| Maximum vertical control | $u_\phi^{\max}$ | 1 |
+| Minimum final time | $T_{\min}$ | 0.1 |
+
+### Qualitative behaviour
+
+The optimal solution exploits the **degrees of freedom of the spherical coordinate system** to minimise the motion time.
+Control inputs $u_\rho$, $u_\theta$, and $u_\phi$ typically follow **bang–bang profiles**, alternating between maximum and minimum values to reach the target efficiently.
+The dynamics are nonlinear due to the dependence of the moments of inertia $I_\theta$ and $I_\phi$ on the radial and angular positions, which affects the acceleration of the arm.
+
+### Characteristics
+
+- Nonlinear, second-order dynamics in spherical coordinates,
+- Free final time with state and control constraints,
+- Bang–bang control structure in the time-optimal solution,
+- Serves as a benchmark for direct transcription and NLP solvers in robotic motion planning.
### References
- Mössner-Beigel, M. (1989). *Thesis*, Heidelberg University.
+ Introduces the time-optimal robotic arm problem and its theoretical formulation.
+
- Vanderbei, R. J. (2001). *Case studies in trajectory optimization: Trains, Planes, and Other Pastimes*.
-- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Argonne National Laboratory. Retrieved from https://www.mcs.anl.gov/~more/cops/cops3.pdf
+ Provides numerical strategies for robotic motion planning and discusses control substitution techniques for time-optimal trajectories.
+
+- BOCOP examples: Robot arm problem. [BOCOP Robot Arm Example](https://github.com/control-toolbox/bocop/tree/main/bocop)
+ Demonstrates the practical implementation using direct transcription and NLP solvers for robotic motion problems.
diff --git a/ext/Descriptions/rocket.md b/ext/Descriptions/rocket.md
index fc6b95a1..0a87d62e 100644
--- a/ext/Descriptions/rocket.md
+++ b/ext/Descriptions/rocket.md
@@ -1,121 +1,94 @@
-This problem models the **maximisation of the final altitude of a vertically launched rocket**, a classical example in optimal control theory with **singular arcs** (see Bryson [7, pp. 392–394]).
-The rocket dynamics involve thrust, drag, gravity, and fuel consumption.
+The **Goddard rocket problem** is a classical optimal control problem that models the ascent of a vertically launched rocket.
+The objective is to **maximise the final altitude** by optimally controlling the thrust, while accounting for atmospheric drag, gravity, and fuel consumption.
+This problem is notable for its inclusion of **singular arcs**, where the control is neither at its maximum nor minimum, complicating the solution process.
+Historically, it was first proposed by **Robert H. Goddard** in 1919 (*A Method of Reaching Extreme Altitudes*), one of the founding works in modern rocketry.
-### System Description
+### Mathematical formulation
-The state variables are:
+The problem can be stated as
```math
-x = (h, v, m)
-```
-
-- $h$: altitude
-- $v$: vertical velocity
-- $m$: mass
-
-The control variable is:
-
-```math
-u = T
-```
-
-- $T$: rocket thrust
-
-### Dynamics
-
-The equations of motion are:
-
-```math
-\dot{h} = v
-```
-
-```math
-\dot{v} = \frac{T - D(h, v) - m g(h)}{m}
-```
-
-```math
-\dot{m} = -\frac{T}{c}
+\begin{aligned}
+\min_{h,v,m,T} \quad & J = -h(T) \\[0.5em]
+\text{s.t.} \quad
+\dot{h}(t) &= v(t), \\[0.25em]
+\dot{v}(t) &= \frac{T(t) - D(h(t), v(t)) - m(t) g(h(t))}{m(t)}, \\[0.25em]
+\dot{m}(t) &= -\frac{T(t)}{c}, \\[0.25em]
+0 \le T(t) &\le T_{\max}, \\[0.25em]
+h(t) &\ge h_0, \\[0.25em]
+v(t) &\ge v_0, \\[0.25em]
+m_f \le m(t) &\le m_0, \\[0.25em]
+h(0) &= h_0, \quad v(0) = v_0, \quad m(0) = m_0, \\[0.25em]
+m(T) &= m_f,
+\end{aligned}
```
where
-- **Drag**:
+- **Drag**:
```math
-D(h, v) = D_c \, v^2 \exp\!\left(-h_c \, \frac{h - h_0}{h_0}\right)
+D(h,v) = D_c \, v^2 \exp\!\left(-h_c \frac{h - h_0}{h_0}\right)
```
-- **Gravity**:
+- **Gravity**:
```math
g(h) = g_0 \left(\frac{h_0}{h}\right)^2
```
-- **Fuel constant**:
-```math
-c = \tfrac{1}{2}\sqrt{g_0 h_0}
-```
-
-### Constraints
-
-- **Control constraints**:
-```math
-0 \leq T(t) \leq T_{\max}
-```
-
-with
+- **Fuel constant**:
```math
-T_{\max} = T_c \, m_0 g_0
+c = \frac{1}{2} \sqrt{g_0 h_0}
```
-- **State constraints**:
+- **Maximum thrust**:
```math
-h(t) \geq h_0
+T_{\max} = T_c m_0 g_0
```
-```math
-v(t) \geq v_0
-```
+### Parameters
-```math
-m_f \leq m(t) \leq m_0
-```
+| Parameter | Symbol | Value |
+|----------------------------|--------|------------------------|
+| Initial altitude | $h_0$ | 1 |
+| Initial velocity | $v_0$ | 0 |
+| Initial mass | $m_0$ | 1 |
+| Gravitational constant | $g_0$ | 1 |
+| Thrust coefficient | $T_c$ | 3.5 |
+| Drag coefficient | $D_c$ | $\frac{1}{2} v_c \frac{m_0}{g_0}$ |
+| Characteristic altitude | $h_c$ | 500 |
+| Characteristic velocity | $v_c$ | 620 |
+| Characteristic mass ratio | $m_c$ | 0.6 |
+| Final mass | $m_f$ | $m_c m_0$ |
-- **Initial conditions**:
-```math
-h(0) = h_0, \quad v(0) = v_0, \quad m(0) = m_0
-```
+### Qualitative behaviour
-- **Final condition**:
-```math
-m(T) = m_f = m_c m_0
-```
+The optimal trajectory often exhibits a **bang–singular–bang structure**, with thrust at maximum, then along a singular arc, and finally at minimum.
+The singular arc arises from the interplay of drag and gravity, producing a non-trivial optimal control law.
+This behaviour illustrates the complexity of the problem and the challenges in numerically resolving singular arcs.
-### Objective
+### Characteristics
-The objective is to **maximise the final altitude** $h(T)$:
+- Nonlinear coupled dynamics with drag and gravity,
+- Singular arcs in the optimal control,
+- State and control constraints,
+- Classical benchmark for testing optimal control algorithms, especially in aerospace applications.
-```math
-J = -h(T) \to \min
-```
+### References
-### Parameters
+- Goddard, R. H. (1919). *A Method of Reaching Extreme Altitudes*. Smithsonian Institution, Washington D.C.
+ The original formulation of the vertical rocket ascent problem, introducing fuel consumption, gravity, and drag effects.
-For the standard nondimensionalised version of the problem:
+- Bryson, A. E. (1999). *Dynamic Optimization*. Addison Wesley Longman. pp. 392–394.
+ Provides an introduction to singular arcs in optimal control with the Goddard rocket as a primary example.
-```math
-h_0 = 1, \quad v_0 = 0, \quad m_0 = 1, \quad g_0 = 1
-```
+- Garfinkel, B. (1963). A Solution of the Goddard Problem. *SIAM Journal on Control*, 1(1), 1–20. [doi:10.1137/0301020](https://epubs.siam.org/doi/10.1137/0301020)
+ Analytical and numerical solutions of the Goddard problem highlighting singular arc structures.
-```math
-T_{\max} = 3.5 g_0 m_0, \quad D_c = \tfrac{1}{2} v_c \tfrac{m_0}{g_0}, \quad c = \tfrac{1}{2}\sqrt{g_0 h_0}
-```
+- Tsiotras, P., & Kelley, H. J. (1992). Goddard Problem with Constrained Time of Flight. *Journal of Guidance, Control, and Dynamics*, 15(2), 394–399. [doi:10.2514/3.20836](https://doi.org/10.2514/3.20836)
+ Study of time-constrained versions of the Goddard problem and associated optimal control strategies.
-with
-
-```math
-h_c = 500, \quad v_c = 620, \quad m_c = 0.6
-```
-
-### References
+- Bonnans, F., Martinon, P., & Trélat, E. (2007). Singular Arcs in the Generalized Goddard's Problem. [arXiv:math/0703911](https://arxiv.org/pdf/math/0703911)
+ Analysis of singular arcs and numerical methods for the Goddard problem.
-- Bryson, A. E. (1999). *Dynamic Optimization*. Addison Wesley Longman. (pp. 392–394)
-- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Argonne National Laboratory. Retrieved from https://www.mcs.anl.gov/~more/cops/cops3.pdf
+- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Argonne National Laboratory. [Link](https://www.mcs.anl.gov/~more/cops/cops3.pdf)
+ Provides a benchmark implementation of the Goddard rocket problem for testing numerical optimization methods.
diff --git a/ext/Descriptions/space_shuttle.md b/ext/Descriptions/space_shuttle.md
index adcb4176..8837af3d 100644
--- a/ext/Descriptions/space_shuttle.md
+++ b/ext/Descriptions/space_shuttle.md
@@ -1,107 +1,80 @@
-The *space shuttle reentry problem* is a classical benchmark in aerospace optimal control, originating from reentry trajectory studies (see Betts 2010, Bulirsch 1971, Dickmanns 1972). It describes the atmospheric descent of the space shuttle from high altitude to the Terminal Area Energy Management (TAEM) interface. The aim is to maximise the crossrange, i.e. the final latitude at TAEM, subject to nonlinear dynamics, control bounds, and path constraints.
+The **Space Shuttle reentry problem** is a classical benchmark in aerospace optimal control.
+It models the atmospheric descent of the space shuttle from high altitude to the Terminal Area Energy Management (TAEM) interface.
+The system includes six state variables: altitude $h(t)$, longitude $\phi(t)$, latitude $\theta(t)$, velocity $v(t)$, flight path angle $\gamma(t)$, and azimuth $\psi(t)$.
+The control variables are the angle of attack $\alpha(t)$ and bank angle $\beta(t)$.
+The goal is to **maximise the final latitude (crossrange) at TAEM** while satisfying aerodynamic, gravitational, and operational constraints [Betts 2010; Bulirsch 1971; Dickmanns 1972].
-The problem can be written as
+### Mathematical formulation
-```math
-\max_{x,\,u} J(x,u) = \theta(t_f),
-```
-
-or equivalently
+The problem can be stated as
-```math
-\min_{x,\,u} J(x,u) = -\theta(t_f),
-```
-
-subject to the dynamics
```math
\begin{aligned}
-\dot{h}(t) &= v(t)\,\sin \gamma(t), \\
-\dot{\phi}(t) &= \tfrac{v(t)}{r(t)} \cos \gamma(t) \sin \psi(t) / \cos \theta(t), \\
-\dot{\theta}(t) &= \tfrac{v(t)}{r(t)} \cos \gamma(t) \cos \psi(t), \\
-\dot{v}(t) &= -\tfrac{D}{m} - g(t)\,\sin \gamma(t), \\
-\dot{\gamma}(t) &= \tfrac{L}{m v(t)} \cos \beta(t)
- + \cos \gamma(t)\Big(\tfrac{v(t)}{r(t)} - \tfrac{g(t)}{v(t)}\Big), \\
-\dot{\psi}(t) &= \tfrac{L}{m v(t)\cos \gamma(t)} \sin \beta(t)
- + \tfrac{v(t)}{r(t)\cos \theta(t)} \cos \gamma(t)\sin \psi(t)\sin \theta(t),
+\min_{h, \phi, \theta, v, \gamma, \psi, \alpha, \beta} \quad & J(h, \theta) = - \theta(T) \\[1em]
+\text{s.t.} \quad &
+\dot{h}(t) = v \sin(\gamma), \\[0.5em]
+& \dot{\phi}(t) = \frac{v}{r} \cos(\gamma) \frac{\sin(\psi)}{\cos(\theta)}, \\[0.5em]
+& \dot{\theta}(t) = \frac{v}{r} \cos(\gamma) \cos(\psi), \\[0.5em]
+& \dot{v}(t) = -\frac{D(h,v,\alpha)}{m} - g(h) \sin(\gamma), \\[0.5em]
+& \dot{\gamma}(t) = \frac{L(h,v,\alpha)}{m v} \cos(\beta) + \cos(\gamma) \left( \frac{v}{r} - \frac{g(h)}{v} \right), \\[0.5em]
+& \dot{\psi}(t) = \frac{L(h,v,\alpha)}{m v \cos(\gamma)} \sin(\beta) + \frac{v}{r \cos(\theta)} \cos(\gamma) \sin(\psi) \sin(\theta), \\[0.5em]
+& \alpha_{\min} \le \alpha(t) \le \alpha_{\max}, \\[0.5em]
+& \beta_{\min} \le \beta(t) \le \beta_{\max}, \\[0.5em]
+& h_{\min} \le h(t) \le h_{\max}, \quad
+v_{\min} \le v(t) \le v_{\max}, \\[0.5em]
+& \gamma_{\min} \le \gamma(t) \le \gamma_{\max}, \quad
+\theta_{\min} \le \theta(t) \le \theta_{\max}, \\[0.5em]
+& h(0) = h_s, \; \phi(0) = \phi_s, \; \theta(0) = \theta_s, \\[0.5em]
+& v(0) = v_s, \; \gamma(0) = \gamma_s, \; \psi(0) = \psi_s, \\[0.5em]
+& h(T) = h_t, \; v(T) = v_t, \; \gamma(T) = \gamma_t.
\end{aligned}
```
-with aerodynamic lift and drag
-```math
-D = \tfrac{1}{2} c_D S \rho v^2, \qquad
-L = \tfrac{1}{2} c_L S \rho v^2,
-```
+The final time $T$ is free but bounded.
-where the coefficients are given by
-```math
-c_D = b_0 + b_1 \alpha^\circ + b_2 (\alpha^\circ)^2, \qquad
-c_L = a_0 + a_1 \alpha^\circ, \qquad
-\alpha^\circ = \tfrac{180}{\pi} \alpha,
-```
-
-and
-```math
-\rho = \rho_0 e^{-h/h_r}, \qquad
-r = R_e + h, \qquad
-g = \tfrac{\mu}{r^2}.
-```
+### Parameters
-### Boundary conditions
-
-At reentry interface ($t=0$):
-```math
-h(0) = 260{,}000 \ \text{ft}, \quad
-v(0) = 25{,}600 \ \text{ft/s}, \quad
-\phi(0)=0, \quad \theta(0)=0, \quad
-\gamma(0)=-1^\circ, \quad \psi(0)=90^\circ.
-```
-
-At TAEM interface ($t=t_f$):
-```math
-h(t_f) = 80{,}000 \ \text{ft}, \quad
-v(t_f) = 2{,}500 \ \text{ft/s}, \quad
-\gamma(t_f) = -5^\circ.
-```
-
-Final time is free within
-```math
-500\Delta t_{\min} \le t_f \le 500\Delta t_{\max}, \qquad
-\Delta t_{\min}=3.5, \quad \Delta t_{\max}=4.5.
-```
-
-### Constraints
-
-- State bounds:
-```math
-h(t) \ge 0, \qquad
--89^\circ \le \theta(t) \le 89^\circ, \qquad
-v(t) \ge 0, \qquad
--89^\circ \le \gamma(t) \le 89^\circ.
-```
-
-- Control bounds:
-```math
--90^\circ \le \alpha(t) \le 90^\circ, \qquad
--89^\circ \le \beta(t) \le 1^\circ.
-```
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Initial altitude | $h_s$ | $2.6 \times 10^5$ ft |
+| Initial velocity | $v_s$ | $2.56 \times 10^4$ ft/s |
+| Initial flight path angle | $\gamma_s$ | $-1^\circ$ |
+| Initial azimuth | $\psi_s$ | $90^\circ$ |
+| Final altitude | $h_t$ | $0.8 \times 10^5$ ft |
+| Final velocity | $v_t$ | $0.25 \times 10^4$ ft/s |
+| Final flight path angle | $\gamma_t$ | $-5^\circ$ |
+| Mass | $m$ | $w/g_0$ |
+| Reference area | $S$ | 2690 |
+| Earth's radius | $R_e$ | 20902900 ft |
+| Gravitational parameter | $\mu$ | $0.14076539 \cdot 10^{17}$ |
### Qualitative behaviour
- The optimal trajectory balances **lift** and **drag** to control heating and deceleration while extending crossrange.
- The **bank angle** $\beta$ determines heading changes and crossrange capability.
-- The solution typically includes a combination of steep reentry to dissipate energy and crossrange manoeuvres to reach the target latitude.
+- The solution typically combines **steep reentry** to dissipate energy and **crossrange manoeuvres** to reach the target latitude.
### Characteristics
- Nonlinear, six–state dynamics with two controls.
- Strongly nonlinear aerodynamic coefficients.
- Free final time with bounded range.
-- Path constraints on both states and controls.
-- Widely used as a benchmark in optimal control and trajectory optimisation.
+- Path constraints on states and controls.
+- Widely used as a benchmark in optimal control and trajectory optimisation.
+
+### References and relevance
+
+- **Betts, J.T. (2010)**. *Practical Methods for Optimal Control and Estimation Using Nonlinear Programming*. SIAM.
+ Provides a comprehensive discussion of trajectory optimisation for aerospace vehicles, including space shuttle reentry problems. It details numerical methods applicable to nonlinear, constrained, free-final-time problems.
+
+- **Bulirsch, R. (1971)**. Numerical solution of optimal control problems with state constraints by direct methods. *Numerische Mathematik*.
+ Introduces early direct methods for optimal control problems with constraints, which are foundational for solving shuttle reentry trajectory problems with bounded states and controls.
+
+- **Dickmanns, E.D. (1972)**. Numerical solution methods for nonlinear optimal control problems with state constraints. *Automatica*.
+ Focuses on numerical techniques for nonlinear optimal control with state constraints, directly relevant to handling the shuttle’s aerodynamic and flight path restrictions.
-### References
+- **Ascher, U.M., Mattheij, R.M.M., & Russell, R.D. (1988)**. *Numerical Solution of Boundary Value Problems for Ordinary Differential Equations*.
+ Provides practical algorithms for solving boundary value problems, which underpin the solution of two-point boundary value problems like the shuttle reentry trajectory with fixed initial and terminal states.
-- Betts, J.T. (2010). *Practical Methods for Optimal Control and Estimation Using Nonlinear Programming*. SIAM.
-- Bulirsch, R. (1971). Numerical solution of optimal control problems with state constraints by direct methods. *Numerische Mathematik*.
-- Dickmanns, E.D. (1972). Numerical solution methods for nonlinear optimal control problems with state constraints. *Automatica*.
-- Ascher, U.M., Mattheij, R.M.M., & Russell, R.D. (1988). *Numerical Solution of Boundary Value Problems for Ordinary Differential Equations*.
+- **Betts, J.T. (2001)**. Survey of numerical methods for trajectory optimization. *Journal of Guidance, Control, and Dynamics*, 24(4), 643–653.
+ Reviews trajectory optimisation methods including direct transcription and collocation, which are commonly applied to the space shuttle reentry benchmark.
diff --git a/ext/Descriptions/steering.md b/ext/Descriptions/steering.md
index a1393ec9..f8903e11 100644
--- a/ext/Descriptions/steering.md
+++ b/ext/Descriptions/steering.md
@@ -1,71 +1,66 @@
-This problem models the **time-optimal steering of a vehicle** with bounded control inputs, inspired by the COPS collection (More et al., 2001).
-The goal is to move the system from an initial state to a specified final state in minimum time while respecting constraints on the steering input.
+The **particle steering problem** is a classical benchmark in time-optimal control.
+It models the task of steering a particle from a given initial state to a specified terminal state while minimising travel time.
+The system includes four state variables representing position and velocity components, and a single scalar control input that determines the steering direction.
+Control bounds and dynamics constraints must be satisfied throughout the trajectory [Athans & Falb 1966; Bryson & Ho 1975].
-### System Dynamics
+### Mathematical formulation
-The system has four states and one control:
-
-- $x_1$ : horizontal position
-- $x_2$ : vertical position
-- $x_3$ : velocity in $x_1$ direction
-- $x_4$ : velocity in $x_2$ direction
-- $u$ : steering angle (control)
-
-The dynamics are expressed as:
+The problem can be stated as
```math
-\dot{x}_1 = x_3
+\begin{aligned}
+\min_{x,u} \quad & J(x,u) = t_f \\[0.5em]
+\text{s.t.} \quad &
+\dot{x}_1(t) = x_3(t), \\[0.5em]
+& \dot{x}_2(t) = x_4(t), \\[0.5em]
+& \dot{x}_3(t) = a \cos(u(t)), \\[0.5em]
+& \dot{x}_4(t) = a \sin(u(t)), \\[0.5em]
+& u_{\min} \le u(t) \le u_{\max}, \\[0.5em]
+& x(0) = x_s, \quad x(t_f) = x_f, \\[0.5em]
+& t_f \ge 0.
+\end{aligned}
```
-```math
-\dot{x}_2 = x_4
-```
+where
-```math
-\dot{x}_3 = a \cos(u)
-```
+- the state $x = (x_1, x_2, x_3, x_4)$ are the particle states (position and velocity components),
+- the control $u$ is the steering control, bounded by $u_{\min} = -\pi/2$ and $u_{\max} = \pi/2$,
+- the parameter $a$ is a constant scaling the acceleration.
-```math
-\dot{x}_4 = a \sin(u)
-```
-
-where $a$ is a constant acceleration parameter, and $u$ is constrained by:
+### Parameter values
-```math
-u_{\min} \le u(t) \le u_{\max}
-```
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Acceleration magnitude | $a$ | 100 |
+| Control lower bound | $u_{\min}$ | $-\pi/2$ |
+| Control upper bound | $u_{\max}$ | $\pi/2$ |
+| Initial state | $x_s$ | $(0,0,0,0)$ |
+| Final state | $x_f$ | $(\text{unspecified},5,45,0)$ |
-### Boundary Conditions
+### Qualitative behaviour
-- **Initial conditions**:
+The optimal trajectory typically exhibits **bang–bang control**, where the steering input switches between its extreme values to minimise time.
+The trajectory balances position and velocity to reach the terminal state in minimal time.
+Analytical solutions exist for simplified cases, while numerical methods are required for general conditions.
-```math
-x(0) = x_s
-```
+### Characteristics
-- **Final conditions**:
+- Four-dimensional state with single scalar control,
+- Time-optimal objective,
+- Control constraints and boundary conditions,
+- Bang–bang control structure with possible singular arcs,
+- Serves as a benchmark for time-optimal control methods.
-```math
-x(T) = x_f, \quad T \ge 0
-```
-
-- **Control constraints**:
-
-```math
-u_{\min} \le u(t) \le u_{\max}
-```
-
-### Objective
-
-The goal is to **minimize the final time** $T$:
+### References
-```math
-J = T \to \min
-```
+- Athans, M., & Falb, P.L. (1966). *Optimal Control: An Introduction to the Theory and Its Applications*. McGraw-Hill.
+ This classic text introduces the theory of optimal control, including time-optimal problems and bang–bang solutions, providing foundational concepts relevant to particle steering problems.
-subject to the dynamics, boundary conditions, and control constraints.
+- Bryson, A.E., & Ho, Y.-C. (1975). *Applied Optimal Control: Optimization, Estimation, and Control*. Hemisphere Publishing.
+ The book presents methods and examples of time-optimal and energy-optimal control problems, including systems with bang–bang and singular arc solutions, directly relevant to the steering problem.
-### References
+- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Mathematics and Computer Science Division, Argonne National Laboratory. [https://www.mcs.anl.gov/~more/cops/cops3.pdf](https://www.mcs.anl.gov/~more/cops/cops3.pdf)
+ This report lists the particle steering problem as a benchmark for nonlinear constrained optimization, providing the original formulation and serving as a reference for numerical experimentation with optimal control solvers.
-- More, J., Garbow, B., Hillstrom, K., & Watson, L. (2001). *COPS: Constrained Optimization Problem Set* (COPS3). Mathematics and Computer Science Division, Argonne National Laboratory. Retrieved from https://www.mcs.anl.gov/~more/cops/cops3.pdf
-- Cesari, L. (1983). *Optimization – Theory and Applications. Problems with Ordinary Differential Equations*. Springer-Verlag.
+- PSOPT Example Set. *Particle Steering Problem*. [https://www.psopt.net/list-of-examples?utm_source=chatgpt.com](https://www.psopt.net/list-of-examples?utm_source=chatgpt.com)
+ A modern implementation of the particle steering problem for testing direct optimal control methods, illustrating numerical approaches for time-optimal trajectory problems.
diff --git a/ext/Descriptions/truck_trailer.md b/ext/Descriptions/truck_trailer.md
deleted file mode 100644
index e564fe21..00000000
--- a/ext/Descriptions/truck_trailer.md
+++ /dev/null
@@ -1,86 +0,0 @@
-This problem models the **minimum-time maneuvering of a truck with two trailers** while respecting steering, velocity, and articulation constraints.
-The goal is to move the vehicle from an initial configuration to a target position and orientation while minimizing final time and reducing excessive trailer articulation.
-
-### System Description
-
-The system has **7 states** and **2 controls**:
-
-- **States**:
- - $x_2, y_2$: position of the rear trailer axle
- - $\theta_0, \theta_1, \theta_2$: orientations of the truck and two trailers
- - $v_0$: longitudinal velocity of the truck
- - $\delta_0$: steering angle of the truck
-
-- **Controls**:
- - $\dot{v}_0$: acceleration of the truck
- - $\dot{\delta}_0$: steering rate
-
-- **Auxiliary variables**:
- - $\beta_{01} = \theta_0 - \theta_1$
- - $\beta_{12} = \theta_1 - \theta_2$
-
-### Constraints
-
-- **Time horizon**: $1 \le t_f \le 1000$
-- **State constraints**:
-
-```math
--\pi/2 \le \theta_0(t), \theta_1(t) \le \pi/2, \quad
--\frac{0.2 v_{\rm max}}{1} \le v_0(t) \le \frac{0.2 v_{\rm max}}{1}, \quad
--\pi/6 \le \delta_0(t) \le \pi/6
-```
-
-- **Control constraints**:
-
-```math
--1 \le \dot{v}_0(t) \le 1, \quad
--\pi/10 \le \dot{\delta}_0(t) \le \pi/10
-```
-
-- **Path constraints**:
-
-```math
--\pi/2 \le \beta_{01}(t), \beta_{12}(t) \le \pi/2
-```
-
-- **Boundary conditions**:
-
-```math
-x_2(0) = x_{2,0}, \quad y_2(0) = y_{2,0}, \quad \theta_0(0) = \theta_{0,0}, \quad \theta_1(0) = \theta_{1,0}, \quad \theta_2(0) = \theta_{2,0}
-```
-
-```math
-x_2(t_f) = x_{2,f}, \quad y_2(t_f) = y_{2,f}, \quad \theta_2(t_f) = \theta_{2,f}, \quad \beta_{01}(t_f) = \theta_{0,f} - \theta_{1,f}, \quad \beta_{12}(t_f) = \theta_{1,f} - \theta_{2,f}
-```
-
-### Dynamics
-
-The truck-trailer kinematics are governed by
-
-```math
-\begin{aligned}
-\dot{\theta}_0 &= \frac{v_0}{L_0} \tan\delta_0, \\
-\dot{\theta}_1 &= \frac{v_0}{L_1} \sin\beta_{01} - \frac{M_0}{L_1} \cos\beta_{01} \, \dot{\theta}_0, \\
-v_1 &= v_0 \cos\beta_{01} + M_0 \sin\beta_{01} \, \dot{\theta}_0, \\
-\dot{\theta}_2 &= \frac{v_1}{L_2} \sin\beta_{12} - \frac{M_1}{L_2} \cos\beta_{12} \, \dot{\theta}_1, \\
-v_2 &= v_1 \cos\beta_{12} + M_1 \sin\beta_{12} \, \dot{\theta}_1, \\
-\dot{x}_2 &= v_2 \cos\theta_2, \quad \dot{y}_2 = v_2 \sin\theta_2, \\
-\dot{v}_0 &= \dot{v}_0^{\rm control}, \quad \dot{\delta}_0 = \dot{\delta}_0^{\rm control}
-\end{aligned}
-```
-
-where $L_i$ and $M_i$ are the vehicle and trailer geometric parameters.
-
-### Objective
-
-The goal is to **minimize the final time** while reducing large trailer articulation angles:
-
-```math
-J = t_f + \int_0^{t_f} (\beta_{01}^2(t) + \beta_{12}^2(t)) \, dt \to \min
-```
-
-### References
-
-- Vanroye, L., Sathya, A., De Schutter, J., & Decré, W. (2023). FATROP: A Fast Constrained Optimal Control Problem Solver for Robot Trajectory Optimization and Control. *arXiv preprint arXiv:2303.16746*. Retrieved from https://arxiv.org/pdf/2303.16746
-- Kretzschmar, H., & Burgard, W. (2019). Optimal motion planning for truck and trailer systems: A review. *IEEE Transactions on Intelligent Vehicles*, 4(3), 256–271.
-- Falcone, P., Borrelli, F., Asgari, J., Tseng, H. E., & Hrovat, D. (2007). Predictive active steering control for autonomous vehicle systems. *IEEE Transactions on Control Systems Technology*, 15(3), 566–580.
\ No newline at end of file
diff --git a/ext/Descriptions/vanderpol.md b/ext/Descriptions/vanderpol.md
index 47c2ab8e..885a0c96 100644
--- a/ext/Descriptions/vanderpol.md
+++ b/ext/Descriptions/vanderpol.md
@@ -1,40 +1,55 @@
-This problem considers the **controlled Van der Pol oscillator**, a classical nonlinear system with non-linear damping:
+The **Van der Pol oscillator control problem** is a benchmark in nonlinear optimal control.
+It models a two-dimensional oscillator with nonlinear damping and a control input.
+The state vector is $x(t) = [x_1(t), x_2(t)]^\top$, and the control variable is $u(t)$.
+The goal is to minimise a quadratic cost functional over a **fixed horizon** $T = 2$ while steering the system from a given initial condition.
-```math
-\ddot{x}(t) - \varepsilon \, \omega \, (1 - x(t)^2) \dot{x}(t) + \omega^2 x(t) = u(t),
-```
+### Mathematical formulation
-with control $u(t)$ and cost functional
-
-```math
-J(x,u) = \frac{1}{2} \int_0^{t_f} \big( x_1(t)^2 + x_2(t)^2 + u(t)^2 \big) \, dt.
-```
-
-The system is represented in first-order form as
+The problem can be stated as
```math
\begin{aligned}
-\dot x_1(t) &= x_2(t), \\
-\dot x_2(t) &= \varepsilon \, \omega \, (1 - x_1(t)^2) x_2(t) - \omega^2 x_1(t) + u(t),
+\min_{x,u} \quad & J(x,u) = \frac{1}{2} \int_0^2 \big( x_1(t)^2 + x_2(t)^2 + u(t)^2 \big) \, dt \\[0.5em]
+\text{s.t.} \quad &
+\dot{x}_1(t) = x_2(t), \\[0.5em]
+& \dot{x}_2(t) = \varepsilon \, \omega \, (1 - x_1(t)^2) x_2(t) - \omega^2 x_1(t) + u(t), \\[0.5em]
+& x(0) = [1, 0], \\[0.5em]
+& u(t) \in \mathbb{R}.
\end{aligned}
```
-with boundary condition
-
-```math
-x(0) = (1,0),
-```
+### Parameters
-and horizon $t_f = 2$. In our simulations we use parameters $\varepsilon=1$, $\omega=1$.
+| Parameter | Symbol | Value |
+|-----------|--------|-------|
+| Oscillator frequency | $\omega$ | 1 |
+| Nonlinearity coefficient | $\varepsilon$ | 1 |
+| Final time | $T$ | 2 |
### Qualitative behaviour
-The Van der Pol oscillator exhibits **limit cycle behaviour** for $\varepsilon>0$, where trajectories converge to a stable periodic orbit.
-The non-linear damping leads to relaxation oscillations: slow accumulation along one branch followed by a rapid transition.
-Control input $u(t)$ allows influencing the amplitude and phase of the oscillations, and the quadratic cost penalises deviations from the origin as well as control effort.
+- The optimal control $u(t)$ regulates the nonlinear oscillations to minimise both state excursions and control effort.
+- Due to the quadratic cost, the control is typically **smooth**, without bang–bang or singular arcs.
+- Serves as a standard test for **direct transcription methods** and numerical solvers for nonlinear optimal control problems.
+
+### Characteristics
+
+- Two-dimensional nonlinear dynamics with cubic damping,
+- Quadratic cost on both states and control,
+- Fixed final time,
+- Unconstrained control (extensions may include bounds),
+- Useful benchmark for testing optimal control algorithms and integration schemes.
### References
-- BOCOP repository: [https://github.com/control-toolbox/bocop](https://github.com/control-toolbox/bocop)
-- Wikipedia: [Van der Pol oscillator](https://en.wikipedia.org/wiki/Van_der_Pol_oscillator)
-- Heitmann, S., Breakspear, M., *BOCOP: User Guide*, 2017–2022.
+- **Dymos Documentation**. [*The Van der Pol Oscillator*.](https://openmdao.github.io/dymos/examples/vanderpol/vanderpol.html).
+ Provides an implementation of the Van der Pol oscillator in an optimal control framework, illustrating problem setup, discretisation, and solution.
+
+- **mintOC**. [*Van der Pol Oscillator*](https://mintoc.de/index.php/Van_der_Pol_Oscillator).
+ Presents the Van der Pol optimal control problem with numerical results, highlighting cost function and control profiles.
+
+- James, E.M. (1974). *Time Optimal Control and the Van der Pol Oscillator*. *IMA Journal of Applied Mathematics*, 13(1), 67–78.
+ Discusses theoretical aspects of time-optimal control for the Van der Pol oscillator, including analysis of control strategies.
+
+- Van Dooren, R. (1987). *Numerical Study of the Controlled Van der Pol Oscillator in Optimal Control*. *Numerische Mathematik*, 51(4), 471–485.
+ Explores numerical methods for solving the Van der Pol optimal control problem, providing insights into solution behavior.
\ No newline at end of file
diff --git a/ext/JuMPModels/bioreactor.jl b/ext/JuMPModels/bioreactor.jl
index 363754f5..68613873 100644
--- a/ext/JuMPModels/bioreactor.jl
+++ b/ext/JuMPModels/bioreactor.jl
@@ -37,7 +37,7 @@ function OptimalControlProblems.bioreactor(
# parameters
β = 1
c = 2
- gamma = 1
+ γ = 1
halfperiod = 5
Ks = 0.05
μ2m = 0.1
@@ -88,7 +88,7 @@ function OptimalControlProblems.bioreactor(
# dynamics
dy[k = 0:N], μ[k] * y[k] / (1 + y[k]) - (r + u[k]) * y[k]
- ds[k = 0:N], -μ2[k] * b[k] + u[k] * β * (gamma * y[k] - s[k])
+ ds[k = 0:N], -μ2[k] * b[k] + u[k] * β * (γ * y[k] - s[k])
db[k = 0:N], (μ2[k] - u[k] * β) * b[k]
# objective
diff --git a/ext/JuMPModels/ducted_fan.jl b/ext/JuMPModels/ducted_fan.jl
index 9fbc16c1..3208486d 100644
--- a/ext/JuMPModels/ducted_fan.jl
+++ b/ext/JuMPModels/ducted_fan.jl
@@ -51,7 +51,7 @@ function OptimalControlProblems.ducted_fan(
@variable(model, vα[0:N], start = 0.1)
@variable(model, -5 <= u₁[0:N] <= 5, start = 0.1) # [N]
@variable(model, 0 <= u₂[0:N] <= 17, start = 1) # [N]
- @variable(model, 0.1 <= tf, start = 1)
+ @variable(model, 0.1 <= tf, start = 1.5)
# Boundary constraints
@constraints(
diff --git a/ext/JuMPModels/truck_trailer.jl b/ext/JuMPModels/truck_trailer.jl
deleted file mode 100644
index a9bce34f..00000000
--- a/ext/JuMPModels/truck_trailer.jl
+++ /dev/null
@@ -1,164 +0,0 @@
-"""
-$(TYPEDSIGNATURES)
-
-Constructs and returns a JuMP model for the **Truck Trailer Problem**.
-The model represents the dynamics of a truck with two trailers, including the truck's velocity, steering angles, and articulation angles between the truck and trailers.
-The objective is to minimise the total time taken to park the truck and trailers aligned vertically at a specified target location while maintaining the physical constraints and vehicle kinematics.
-
-# Arguments
-
-- `::JuMPBackend`: Specifies the backend for building the JuMP model.
-- `N::Int=200`: (Keyword) Number of discretisation steps for the time horizon.
-
-# Returns
-
-- `model::JuMP.Model`: A JuMP model representing the truck trailer optimal control problem.
-
-# Example
-
-```julia-repl
-julia> using OptimalControlProblems
-julia> using JuMP
-
-julia> model = OptimalControlProblems.truck_trailer(JuMPBackend(); N=100)
-```
-
-# References
-
-- Problem formulation available at: https://arxiv.org/pdf/2303.16746
-"""
-function OptimalControlProblems.truck_trailer(
- ::JuMPBackend, args...; N::Int=steps_number_data(:truck_trailer), kwargs...
-)
-
- # parameters
- data=[0.4 0.1 0.2; 1.1 0.2 0.2; 0.8 0.1 0.2]
- L0 = data[1, 1]
- M0 = data[1, 2]
- W0 = data[1, 3]
- L1 = data[2, 1]
- M1 = data[2, 2]
- W1 = data[2, 3]
- L2 = data[3, 1]
- M2 = data[3, 2]
- W2 = data[3, 3]
- speedf = 1
- x2_t0 = 0
- y2_t0 = 0
- θ2_t0 = 0
- θ1_t0 = 0
- θ0_t0 = 0
- x2_tf = 0
- y2_tf = -2
- θ2_tf = π / 2
- θ1_tf = π / 2
- θ0_tf = π / 2
-
- # model
- model = JuMP.Model(args...; kwargs...)
-
- # state, control, variable (final time) and initial guess
- @variables(
- model,
- begin
-
- # Final time
- 1 <= tf <= 1000, (start = 10)
-
- # State variables
- x2[0:N], (start = 0.1)
- y2[0:N], (start = 0.1)
- -π / 2 <= θ0[0:N] <= π / 2, (start = 0.1)
- -π / 2 <= θ1[0:N] <= π / 2, (start = 0.1)
- θ2[0:N], (start = 0.1)
- -0.2 * speedf <= v0[0:N] <= 0.2 * speedf, (start = 0.1)
- -π / 6 <= δ0[0:N] <= π / 6, (start = 0.1)
-
- # Control variables
- -1 <= dv0[0:N] <= 1, (start = 0.1)
- -π / 10 <= dδ0[0:N] <= π / 10, (start = 0.1)
- end
- )
-
- # # positions
- # @expressions(
- # model,
- # begin
- # x1[i=0:N], x2[i] + L2 * cos(θ2[i]) + M1 * cos(θ1[i])
- # y1[i=0:N], y2[i] + L2 * sin(θ2[i]) + M1 * sin(θ1[i])
- # x0[i=0:N], x1[i] + L1 * cos(θ1[i]) + M0 * cos(θ0[i])
- # y0[i=0:N], y1[i] + L1 * sin(θ1[i]) + M0 * sin(θ0[i])
- # end
- # )
-
- # intermediate variables
- @expressions(
- model,
- begin
- β01[i = 0:N], θ0[i] - θ1[i]
- β12[i = 0:N], θ1[i] - θ2[i]
- step, tf / N
- end
- )
-
- @constraints(
- model,
- begin
- β01_con[i = 0:N], -π / 2 <= β01[i] <= π / 2
- β12_con[i = 0:N], -π / 2 <= β12[i] <= π / 2
- end
- )
-
- # boundary conditions
- @constraints(
- model,
- begin
-
- # initial constraints
- x2[0] == x2_t0
- y2[0] == y2_t0
- θ0[0] == θ0_t0
- θ1[0] == θ1_t0
- θ2[0] == θ2_t0
-
- # final constraint
- x2[N] == x2_tf
- y2[N] == y2_tf
- θ2[N] == θ2_tf
- β01[N] == θ0_tf - θ1_tf
- β12[N] == θ1_tf - θ2_tf
- end
- )
-
- @expression(model, dθ0[i = 0:N], v0[i] / L0 * tan(δ0[i]))
- @expression(
- model, dθ1[i = 0:N], v0[i] / L1 * sin(β01[i]) - M0 / L1 * cos(β01[i]) * dθ0[i]
- )
- @expression(model, v1[i = 0:N], v0[i] * cos(β01[i]) + M0 * sin(β01[i]) * dθ0[i])
- @expression(
- model, dθ2[i = 0:N], v1[i] / L2 * sin(β12[i]) - M1 / L2 * cos(β12[i]) * dθ1[i]
- )
- @expression(model, v2[i = 0:N], v1[i] * cos(β12[i]) + M1 * sin(β12[i]) * dθ1[i])
- @expression(model, dx2[i = 0:N], v2[i] * cos(θ2[i]))
- @expression(model, dy2[i = 0:N], v2[i] * sin(θ2[i]))
-
- # Dynamics
- @constraints(
- model,
- begin
- ∂x2[i = 1:N], x2[i] == x2[i - 1] + 0.5 * step * (dx2[i] + dx2[i - 1])
- ∂y2[i = 1:N], y2[i] == y2[i - 1] + 0.5 * step * (dy2[i] + dy2[i - 1])
- ∂θ0[i = 1:N], θ0[i] == θ0[i - 1] + 0.5 * step * (dθ0[i] + dθ0[i - 1])
- ∂θ1[i = 1:N], θ1[i] == θ1[i - 1] + 0.5 * step * (dθ1[i] + dθ1[i - 1])
- ∂θ2[i = 1:N], θ2[i] == θ2[i - 1] + 0.5 * step * (dθ2[i] + dθ2[i - 1])
- ∂v0[i = 1:N], v0[i] == v0[i - 1] + 0.5 * step * (dv0[i] + dv0[i - 1])
- ∂δ0[i = 1:N], δ0[i] == δ0[i - 1] + 0.5 * step * (dδ0[i] + dδ0[i - 1])
- end
- )
-
- # objective
- @expression(model, dc[i = 0:N], β01[i]^2 + β12[i]^2)
- @objective(model, Min, tf + 0.5 * step * sum(dc[i] + dc[i - 1] for i in 1:N))
-
- return model
-end
diff --git a/ext/MetaData/truck_trailer.jl b/ext/MetaData/truck_trailer.jl
deleted file mode 100644
index ea63ff38..00000000
--- a/ext/MetaData/truck_trailer.jl
+++ /dev/null
@@ -1,10 +0,0 @@
-truck_trailer_meta = OrderedDict(
- :name => "truck_trailer",
- :N => 200,
- :minimise => true,
- :state_name => ["x2", "y2", "θ0", "θ1", "θ2", "v0", "δ0"],
- :costate_name => ["∂x2", "∂y2", "∂θ0", "∂θ1", "∂θ2", "∂v0", "∂δ0"],
- :control_name => ["dv0", "dδ0"],
- :variable_name => ["tf"],
- :final_time => (:free, 1), # first component of the variable
-)
diff --git a/ext/OptimalControlModels/bioreactor.jl b/ext/OptimalControlModels/bioreactor.jl
index 4560d42c..99e5bd28 100644
--- a/ext/OptimalControlModels/bioreactor.jl
+++ b/ext/OptimalControlModels/bioreactor.jl
@@ -55,7 +55,7 @@ function OptimalControlProblems.bioreactor(
# parameters
β = 1
c = 2
- gamma = 1
+ γ = 1
halfperiod = 5
Ks = 0.05
μ2m = 0.1
@@ -81,7 +81,7 @@ function OptimalControlProblems.bioreactor(
ẋ(t) == [
μ * y(t) / (1 + y(t)) - (r + u(t)) * y(t),
- -μ2 * b(t) + u(t) * β * (gamma * y(t) - s(t)),
+ -μ2 * b(t) + u(t) * β * (γ * y(t) - s(t)),
(μ2 - u(t) * β) * b(t),
]
diff --git a/ext/OptimalControlModels/ducted_fan.jl b/ext/OptimalControlModels/ducted_fan.jl
index 81bebc22..ba077e2e 100644
--- a/ext/OptimalControlModels/ducted_fan.jl
+++ b/ext/OptimalControlModels/ducted_fan.jl
@@ -95,7 +95,7 @@ function OptimalControlProblems.ducted_fan(
# initial guess
xinit = [0.1, 0.1, -0.1, 0.1, 0.1, 0.1] # [x₁, v₁, x₂, v₂, α, vα]
uinit = [0.1, 1] # [u₁, u₂]
- varinit = [1] # [tf]
+ varinit = [1.5] # [tf]
init = (state=xinit, control=uinit, variable=varinit)
# discretise the optimal control problem
diff --git a/ext/OptimalControlModels/truck_trailer.jl b/ext/OptimalControlModels/truck_trailer.jl
deleted file mode 100644
index 668642a1..00000000
--- a/ext/OptimalControlModels/truck_trailer.jl
+++ /dev/null
@@ -1,147 +0,0 @@
-"""
-$(TYPEDSIGNATURES)
-
-Constructs an **OptimalControl problem** for a truck with two trailers, starting horizontally aligned.
-The objective is to minimise the total time required to park the truck and trailers such that they are aligned vertically at a specified target location, while respecting vehicle dynamics and control constraints.
-The problem includes path constraints for articulation angles between the trailers.
-
-# Arguments
-
-- `::OptimalControlBackend`: Placeholder type specifying the OptimalControl backend or solver interface.
-- `N::Int=200`: (Keyword) Number of discretisation points for the direct transcription grid.
-
-# Returns
-
-- `docp`: The direct optimal control problem object representing the truck-trailer parking problem.
-- `nlp`: The corresponding nonlinear programming model obtained from the DOCP, suitable for numerical optimisation.
-
-# Example
-
-```julia-repl
-julia> using OptimalControlProblems
-
-julia> docp = OptimalControlProblems.truck_trailer(OptimalControlBackend(); N=200);
-```
-"""
-function OptimalControlProblems.truck_trailer(
- ::OptimalControlBackend,
- description::Symbol...;
- N::Int=steps_number_data(:truck_trailer),
- kwargs...,
-)
-
- # parameters
- data=[0.4 0.1 0.2; 1.1 0.2 0.2; 0.8 0.1 0.2]
- L0 = data[1, 1]
- M0 = data[1, 2]
- W0 = data[1, 3]
- L1 = data[2, 1]
- M1 = data[2, 2]
- W1 = data[2, 3]
- L2 = data[3, 1]
- M2 = data[3, 2]
- W2 = data[3, 3]
- speedf = 1
- x2_t0 = 0
- y2_t0 = 0
- θ2_t0 = 0
- θ1_t0 = 0
- θ0_t0 = 0
- x2_tf = 0
- y2_tf = -2
- θ2_tf = π / 2
- θ1_tf = π / 2
- θ0_tf = π / 2
-
- # model
- ocp = @def begin
- tf ∈ R, variable
- t ∈ [0, tf], time
- x = (x2, y2, θ0, θ1, θ2, v0, δ0) ∈ R⁷, state
- u = (dv0, dδ0) ∈ R², control
-
- # auxiliary variables
- β01 = θ0 - θ1
- β12 = θ1 - θ2
- # x1 = x2 + L2 * cos(θ2) + M1 * cos(θ1)
- # y1 = y2 + L2 * sin(θ2) + M1 * sin(θ1)
- # x0 = x1 + L1 * cos(θ1) + M0 * cos(θ0)
- # y0 = y1 + L1 * sin(θ1) + M0 * sin(θ0)
-
- # final time constraints
- 1 ≤ tf ≤ 1000
-
- # state constraints
- -π / 2 ≤ θ0(t) ≤ π / 2, (θ0_con)
- -π / 2 ≤ θ1(t) ≤ π / 2, (θ1_con)
- -0.2 * speedf ≤ v0(t) ≤ 0.2 * speedf, (v0_con)
- -π / 6 ≤ δ0(t) ≤ π / 6, (δ0_con)
-
- # control constraints
- -1 ≤ dv0(t) ≤ 1, (v0_dot_con)
- -π / 10 ≤ dδ0(t) ≤ π / 10, (δ0_dot_con)
-
- # path constraints
- -π / 2 ≤ β01(t) ≤ π / 2, (β01_con)
- -π / 2 ≤ β12(t) ≤ π / 2, (β12_con)
-
- # initial conditions
- x2(0) == x2_t0, (x2_t0_con)
- y2(0) == y2_t0, (y2_t0_con)
- θ0(0) == θ0_t0, (θ0_t0_con)
- θ1(0) == θ1_t0, (θ1_t0_con)
- θ2(0) == θ2_t0, (θ2_t0_con)
-
- # final conditions
- x2(tf) == x2_tf, (x2_tf_con)
- y2(tf) == y2_tf, (y2_tf_con)
- θ2(tf) == θ2_tf, (θ2_tf_con)
- β01(tf) == θ0_tf - θ1_tf, (β01_tf_con)
- β12(tf) == θ1_tf - θ2_tf, (β12_tf_con)
-
- # dynamics
- ẋ(t) == dynamics(x(t), u(t))
-
- # objective
- tf + ∫(β01(t)^2 + β12(t)^2) → min
- end
-
- function dynamics(x, u)
- x2, y2, θ0, θ1, θ2, v0, δ0 = x
- dv0, dδ0 = u
-
- β01 = θ0 - θ1
- β12 = θ1 - θ2
-
- dθ0 = v0 / L0 * tan(δ0)
- dθ1 = v0 / L1 * sin(β01) - M0 / L1 * cos(β01) * dθ0
-
- v1 = v0 * cos(β01) + M0 * sin(β01) * dθ0
- dθ2 = v1 / L2 * sin(β12) - M1 / L2 * cos(β12) * dθ1
- v2 = v1 * cos(β12) + M1 * sin(β12) * dθ1
-
- dx2 = v2 * cos(θ2)
- dy2 = v2 * sin(θ2)
-
- return [dx2, dy2, dθ0, dθ1, dθ2, dv0, dδ0]
- end
-
- # initial guess
- xinit = [0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1] # [x2, y2, θ0, θ1, θ2, v0, δ0]
- uinit = [0.1, 0.1] # [dv0, dδ0]
- varinit = [10] # [tf]
- init = (state=xinit, control=uinit, variable=varinit)
-
- # discretise the optimal control problem
- docp = direct_transcription(
- ocp,
- description...;
- lagrange_to_mayer=false,
- init=init,
- grid_size=N,
- disc_method=:trapeze,
- kwargs...,
- )
-
- return docp
-end
diff --git a/src/OptimalControlProblems.jl b/src/OptimalControlProblems.jl
index 4c063725..08335ffa 100644
--- a/src/OptimalControlProblems.jl
+++ b/src/OptimalControlProblems.jl
@@ -27,8 +27,7 @@ function build_ocp_solution(
end
export nlp_model, ocp_model, build_ocp_solution
-
-# -----------------
+#
"""
$(TYPEDEF)
diff --git a/test/Project.toml b/test/Project.toml
index 89cf12fd..2e076448 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -5,6 +5,7 @@ CTDirect = "790bbbee-bee9-49ee-8912-a9de031322d5"
Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59"
Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
+NLPModels = "a4795742-8479-5a88-8948-cc11e1c8c1a6"
NLPModelsIpopt = "f4238b75-b362-5c4c-b852-0801c9a21d71"
OptimalControl = "5f98b655-cc9a-415a-b60e-744165666948"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
@@ -18,6 +19,7 @@ CTDirect = "0.16"
Interpolations = "0.15, 0.16"
Ipopt = "1"
JuMP = "1"
+NLPModels = "0.21"
NLPModelsIpopt = "0.10"
OptimalControl = "1"
Plots = "1"
diff --git a/test/extras/beam.jl b/test/extras/beam.jl
new file mode 100644
index 00000000..d2a80da0
--- /dev/null
+++ b/test/extras/beam.jl
@@ -0,0 +1,28 @@
+using Pkg
+Pkg.activate("test")
+using OptimalControl
+using NLPModelsIpopt
+using Plots
+
+function Beam(a)
+
+ tf = 1
+
+ ocp = @def begin
+ t ∈ [0, tf], time
+ x ∈ R², state
+ u ∈ R, control
+ x(0) == [0, 1]
+ x(tf) == [0, -1]
+ ẋ(t) == [x₂(t), u(t)]
+ 0 ≤ x₁(t) ≤ a
+ ∫(u(t)^2) → min
+ end
+
+ return ocp
+
+end
+
+plot(solve(Beam(1/2))) # state constraint inactive
+plot(solve(Beam(1/5))) # state constraint active: touch point
+plot(solve(Beam(1/8))) # state constraint active: arc
\ No newline at end of file
diff --git a/test/figures/init/beam.pdf b/test/figures/init/beam.pdf
index 4631f547..2280a2a2 100644
Binary files a/test/figures/init/beam.pdf and b/test/figures/init/beam.pdf differ
diff --git a/test/figures/init/bioreactor.pdf b/test/figures/init/bioreactor.pdf
index 91249ce4..f638d896 100644
Binary files a/test/figures/init/bioreactor.pdf and b/test/figures/init/bioreactor.pdf differ
diff --git a/test/figures/init/cart_pendulum.pdf b/test/figures/init/cart_pendulum.pdf
index 17140a66..2b4ca9ca 100644
Binary files a/test/figures/init/cart_pendulum.pdf and b/test/figures/init/cart_pendulum.pdf differ
diff --git a/test/figures/init/chain.pdf b/test/figures/init/chain.pdf
index 195d6bca..d8981612 100644
Binary files a/test/figures/init/chain.pdf and b/test/figures/init/chain.pdf differ
diff --git a/test/figures/init/dielectrophoretic_particle.pdf b/test/figures/init/dielectrophoretic_particle.pdf
index 2a8f3270..5b94cbad 100644
Binary files a/test/figures/init/dielectrophoretic_particle.pdf and b/test/figures/init/dielectrophoretic_particle.pdf differ
diff --git a/test/figures/init/double_oscillator.pdf b/test/figures/init/double_oscillator.pdf
index 220597ba..e862b180 100644
Binary files a/test/figures/init/double_oscillator.pdf and b/test/figures/init/double_oscillator.pdf differ
diff --git a/test/figures/init/ducted_fan.pdf b/test/figures/init/ducted_fan.pdf
index 9356ce84..affcf175 100644
Binary files a/test/figures/init/ducted_fan.pdf and b/test/figures/init/ducted_fan.pdf differ
diff --git a/test/figures/init/electric_vehicle.pdf b/test/figures/init/electric_vehicle.pdf
index d9544d82..376ce411 100644
Binary files a/test/figures/init/electric_vehicle.pdf and b/test/figures/init/electric_vehicle.pdf differ
diff --git a/test/figures/init/glider.pdf b/test/figures/init/glider.pdf
index 844451ca..12cfa497 100644
Binary files a/test/figures/init/glider.pdf and b/test/figures/init/glider.pdf differ
diff --git a/test/figures/init/insurance.pdf b/test/figures/init/insurance.pdf
index cd3c3e11..6acc78d2 100644
Binary files a/test/figures/init/insurance.pdf and b/test/figures/init/insurance.pdf differ
diff --git a/test/figures/init/jackson.pdf b/test/figures/init/jackson.pdf
index b90a9e1c..b98e0192 100644
Binary files a/test/figures/init/jackson.pdf and b/test/figures/init/jackson.pdf differ
diff --git a/test/figures/init/moonlander.pdf b/test/figures/init/moonlander.pdf
index bac2e36d..4361a558 100644
Binary files a/test/figures/init/moonlander.pdf and b/test/figures/init/moonlander.pdf differ
diff --git a/test/figures/init/quadrotor.pdf b/test/figures/init/quadrotor.pdf
deleted file mode 100644
index b3a37e89..00000000
Binary files a/test/figures/init/quadrotor.pdf and /dev/null differ
diff --git a/test/figures/init/robbins.pdf b/test/figures/init/robbins.pdf
index 299a93ae..2936a89c 100644
Binary files a/test/figures/init/robbins.pdf and b/test/figures/init/robbins.pdf differ
diff --git a/test/figures/init/robot.pdf b/test/figures/init/robot.pdf
index fd55b787..4cc3b6c1 100644
Binary files a/test/figures/init/robot.pdf and b/test/figures/init/robot.pdf differ
diff --git a/test/figures/init/rocket.pdf b/test/figures/init/rocket.pdf
index ae7dd14b..7f1a025e 100644
Binary files a/test/figures/init/rocket.pdf and b/test/figures/init/rocket.pdf differ
diff --git a/test/figures/init/space_shuttle.pdf b/test/figures/init/space_shuttle.pdf
index 6149a5d1..7dfa0923 100644
Binary files a/test/figures/init/space_shuttle.pdf and b/test/figures/init/space_shuttle.pdf differ
diff --git a/test/figures/init/steering.pdf b/test/figures/init/steering.pdf
index 95e8da72..5d798268 100644
Binary files a/test/figures/init/steering.pdf and b/test/figures/init/steering.pdf differ
diff --git a/test/figures/init/truck_trailer.pdf b/test/figures/init/truck_trailer.pdf
deleted file mode 100644
index 7e3a3522..00000000
Binary files a/test/figures/init/truck_trailer.pdf and /dev/null differ
diff --git a/test/figures/init/vanderpol.pdf b/test/figures/init/vanderpol.pdf
index 96ea7506..ed3e83e9 100644
Binary files a/test/figures/init/vanderpol.pdf and b/test/figures/init/vanderpol.pdf differ
diff --git a/test/figures/solution/beam.pdf b/test/figures/solution/beam.pdf
index 311dcfeb..50bedd09 100644
Binary files a/test/figures/solution/beam.pdf and b/test/figures/solution/beam.pdf differ
diff --git a/test/figures/solution/bioreactor.pdf b/test/figures/solution/bioreactor.pdf
index d5dedaeb..639807c9 100644
Binary files a/test/figures/solution/bioreactor.pdf and b/test/figures/solution/bioreactor.pdf differ
diff --git a/test/figures/solution/cart_pendulum.pdf b/test/figures/solution/cart_pendulum.pdf
index a8796c9f..a01b8f4b 100644
Binary files a/test/figures/solution/cart_pendulum.pdf and b/test/figures/solution/cart_pendulum.pdf differ
diff --git a/test/figures/solution/chain.pdf b/test/figures/solution/chain.pdf
index 6aeaa820..b60dc74f 100644
Binary files a/test/figures/solution/chain.pdf and b/test/figures/solution/chain.pdf differ
diff --git a/test/figures/solution/dielectrophoretic_particle.pdf b/test/figures/solution/dielectrophoretic_particle.pdf
index 6fe975a5..01ab9a66 100644
Binary files a/test/figures/solution/dielectrophoretic_particle.pdf and b/test/figures/solution/dielectrophoretic_particle.pdf differ
diff --git a/test/figures/solution/double_oscillator.pdf b/test/figures/solution/double_oscillator.pdf
index 7744a0e6..1fadf310 100644
Binary files a/test/figures/solution/double_oscillator.pdf and b/test/figures/solution/double_oscillator.pdf differ
diff --git a/test/figures/solution/ducted_fan.pdf b/test/figures/solution/ducted_fan.pdf
index e8304c7d..8b20d6dc 100644
Binary files a/test/figures/solution/ducted_fan.pdf and b/test/figures/solution/ducted_fan.pdf differ
diff --git a/test/figures/solution/electric_vehicle.pdf b/test/figures/solution/electric_vehicle.pdf
index a3151548..75ed5771 100644
Binary files a/test/figures/solution/electric_vehicle.pdf and b/test/figures/solution/electric_vehicle.pdf differ
diff --git a/test/figures/solution/glider.pdf b/test/figures/solution/glider.pdf
index 4be1a221..ada9b5f6 100644
Binary files a/test/figures/solution/glider.pdf and b/test/figures/solution/glider.pdf differ
diff --git a/test/figures/solution/insurance.pdf b/test/figures/solution/insurance.pdf
index cbe6996a..503d82e1 100644
Binary files a/test/figures/solution/insurance.pdf and b/test/figures/solution/insurance.pdf differ
diff --git a/test/figures/solution/jackson.pdf b/test/figures/solution/jackson.pdf
index 6fc6c93d..19880023 100644
Binary files a/test/figures/solution/jackson.pdf and b/test/figures/solution/jackson.pdf differ
diff --git a/test/figures/solution/moonlander.pdf b/test/figures/solution/moonlander.pdf
index 62b7c641..cd118905 100644
Binary files a/test/figures/solution/moonlander.pdf and b/test/figures/solution/moonlander.pdf differ
diff --git a/test/figures/solution/quadrotor.pdf b/test/figures/solution/quadrotor.pdf
deleted file mode 100644
index 2f009bc0..00000000
Binary files a/test/figures/solution/quadrotor.pdf and /dev/null differ
diff --git a/test/figures/solution/robbins.pdf b/test/figures/solution/robbins.pdf
index e752cb0d..7d64f510 100644
Binary files a/test/figures/solution/robbins.pdf and b/test/figures/solution/robbins.pdf differ
diff --git a/test/figures/solution/robot.pdf b/test/figures/solution/robot.pdf
index 883df745..d88efe26 100644
Binary files a/test/figures/solution/robot.pdf and b/test/figures/solution/robot.pdf differ
diff --git a/test/figures/solution/rocket.pdf b/test/figures/solution/rocket.pdf
index ccdd6363..f05ba447 100644
Binary files a/test/figures/solution/rocket.pdf and b/test/figures/solution/rocket.pdf differ
diff --git a/test/figures/solution/space_shuttle.pdf b/test/figures/solution/space_shuttle.pdf
index 70d8f3e5..453ba5e5 100644
Binary files a/test/figures/solution/space_shuttle.pdf and b/test/figures/solution/space_shuttle.pdf differ
diff --git a/test/figures/solution/steering.pdf b/test/figures/solution/steering.pdf
index 4be292ca..dd15cf2f 100644
Binary files a/test/figures/solution/steering.pdf and b/test/figures/solution/steering.pdf differ
diff --git a/test/figures/solution/truck_trailer.pdf b/test/figures/solution/truck_trailer.pdf
deleted file mode 100644
index ea5788f5..00000000
Binary files a/test/figures/solution/truck_trailer.pdf and /dev/null differ
diff --git a/test/figures/solution/vanderpol.pdf b/test/figures/solution/vanderpol.pdf
index f3bebda9..b28704f9 100644
Binary files a/test/figures/solution/vanderpol.pdf and b/test/figures/solution/vanderpol.pdf differ
diff --git a/test/runtests.jl b/test/runtests.jl
index a4df4153..aa1491bf 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -3,12 +3,14 @@ using CTBase
using CTDirect
using Ipopt
using JuMP
+using NLPModels
using NLPModelsIpopt
using OptimalControl
using OptimalControlProblems
using Test
using Plots
using Plots.PlotMeasures # for leftmargin, bottommargin
+using Printf
using Interpolations
include("utils.jl")
@@ -35,7 +37,7 @@ end
# list_of_problems = setdiff(list_of_problems, problems_to_exclude)
# list_of_problems = [
-# :quadrotor
+# :truck_trailer
# ]
# The list of all the problems to test
diff --git a/test/test_JuMP.jl b/test/test_JuMP.jl
index c53c8dbe..2f0ac509 100644
--- a/test/test_JuMP.jl
+++ b/test/test_JuMP.jl
@@ -41,13 +41,9 @@ function test_JuMP()
# Test
res = @my_test_broken termination_status(nlp) == MOI.LOCALLY_SOLVED
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
+ keep_problem = keep_problem && res
+ DEBUG && res && println("│ \033[1;32mTest Passed\033[0m")
+ DEBUG && !res && println("│ \033[1;31mTest Failed\033[0m")
DEBUG && println("│")
DEBUG && println("└─")
diff --git a/test/test_OptimalControl.jl b/test/test_OptimalControl.jl
index b7532b57..0596c966 100644
--- a/test/test_OptimalControl.jl
+++ b/test/test_OptimalControl.jl
@@ -43,13 +43,9 @@ function test_OptimalControl()
# Test
res = @my_test_broken (sol.status == :first_order || sol.status == :acceptable)
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
+ keep_problem = keep_problem && res
+ DEBUG && res && println("│ \033[1;32mTest Passed\033[0m")
+ DEBUG && !res && println("│ \033[1;31mTest Failed\033[0m")
DEBUG && println("│")
DEBUG && println("└─")
diff --git a/test/test_quick.jl b/test/test_quick.jl
index 5760ab26..d83959b0 100644
--- a/test/test_quick.jl
+++ b/test/test_quick.jl
@@ -60,12 +60,8 @@ function test_quick()
res = @my_test_broken o_di < o_bd
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
+ DEBUG && res && println("│ \033[1;32mTest Passed\033[0m")
+ DEBUG && !res && println("│ \033[1;31mTest Failed\033[0m")
DEBUG && println("│")
max_r_err = max(max_r_err, o_di/(0.5*(abs(o_oc) + abs(o_jp))))
diff --git a/test/utils.jl b/test/utils.jl
index 4355ae31..7dd4e43e 100644
--- a/test/utils.jl
+++ b/test/utils.jl
@@ -51,7 +51,8 @@ julia> @macroexpand @my_test_broken 1 == 2
"""
macro my_test_broken(e)
return esc(quote
- @test $e broken=!$e
+ res = @test $e broken=!$e
+ typeof(res) == Test.Pass
end)
end
@@ -115,272 +116,185 @@ function comparison(; max_iter, test_name)
v_vars = metadata[f][:variable_name]
@testset "$(string(f)) ($(string(test_name)))" verbose=VERBOSE begin
- DEBUG && println("\n", "┌─ ", string(f), " (", string(test_name), ")")
+ DEBUG && println("\n┌─ ", string(f), " (", string(test_name), ")")
DEBUG && println("│")
########## OptimalControl ##########
-
- # set up the OptimalControl model
docp = OptimalControlProblems.eval(f)(OptimalControlBackend(); N=N)
- nlp = nlp_model(docp)
-
- # solve the problem
- nlp_sol = NLPModelsIpopt.ipopt(nlp; Options...)
-
- # build the solution
+ nlp_oc = nlp_model(docp)
+ nlp_sol = NLPModelsIpopt.ipopt(nlp_oc; Options...)
sol = build_ocp_solution(docp, nlp_sol)
+ sol_oc = deepcopy(sol) # for plotting
- sol_oc = deepcopy(sol) # for plotting
-
- # retrieve values of variables that we compare
t_oc = time_grid(sol)
x_oc = state(sol).(t_oc)
u_oc = control(sol).(t_oc)
o_oc = objective(sol)
- i_oc = nlp_sol.iter # iterations(sol) returns 0!
+ i_oc = iterations(sol)
v_oc = variable(sol)
############### JuMP ###############
+ nlp_jp = OptimalControlProblems.eval(f)(JuMPBackend(); N=N)
+ set_optimizer(nlp_jp, Ipopt.Optimizer)
+ set_silent(nlp_jp)
+ set_optimizer_attribute(nlp_jp, "tol", Options[:tol])
+ set_optimizer_attribute(nlp_jp, "max_iter", Options[:max_iter])
+ set_optimizer_attribute(nlp_jp, "mu_strategy", Options[:mu_strategy])
+ set_optimizer_attribute(nlp_jp, "linear_solver", "mumps")
+ set_optimizer_attribute(nlp_jp, "max_wall_time", Options[:max_wall_time])
+ set_optimizer_attribute(nlp_jp, "sb", Options[:sb])
+ optimize!(nlp_jp)
+
+ t_jp = time_grid(f, nlp_jp)
+ x_jp = state(f, nlp_jp).(t_jp)
+ u_jp = control(f, nlp_jp).(t_jp)
+ o_jp = objective(f, nlp_jp)
+ i_jp = iterations(f, nlp_jp)
+ v_jp = variable(f, nlp_jp)
+ p_jp = costate(f, nlp_jp).(t_jp)
+
+ ########## Iterations ##########
+ DEBUG && println("├─ Iterations → OC: ", i_oc, ", JP: ", i_jp)
- # set up the JuMP model
- nlp = OptimalControlProblems.eval(f)(JuMPBackend(); N=N)
- set_optimizer(nlp, Ipopt.Optimizer)
- set_silent(nlp)
- set_optimizer_attribute(nlp, "tol", Options[:tol])
- set_optimizer_attribute(nlp, "max_iter", Options[:max_iter])
- set_optimizer_attribute(nlp, "mu_strategy", Options[:mu_strategy])
- set_optimizer_attribute(nlp, "linear_solver", "mumps")
- set_optimizer_attribute(nlp, "max_wall_time", Options[:max_wall_time])
- set_optimizer_attribute(nlp, "sb", Options[:sb])
-
- # solve the model
- optimize!(nlp)
-
- # retrieve values
- t_jp = time_grid(f, nlp)
- x_jp = state(f, nlp).(t_jp)
- u_jp = control(f, nlp).(t_jp)
- o_jp = objective_value(nlp)
- i_jp = barrier_iterations(nlp)
- p_jp = costate(f, nlp).(t_jp)
- v_jp = variable(f, nlp)
-
- ############ TEST ############
-
- DEBUG && println("├─ iterations")
- DEBUG && println("│")
- DEBUG && println("│ i_oc = ", i_oc)
- DEBUG && println("│ i_jp = ", i_jp)
- DEBUG && println("│")
-
- # do we keep or remove the problem from the list
keep_problem = true
- # time grids
+ ########## Variables / Constraints ##########
+ @testset "nlp" verbose=VERBOSE begin
+ nb_var_oc, nb_var_jp = get_nvar(nlp_oc), num_variables(nlp_jp)
+ res = @my_test_broken nb_var_oc == nb_var_jp
+ keep_problem = keep_problem && res
+ DEBUG && @printf("├─ Variables → OC: %d JP: %d %s\n", nb_var_oc, nb_var_jp,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
+
+ nb_con_oc = get_ncon(nlp_oc)
+ nb_con_jp = num_constraints(nlp_jp; count_variable_in_set_constraints=false)
+ res = @my_test_broken nb_con_oc == nb_con_jp
+ keep_problem = keep_problem && res
+ DEBUG && @printf("├─ Constraints → OC: %d JP: %d %s\n", nb_con_oc, nb_con_jp,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
+ end
+
+ ########## Time Grid ##########
test_grid_ok = true
@testset "grid" verbose=VERBOSE begin
-
# final time
- DEBUG && println("├─ final time")
- DEBUG && println("│")
- DEBUG && println("│ tf oc = ", t_oc[end])
- DEBUG && println("│ tf jp = ", t_jp[end])
- DEBUG && println("│")
-
- # length
+ tf_di = abs(t_oc[end] - t_jp[end])
+ tf_bd = max(0.5*(t_oc[end]+t_jp[end])*ε_rel_grid, ε_abs_grid)
+ res = @my_test_broken tf_di < tf_bd
+ r_err = tf_di / (0.5*(t_oc[end]+t_jp[end]))
+ DEBUG && @printf("├─ Final time → OC: %.3e JP: %.3e\n", t_oc[end], t_jp[end])
+ DEBUG && @printf("│ r_err=%.3e a_err=%.3e bound=%.3e %s\n",
+ r_err, tf_di, tf_bd,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
+
+ # length of the grids
res = @my_test_broken length(t_oc) == length(t_jp)
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
- test_grid_ok = test_grid_ok && (typeof(res) == Test.Pass)
-
- DEBUG && println("├─ grid length")
- DEBUG && println("│")
- DEBUG && println("│ length(t_oc) = ", length(t_oc))
- DEBUG && println("│ length(t_jp) = ", length(t_jp))
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
- DEBUG && println("│")
-
- # values
- ti_oc_max = NaN
- ti_jp_max = NaN
- ti_di_max = NaN
- ti_bd_max = NaN
- ti_te_max = Inf
- itera_max = NaN
+ keep_problem = keep_problem && res
+ test_grid_ok = test_grid_ok && res
+ DEBUG && @printf("├─ Grid length → OC: %d JP: %d %s\n", length(t_oc), length(t_jp),
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
+
+ # max error
if test_grid_ok
+ ti_di_max, ti_bd_max, itera_max = 0, NaN, 0
for i in eachindex(t_oc)
ti_di = t_oc[i] - t_jp[i]
- ti_bd = max(
- 0.5*(abs(t_oc[i]) + abs(t_jp[i]))*ε_rel_grid, ε_abs_grid
- )
- res = @my_test_broken ti_di < ti_bd
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
- test_grid_ok = test_grid_ok && (typeof(res) == Test.Pass)
- if ti_bd - ti_di < ti_te_max
- ti_te_max = ti_bd - ti_di
- itera_max = i
- ti_oc_max = t_oc[i]
- ti_jp_max = t_jp[i]
- ti_di_max = ti_di
- ti_bd_max = ti_bd
+ ti_bd = max(0.5*(abs(t_oc[i])+abs(t_jp[i]))*ε_rel_grid, ε_abs_grid)
+ res = @my_test_broken abs(ti_di) < ti_bd
+ keep_problem = keep_problem && res
+ test_grid_ok = test_grid_ok && res
+ if abs(ti_di) ≥ abs(ti_di_max)
+ ti_di_max, ti_bd_max, itera_max = ti_di, ti_bd, i
end
end
+ r_err = abs(ti_di_max)/(0.5*(abs(t_oc[itera_max])+abs(t_jp[itera_max])))
+ DEBUG && @printf("├─ Grid max error → iter=%d r_err=%.3e a_err=%.3e bound=%.3e %s\n",
+ itera_max, r_err, abs(ti_di_max), ti_bd_max,
+ r_err<1.0 ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
end
-
- DEBUG && println("├─ grid values (max error)")
- DEBUG && println("│")
- DEBUG && println("│ iter = ", itera_max)
- DEBUG && println("│ ti oc = ", ti_oc_max)
- DEBUG && println("│ ti jp = ", ti_jp_max)
- DEBUG && println(
- "│ r_err = ", ti_di_max/(0.5*(abs(ti_oc_max) + abs(ti_jp_max)))
- )
- DEBUG && println("│ a_err = ", ti_di_max)
- DEBUG && println("│ bound = ", ti_bd_max)
- DEBUG && test_grid_ok && println("│ \033[1;32mTest Passed\033[0m")
- DEBUG && !test_grid_ok && println("│ \033[1;31mTest Failed\033[0m")
- DEBUG && println("│")
end
- # state
+ ########## States ##########
if test_grid_ok
@testset "state" verbose=VERBOSE begin
+ DEBUG && println("├─ States")
for i in eachindex(x_vars)
@testset "$(x_vars[i])" verbose=VERBOSE begin
- xi_oc = [x_oc[k][i] for k in eachindex(t_oc)]
- xi_jp = [x_jp[k][i] for k in eachindex(t_jp)]
- L2_di = L2_norm(t_oc, xi_oc-xi_jp)
- L2_oc = L2_norm(t_oc, xi_oc)
- L2_jp = L2_norm(t_oc, xi_jp)
- L2_bd = max(0.5*(L2_oc + L2_jp)*ε_rel_state, ε_abs_state)
+ xi_oc, xi_jp = [x_oc[k][i] for k in eachindex(t_oc)], [x_jp[k][i] for k in eachindex(t_jp)]
+ L2_di = L2_norm(t_oc, xi_oc - xi_jp)
+ L2_bd = max(0.5*(L2_norm(t_oc, xi_oc)+L2_norm(t_oc, xi_jp))*ε_rel_state, ε_abs_state)
res = @my_test_broken L2_di < L2_bd
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
-
- DEBUG && println("├─ state $(x_vars[i])")
- DEBUG && println("│")
- DEBUG && println("│ L2 oc = ", L2_oc)
- DEBUG && println("│ L2 jp = ", L2_jp)
- DEBUG && println("│ r_err = ", L2_di/(0.5*(L2_oc + L2_jp)))
- DEBUG && println("│ a_err = ", L2_di)
- DEBUG && println("│ bound = ", L2_bd)
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
- DEBUG && println("│")
+ keep_problem = keep_problem && res
+ r_err = L2_di / (0.5*(L2_norm(t_oc, xi_oc)+L2_norm(t_oc, xi_jp)))
+ DEBUG && @printf("│ %-6s r_err=%.3e a_err=%.3e bound=%.3e %s\n",
+ x_vars[i], r_err, L2_di, L2_bd,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
end
end
end
end
- # control
+ ########## Controls ##########
if test_grid_ok
@testset "control" verbose=VERBOSE begin
+ DEBUG && println("├─ Controls")
for i in eachindex(u_vars)
@testset "$(u_vars[i])" verbose=VERBOSE begin
- ui_oc = [u_oc[k][i] for k in eachindex(t_oc)]
- ui_jp = [u_jp[k][i] for k in eachindex(t_jp)]
- L2_di = L2_norm(t_oc, ui_oc-ui_jp)
- L2_oc = L2_norm(t_oc, ui_oc)
- L2_jp = L2_norm(t_oc, ui_jp)
- L2_bd = max(0.5*(L2_oc + L2_jp)*ε_rel_control, ε_abs_control)
+ ui_oc, ui_jp = [u_oc[k][i] for k in eachindex(t_oc)], [u_jp[k][i] for k in eachindex(t_jp)]
+ L2_di = L2_norm(t_oc, ui_oc - ui_jp)
+ L2_bd = max(0.5*(L2_norm(t_oc, ui_oc)+L2_norm(t_oc, ui_jp))*ε_rel_control, ε_abs_control)
res = @my_test_broken L2_di < L2_bd
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
-
- DEBUG && println("├─ control $(u_vars[i])")
- DEBUG && println("│")
- DEBUG && println("│ L2 oc = ", L2_oc)
- DEBUG && println("│ L2 jp = ", L2_jp)
- DEBUG && println("│ r_err = ", L2_di/(0.5*(L2_oc + L2_jp)))
- DEBUG && println("│ a_err = ", L2_di)
- DEBUG && println("│ bound = ", L2_bd)
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
- DEBUG && println("│")
+ keep_problem = keep_problem && res
+ r_err = L2_di / (0.5*(L2_norm(t_oc, ui_oc)+L2_norm(t_oc, ui_jp)))
+ DEBUG && @printf("│ %-6s r_err=%.3e a_err=%.3e bound=%.3e %s\n",
+ u_vars[i], r_err, L2_di, L2_bd,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
end
end
end
end
- # variable
+ ########## Variables ##########
if test_grid_ok && !isnothing(v_vars)
@testset "variable" verbose=VERBOSE begin
+ DEBUG && println("├─ Variables")
for i in eachindex(v_vars)
@testset "$(v_vars[i])" verbose=VERBOSE begin
- vi_oc = v_oc[i]
- vi_jp = v_jp[i]
+ vi_oc, vi_jp = v_oc[i], v_jp[i]
vi_di = abs(vi_oc-vi_jp)
- vi_bd = max(
- 0.5*(abs(vi_oc) + abs(vi_jp))*ε_rel_control, ε_abs_control
- )
+ vi_bd = max(0.5*(abs(vi_oc)+abs(vi_jp))*ε_rel_control, ε_abs_control)
res = @my_test_broken vi_di < vi_bd
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
-
- DEBUG && println("├─ variable $(v_vars[i])")
- DEBUG && println("│")
- DEBUG && println("│ vi oc = ", vi_oc)
- DEBUG && println("│ vi jp = ", vi_jp)
- DEBUG && println(
- "│ r_err = ", vi_di/(0.5*(abs(vi_oc) + abs(vi_jp)))
- )
- DEBUG && println("│ a_err = ", vi_di)
- DEBUG && println("│ bound = ", vi_bd)
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
- DEBUG && println("│")
+ keep_problem = keep_problem && res
+ r_err = vi_di / (0.5*(abs(vi_oc)+abs(vi_jp)))
+ DEBUG && @printf("│ %-6s r_err=%.3e a_err=%.3e bound=%.3e %s\n",
+ v_vars[i], r_err, vi_di, vi_bd,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
end
end
end
end
- # objective
+ ########## Objective ##########
@testset "objective" verbose=VERBOSE begin
o_di = abs(o_oc-o_jp)
- o_bd = max(0.5*(abs(o_oc) + abs(o_jp))*ε_rel_objective, ε_abs_objective)
+ o_bd = max(0.5*(abs(o_oc)+abs(o_jp))*ε_rel_objective, ε_abs_objective)
res = @my_test_broken o_di < o_bd
- keep_problem = keep_problem && (typeof(res) == Test.Pass)
-
- DEBUG && println("├─ objective")
- DEBUG && println("│")
- DEBUG && println("│ o_oc = ", o_oc)
- DEBUG && println("│ o_jp = ", o_jp)
- DEBUG && println("│ r_err = ", o_di/(0.5*(abs(o_oc) + abs(o_jp))))
- DEBUG && println("│ a_err = ", o_di)
- DEBUG && println("│ bound = ", o_bd)
- DEBUG &&
- (typeof(res) == Test.Pass) &&
- println("│ \033[1;32mTest Passed\033[0m")
- DEBUG &&
- (typeof(res) != Test.Pass) &&
- println("│ \033[1;31mTest Failed\033[0m")
- DEBUG && println("│")
+ keep_problem = keep_problem && res
+ r_err = o_di / (0.5*(abs(o_oc)+abs(o_jp)))
+ DEBUG && println("├─ Objective")
+ DEBUG && @printf("│ r_err=%.3e a_err=%.3e bound=%.3e %s\n",
+ r_err, o_di, o_bd,
+ res ? "\033[1;32mPASS\033[0m" : "\033[1;31mFAIL\033[0m")
end
- #
DEBUG && println("└─")
- # do we keep or remove the problem from the list
if !keep_problem
global LIST_OF_PROBLEMS_FINAL
LIST_OF_PROBLEMS_FINAL = setdiff(LIST_OF_PROBLEMS_FINAL, [f])
end
############ PLOT ############
-
figdir = joinpath(@__DIR__, "figures", string(test_name))
isdir(figdir) || mkpath(figdir)
@@ -429,6 +343,7 @@ function comparison(; max_iter, test_name)
# save figure
savefig(plt, joinpath(figdir, "$f" * ".pdf"))
- end
- end
+
+ end# end testset
+ end # end for
end