Public API

This page lists exported symbols of CTModels.OCP.


From CTModels.OCP

_serialize_solution

CTModels.OCP._serialize_solutionFunction
_serialize_solution(
    sol::CTModels.OCP.Solution
) -> Dict{String, Any}

Serialize a solution into discrete data for export to persistent storage (JLD2, JSON, etc.).

This function converts a Solution object (which may contain interpolated functions) into a fully discrete, serializable representation. All trajectory functions are evaluated on their respective time grids and stored as matrices. The serialization format automatically adapts based on whether the solution uses unified or multiple time grids.

Serialization Formats

The function produces two different formats depending on the solution's time grid model:

Unified Time Grid Format (Legacy)

When all grids are identical (UnifiedTimeGridModel), produces:

Dict(
    "time_grid" => T,                    # Single grid for all components
    "state" => Matrix,                   # Discretized on T
    "control" => Matrix,                 # Discretized on T
    "costate" => Matrix,                 # Discretized on T
    "path_constraints_dual" => Matrix,   # Discretized on T
    # ... other fields
)

Multiple Time Grids Format

When grids differ (MultipleTimeGridModel), produces:

Dict(
    "time_grid_state" => T_state,        # State-specific grid
    "time_grid_control" => T_control,    # Control-specific grid
    "time_grid_costate" => T_costate,    # Costate-specific grid
    "time_grid_path" => T_path,          # Path constraints grid
    "state" => Matrix,                   # Discretized on T_state
    "control" => Matrix,                 # Discretized on T_control
    "costate" => Matrix,                 # Discretized on T_costate
    "path_constraints_dual" => Matrix,   # Discretized on T_path
    # ... other fields
)

Arguments

  • sol::Solution: Solution object to serialize (may contain functions or matrices)

Returns

  • Dict{String, Any}: Complete serializable dictionary containing:
    • Time grids: Either single "time_grid" or four separate grids
    • Trajectories: "state", "control", "costate" as Matrix{Float64}
    • Variable: "variable" as Vector{Float64} (time-independent)
    • Objective: "objective" as Float64
    • Dual variables: All constraint duals (can be nothing if not present)
      • "path_constraints_dual": Path constraint duals on path grid
      • "state_constraints_lb_dual", "state_constraints_ub_dual": State box duals on state grid
      • "control_constraints_lb_dual", "control_constraints_ub_dual": Control box duals on control grid
      • "boundary_constraints_dual": Boundary duals (time-independent vector)
      • "variable_constraints_lb_dual", "variable_constraints_ub_dual": Variable duals (vectors)
    • Solver info: "iterations", "message", "status", "successful", "constraints_violation", "infos"

Discretization Behavior

  • Function trajectories: Evaluated at each point of their associated time grid
  • Matrix trajectories: Copied as-is (already discrete)
  • Nothing duals: Preserved as nothing in the dictionary
  • Grid association: Each component is discretized on its correct grid:
    • State and state box duals → T_state
    • Control and control box duals → T_control
    • Costate → T_costate
    • Path constraint duals → T_path

Example

using CTModels

# Solve OCP with multiple grids
sol = solve(ocp, strategy=MyStrategy())

# Serialize to dictionary
data = _serialize_solution(sol)

# Check format
if haskey(data, "time_grid_state")
    # Multiple grids format
    println("State grid: ", length(data["time_grid_state"]), " points")
    println("Control grid: ", length(data["time_grid_control"]), " points")
    println("Costate grid: ", length(data["time_grid_costate"]), " points")
else
    # Unified grid format
    println("Unified grid: ", length(data["time_grid"]), " points")
end

# Export to file (handled by extensions)
export_ocp_solution(sol; filename="solution", format=:JLD)

# Reconstruct from data
sol_reconstructed = _reconstruct_solution_from_data(ocp, data)

Notes

Backward Compatibility

The serialization format is designed for backward compatibility:

  • Old files with single "time_grid" can be read (costate defaults to state grid)
  • New files with four grids are forward-compatible with updated readers
  • The _reconstruct_solution_from_data function handles both formats automatically

Memory Efficiency

When all grids are identical, the unified format avoids storing redundant grid data, reducing file size and memory usage.

Round-Trip Guarantee

The serialized data is fully compatible with build_solution for exact reconstruction:

data = _serialize_solution(sol)
sol_new = build_solution(ocp, data["time_grid_state"], ...; objective=data["objective"], ...)

See also: build_solution, _reconstruct_solution_from_data, export_ocp_solution, import_ocp_solution

_serialize_solution(
    _::CTModels.OCP.UnifiedTimeGridModel,
    sol::CTModels.OCP.Solution,
    dim_x::Int64,
    dim_u::Int64
) -> Dict{String, Any}

Serialize solution with unified time grid (legacy single-grid format).

This method handles solutions where all components share the same time grid. It produces the legacy format with a single "time_grid" key, which is backward-compatible with older versions of the package.

Format Produced

Dict(
    "time_grid" => T,                    # Single unified grid
    "state" => Matrix,                   # All components discretized on T
    "control" => Matrix,
    "costate" => Matrix,
    # ... all other fields
)

Arguments

  • ::UnifiedTimeGridModel: Time grid model type (dispatch parameter)
  • sol::Solution: Solution to serialize
  • dim_x::Int: State dimension
  • dim_u::Int: Control dimension

Returns

  • Dict{String, Any}: Serialized data with single time grid

Notes

This format is used when build_solution is called with identical grids for all components, or when using the legacy single-grid signature. It ensures backward compatibility with files created before the multi-grid feature was introduced.

See also: _serialize_solution(::MultipleTimeGridModel, ...)

_serialize_solution(
    _::CTModels.OCP.MultipleTimeGridModel,
    sol::CTModels.OCP.Solution,
    dim_x::Int64,
    dim_u::Int64
) -> Dict{String, Any}

Serialize solution with multiple independent time grids (modern format).

This method handles solutions where different components use different time grids. It produces the modern format with four separate grid keys (time_grid_state, time_grid_control, time_grid_costate, time_grid_path), preserving the independent discretizations.

Format Produced

Dict(
    "time_grid_state" => T_state,        # State-specific grid
    "time_grid_control" => T_control,    # Control-specific grid
    "time_grid_costate" => T_costate,    # Costate-specific grid
    "time_grid_path" => T_path,          # Path constraints grid
    "state" => Matrix,                   # Discretized on T_state
    "control" => Matrix,                 # Discretized on T_control
    "costate" => Matrix,                 # Discretized on T_costate
    "path_constraints_dual" => Matrix,   # Discretized on T_path
    # ... all other fields
)

Arguments

  • ::MultipleTimeGridModel: Time grid model type (dispatch parameter)
  • sol::Solution: Solution to serialize
  • dim_x::Int: State dimension
  • dim_u::Int: Control dimension

Returns

  • Dict{String, Any}: Serialized data with four independent time grids

Notes

This format is used when build_solution is called with different grids for different components. It allows numerical schemes to use optimal discretizations for each component (e.g., finer grid for state, coarser for control, custom for costate).

The reconstruction function _reconstruct_solution_from_data detects this format by checking for the presence of "time_grid_state" key and handles it appropriately.

See also: _serialize_solution(::UnifiedTimeGridModel, ...), build_solution

append_box_constraints!

CTModels.OCP.append_box_constraints!Function
append_box_constraints!(
    inds,
    lbs,
    ubs,
    labels,
    rg,
    lb,
    ub,
    label
)

Appends box constraint data to the provided vectors.

Arguments

  • inds::Vector{Int}: Vector of indices to which the range rg will be appended.
  • lbs::Vector{<:Real}: Vector of lower bounds to which lb will be appended.
  • ubs::Vector{<:Real}: Vector of upper bounds to which ub will be appended.
  • labels::Vector{String}: Vector of labels to which the label will be repeated and appended.
  • rg::AbstractVector{Int}: Index range corresponding to the constraint variables.
  • lb::AbstractVector{<:Real}: Lower bounds associated with rg.
  • ub::AbstractVector{<:Real}: Upper bounds associated with rg.
  • label::String: Label describing the constraint block (e.g., "state", "control").

Notes

  • All input vectors (rg, lb, ub) must have the same length.
  • The function modifies the inds, lbs, ubs, and labels vectors in-place.
  • If a component index already exists in inds, a warning is emitted indicating that the previous bound will be overwritten by the new constraint. The dual variable dimension remains equal to the state/control/variable dimension, not the number of constraint declarations.

boundary_constraints_dual

CTModels.OCP.boundary_constraints_dualFunction
boundary_constraints_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, BC_Dual<:Union{Nothing, AbstractVector{<:Real}}}
) -> Union{Nothing, AbstractVector{<:Real}}

Return the dual vector associated with the boundary constraints.

Arguments

  • model::DualModel: A model including dual variables for boundary constraints.

Returns

A vector of dual values, or nothing if not set.

boundary_constraints_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, AbstractVector{<:Real}}

Return the dual of the boundary constraints.

boundary_constraints_nl

CTModels.OCP.boundary_constraints_nlFunction
boundary_constraints_nl(
    model::CTModels.OCP.ConstraintsModel{<:Tuple, TB}
) -> Any

Get the nonlinear boundary constraints from the model.

Arguments

  • model: The constraints model from which to retrieve the boundary constraints.

Returns

  • The nonlinear boundary constraints.

Example

# Example of retrieving nonlinear boundary constraints
julia> model = ConstraintsModel(...)
julia> boundary_constraints = boundary_constraints_nl(model)
boundary_constraints_nl(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.AbstractObjectiveModel, <:CTModels.OCP.ConstraintsModel{<:Tuple, TB<:Tuple}}
) -> Any

Return the nonlinear boundary constraints.

build

CTModels.OCP.buildFunction
build(
    constraints::OrderedCollections.OrderedDict{Symbol, Tuple{Symbol, Union{Function, OrdinalRange{<:Int64}}, AbstractVector{<:Real}, AbstractVector{<:Real}}}
) -> CTModels.OCP.ConstraintsModel{TP, TB, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}} where {TP<:Tuple{Vector{Float64}, Any, Vector{Float64}, Vector{Symbol}}, TB<:Tuple{Vector{Float64}, Any, Vector{Float64}, Vector{Symbol}}}

Constructs a ConstraintsModel from a dictionary of constraints.

This function processes a dictionary where each entry defines a constraint with its type, function or index range, lower and upper bounds, and label. It categorizes constraints into path, boundary, state, control, and variable constraints, assembling them into a structured ConstraintsModel.

Arguments

  • constraints::ConstraintsDictType: A dictionary mapping constraint labels to tuples of the form (type, function_or_range, lower_bound, upper_bound).

Returns

  • ConstraintsModel: A structured model encapsulating all provided constraints.

Example

julia> constraints = OrderedDict(
    :c1 => (:path, f1, [0.0], [1.0]),
    :c2 => (:state, 1:2, [-1.0, -1.0], [1.0, 1.0])
)
julia> model = build(constraints)
build(
    pre_ocp::CTModels.OCP.PreModel;
    build_examodel
) -> CTModels.OCP.Model{TD, var"#s179", var"#s1791", var"#s1792", var"#s1793", var"#s1794", var"#s1795", CTModels.OCP.ConstraintsModel{TP, TB, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}}, Nothing} where {TD<:CTModels.OCP.TimeDependence, var"#s179"<:CTModels.OCP.AbstractTimesModel, var"#s1791"<:CTModels.OCP.AbstractStateModel, var"#s1792"<:CTModels.OCP.AbstractControlModel, var"#s1793"<:CTModels.OCP.AbstractVariableModel, var"#s1794"<:Function, var"#s1795"<:CTModels.OCP.AbstractObjectiveModel, TP<:Tuple{Vector{Float64}, Any, Vector{Float64}, Vector{Symbol}}, TB<:Tuple{Vector{Float64}, Any, Vector{Float64}, Vector{Symbol}}}

Converts a mutable PreModel into an immutable Model.

This function finalizes a pre-defined optimal control problem (PreModel) by verifying that all necessary components (times, state, dynamics, objective) are set. It then constructs a Model instance, incorporating optional components like control, variable, and constraints.

Note

Control is optional: calling control! is not required. When omitted, the model is built with control_dimension == 0 (an EmptyControlModel). This is useful for problems where the dynamics depend only on the state, such as pure state-space systems.

Arguments

  • pre_ocp::PreModel: The pre-defined optimal control problem to be finalized.

Returns

  • Model: A fully constructed model ready for solving.

Example without control

julia> pre_ocp = PreModel()
julia> times!(pre_ocp, 0.0, 1.0, 100)
julia> state!(pre_ocp, 2, "x", ["x1", "x2"])
julia> dynamics!(pre_ocp, (t, x, u) -> [-x[2], x[1]])
julia> objective!(pre_ocp, :min, mayer=(x0, xf) -> xf[1]^2)
julia> model = build(pre_ocp)
julia> control_dimension(model)  # 0

Example with control

julia> pre_ocp = PreModel()
julia> times!(pre_ocp, 0.0, 1.0, 100)
julia> state!(pre_ocp, 2, "x", ["x1", "x2"])
julia> control!(pre_ocp, 1, "u", ["u1"])
julia> dynamics!(pre_ocp, (dx, t, x, u, v) -> dx .= x + u)
julia> model = build(pre_ocp)

build_model

CTModels.OCP.build_modelFunction
build_model(
    pre_ocp::CTModels.OCP.PreModel;
    build_examodel
) -> CTModels.OCP.Model{TD, var"#s179", var"#s1791", var"#s1792", var"#s1793", var"#s1794", var"#s1795", CTModels.OCP.ConstraintsModel{TP, TB, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}, Tuple{Vector{Float64}, Vector{Int64}, Vector{Float64}, Vector{Symbol}}}, Nothing} where {TD<:CTModels.OCP.TimeDependence, var"#s179"<:CTModels.OCP.AbstractTimesModel, var"#s1791"<:CTModels.OCP.AbstractStateModel, var"#s1792"<:CTModels.OCP.AbstractControlModel, var"#s1793"<:CTModels.OCP.AbstractVariableModel, var"#s1794"<:Function, var"#s1795"<:CTModels.OCP.AbstractObjectiveModel, TP<:Tuple{Vector{Float64}, Any, Vector{Float64}, Vector{Symbol}}, TB<:Tuple{Vector{Float64}, Any, Vector{Float64}, Vector{Symbol}}}

Build a complete optimal control problem model from a pre-model.

This function is an alias for build(pre_ocp; build_examodel=build_examodel) and constructs a fully validated Model from a PreModel by extracting and organizing all components (times, state, control, variable, dynamics, objective, constraints).

Arguments

  • pre_ocp::PreModel: The pre-model containing all problem components
  • build_examodel=nothing: Optional ExaModel builder function for GPU acceleration

Returns

  • Model: A complete, validated optimal control problem model

Throws

  • Exceptions.PreconditionError: If time dependence has not been set via time_dependence!

Example

using CTModels

# Create and configure a pre-model
pre_ocp = PreModel()
time_dependence!(pre_ocp, autonomous=true)
state!(pre_ocp, 2)
control!(pre_ocp, 1)
dynamics!(pre_ocp, (x, u) -> [x[2], u[1]])
objective!(pre_ocp, :mayer, (x0, xf) -> xf[1]^2)

# Build the model
ocp = build_model(pre_ocp)

See also: build, PreModel, Model, time_dependence!

build_solution

CTModels.OCP.build_solutionFunction
build_solution(
    ocp::CTModels.OCP.Model,
    T_state::Vector{Float64},
    T_control::Vector{Float64},
    T_costate::Vector{Float64},
    T_path::Union{Nothing, Vector{Float64}},
    X::Union{Function, Matrix{Float64}},
    U::Union{Function, Matrix{Float64}},
    v::Vector{Float64},
    P::Union{Function, Matrix{Float64}};
    objective,
    iterations,
    constraints_violation,
    message,
    status,
    successful,
    path_constraints_dual,
    boundary_constraints_dual,
    state_constraints_lb_dual,
    state_constraints_ub_dual,
    control_constraints_lb_dual,
    control_constraints_ub_dual,
    variable_constraints_lb_dual,
    variable_constraints_ub_dual,
    infos
)

Build a solution from an optimal control problem with independent time grids for each component.

This function constructs a Solution object by assembling trajectory data (state, control, costate, path constraint duals) defined on potentially different time discretizations. The solution automatically creates interpolated functions to evaluate trajectories at arbitrary time points, and optimizes storage when all grids are identical.

Time Grid Semantics

The solution supports four independent time grids, each associated with a specific trajectory component:

  • T_state: Time grid for the state trajectory X and state box constraint duals
    • Defines discretization points where state values are known
    • State box constraint duals (state_constraints_lb_dual, state_constraints_ub_dual) share this grid
  • T_control: Time grid for the control trajectory U and control box constraint duals
    • Defines discretization points where control values are known
    • Control box constraint duals (control_constraints_lb_dual, control_constraints_ub_dual) share this grid
    • May differ from T_state (e.g., coarser discretization for piecewise constant controls)
  • T_costate: Time grid for the costate (adjoint) trajectory P
    • Defines discretization points where costate values are known
    • Independent from state grid to accommodate different numerical schemes
    • Example: symplectic integrators may use different grids for state and costate
  • T_path: Time grid for path constraint duals (can be nothing)
    • Defines discretization points for path constraint dual variables
    • Set to nothing if no path constraints exist
    • When nothing, internally defaults to T_state for consistency

Grid Optimization: If all non-nothing grids are identical, the solution uses UnifiedTimeGridModel for memory efficiency. Otherwise, it uses MultipleTimeGridModel to store each grid separately.

Trajectory Data Formats

Trajectory data (X, U, P, path_constraints_dual) can be provided in two formats:

  1. Matrix format: Matrix{Float64} with dimensions (n_points, n_dim)

    • Each row corresponds to a time point in the associated grid
    • Each column corresponds to a component dimension
    • Example: X is (length(T_state), state_dimension(ocp))
  2. Function format: Function that takes time t::Float64 and returns a vector

    • Allows analytical or pre-interpolated trajectories
    • Function signature: t -> Vector{Float64} of appropriate dimension
    • Useful for exact solutions or when data is already interpolated

Arguments

Required Positional Arguments

  • ocp::Model: The optimal control problem model defining dimensions and structure
  • T_state::Vector{Float64}: Time grid for state trajectory (must be strictly increasing)
  • T_control::Vector{Float64}: Time grid for control trajectory (must be strictly increasing)
  • T_costate::Vector{Float64}: Time grid for costate trajectory (must be strictly increasing)
  • T_path::Union{Vector{Float64},Nothing}: Time grid for path constraint duals (or nothing)
  • X::Union{Matrix{Float64},Function}: State trajectory data
  • U::Union{Matrix{Float64},Function}: Control trajectory data
  • v::Vector{Float64}: Variable values (static optimization variables, not time-dependent)
  • P::Union{Matrix{Float64},Function}: Costate (adjoint) trajectory data

Required Keyword Arguments

  • objective::Float64: Optimal objective function value
  • iterations::Int: Number of solver iterations performed
  • constraints_violation::Float64: Maximum constraint violation (feasibility measure)
  • message::String: Solver status message (e.g., "SolveSucceeded", "IterationLimit")
  • status::Symbol: Solver termination status (e.g., :Solve_Succeeded, :Iteration_Limit)
  • successful::Bool: Whether the solve was successful (true/false)

Optional Keyword Arguments (Dual Variables)

All dual variable arguments default to nothing if not provided:

  • path_constraints_dual::Union{Matrix{Float64},Function,Nothing}: Path constraint duals on T_path grid
  • boundary_constraints_dual::Union{Vector{Float64},Nothing}: Boundary constraint duals (time-independent)
  • state_constraints_lb_dual::Union{Matrix{Float64},Nothing}: State lower bound duals on T_state grid
  • state_constraints_ub_dual::Union{Matrix{Float64},Nothing}: State upper bound duals on T_state grid
  • control_constraints_lb_dual::Union{Matrix{Float64},Nothing}: Control lower bound duals on T_control grid
  • control_constraints_ub_dual::Union{Matrix{Float64},Nothing}: Control upper bound duals on T_control grid
  • variable_constraints_lb_dual::Union{Vector{Float64},Nothing}: Variable lower bound duals (time-independent)
  • variable_constraints_ub_dual::Union{Vector{Float64},Nothing}: Variable upper bound duals (time-independent)
  • infos::Dict{Symbol,Any}: Additional solver-specific information (default: empty dict)

Returns

  • sol::Solution: Complete solution object with interpolated trajectory functions and metadata

Example

using CTModels

# Build OCP
ocp = Model(...)
state!(ocp, 2)
control!(ocp, 1)
# ... define dynamics, objective, etc.

# Define independent time grids
T_state = collect(LinRange(0.0, 1.0, 101))    # Fine state grid (101 points)
T_control = collect(LinRange(0.0, 1.0, 51))   # Coarser control grid (51 points)
T_costate = collect(LinRange(0.0, 1.0, 76))   # Custom costate grid (76 points)
T_path = collect(LinRange(0.0, 1.0, 61))      # Path constraint grid (61 points)

# Trajectory data (matrix format)
X = rand(101, 2)  # State on T_state grid
U = rand(51, 1)   # Control on T_control grid
P = rand(76, 2)   # Costate on T_costate grid
v = [0.5, 1.2]    # Static variables

# Build solution
sol = build_solution(
    ocp,
    T_state, T_control, T_costate, T_path,
    X, U, v, P;
    objective=1.23,
    iterations=50,
    constraints_violation=1e-8,
    message="Optimal",
    status=:first_order,
    successful=true
)

# Access trajectories (automatically interpolated)
x_at_t = state(sol)(0.5)      # Interpolated from T_state grid
u_at_t = control(sol)(0.5)    # Interpolated from T_control grid
p_at_t = costate(sol)(0.5)    # Interpolated from T_costate grid

# Query time grids
time_grid(sol, :state)    # Returns T_state
time_grid(sol, :control)  # Returns T_control
time_grid(sol, :costate)  # Returns T_costate

Notes

Box Constraint Dual Dimensions

The dimensions of box constraint dual variables correspond to the component dimension, not the number of constraint declarations:

  • state_constraints_*_dual: Dimension (length(T_state), state_dimension(ocp))
  • control_constraints_*_dual: Dimension (length(T_control), control_dimension(ocp))
  • variable_constraints_*_dual: Dimension variable_dimension(ocp)

If multiple constraints are declared on the same component (e.g., x₂(t) ≤ 1.2 and x₂(t) ≤ 2.0), only the last bound value is retained, and a warning is emitted during model construction.

Grid Validation

All time grids must be:

  • Strictly increasing: T[i] < T[i+1] for all i
  • Non-empty: At least one time point
  • Finite: No Inf or NaN values

The function automatically validates and fixes grids (e.g., converts ranges to vectors).

Memory Optimization

When all grids are identical, the solution uses UnifiedTimeGridModel to store a single grid, reducing memory overhead. This is detected automatically.

Backward Compatibility

A legacy signature build_solution(ocp, T, X, U, v, P; ...) exists for single-grid solutions. It internally calls this multi-grid version with T_state = T_control = T_costate = T_path = T.

See also: Solution, UnifiedTimeGridModel, MultipleTimeGridModel, time_grid, state, control, costate

build_solution(
    ocp::CTModels.OCP.Model,
    T::Vector{Float64},
    X::Union{Function, Matrix{Float64}},
    U::Union{Function, Matrix{Float64}},
    v::Vector{Float64},
    P::Union{Function, Matrix{Float64}};
    objective,
    iterations,
    constraints_violation,
    message,
    status,
    successful,
    path_constraints_dual,
    boundary_constraints_dual,
    state_constraints_lb_dual,
    state_constraints_ub_dual,
    control_constraints_lb_dual,
    control_constraints_ub_dual,
    variable_constraints_lb_dual,
    variable_constraints_ub_dual,
    infos
)

Build a solution from the optimal control problem, the time grid, the state, control, variable, and dual variables.

Arguments

  • ocp::Model: the optimal control problem.
  • T::Vector{Float64}: the time grid.
  • X::Matrix{Float64}: the state trajectory.
  • U::Matrix{Float64}: the control trajectory.
  • v::Vector{Float64}: the variable trajectory.
  • P::Matrix{Float64}: the costate trajectory.
  • objective::Float64: the objective value.
  • iterations::Int: the number of iterations.
  • constraints_violation::Float64: the constraints violation.
  • message::String: the message associated to the status criterion.
  • status::Symbol: the status criterion.
  • successful::Bool: the successful status.
  • path_constraints_dual::Matrix{Float64}: the dual of the path constraints.
  • boundary_constraints_dual::Vector{Float64}: the dual of the boundary constraints.
  • state_constraints_lb_dual::Matrix{Float64}: the lower bound dual of the state constraints.
  • state_constraints_ub_dual::Matrix{Float64}: the upper bound dual of the state constraints.
  • control_constraints_lb_dual::Matrix{Float64}: the lower bound dual of the control constraints.
  • control_constraints_ub_dual::Matrix{Float64}: the upper bound dual of the control constraints.
  • variable_constraints_lb_dual::Vector{Float64}: the lower bound dual of the variable constraints.
  • variable_constraints_ub_dual::Vector{Float64}: the upper bound dual of the variable constraints.
  • infos::Dict{Symbol,Any}: additional solver information dictionary.

Returns

  • sol::Solution: the optimal control solution.

Notes

The dimensions of box constraint dual variables (state_constraints_*_dual, control_constraints_*_dual, variable_constraints_*_dual) correspond to the state/control/variable dimension, not the number of constraint declarations. If multiple constraints are declared on the same component (e.g., x₂(t) ≤ 1.2 and x₂(t) ≤ 2.0), only the last bound value is retained, and a warning is emitted during model construction.

constraints

CTModels.OCP.constraintsFunction
constraints(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.AbstractObjectiveModel, C<:CTModels.OCP.AbstractConstraintsModel}
) -> CTModels.OCP.AbstractConstraintsModel

Return the constraints struct.

constraints_violation

control

CTModels.OCP.controlFunction
control(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel, <:CTModels.OCP.AbstractStateModel, T<:CTModels.OCP.AbstractControlModel}
) -> CTModels.OCP.AbstractControlModel

Return the control struct.

control(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.ControlModelSolution{TS<:Function}}
) -> Function

Return the control as a function of time.

julia> u  = control(sol)
julia> t0 = time_grid(sol)[1]
julia> u0 = u(t0) # control at the initial time
control(init::CTModels.Init.AbstractInitialGuess) -> Any

Return the control trajectory from an initial guess.

control(sol::CTModels.OCP.AbstractSolution) -> Function

Return the control trajectory from a solution.

control_components

CTModels.OCP.control_componentsFunction
control_components(
    ocp::CTModels.OCP.Model
) -> Vector{String}

Return the names of the components of the control.

control_components(
    sol::CTModels.OCP.Solution
) -> Vector{String}

Return the names of the components of the control.

control_constraints_lb_dual

CTModels.OCP.control_constraints_lb_dualFunction
control_constraints_lb_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, CC_LB_Dual<:Union{Nothing, Function}}
) -> Union{Nothing, Function}

Return the dual function associated with the lower bounds of control constraints.

Arguments

  • model::DualModel: A model including dual variables for control lower bounds.

Returns

A function mapping time t to a vector of dual values, or nothing if not set.

control_constraints_lb_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, Function}

Return the lower bound dual of the control constraints.

control_constraints_ub_dual

CTModels.OCP.control_constraints_ub_dualFunction
control_constraints_ub_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, CC_UB_Dual<:Union{Nothing, Function}}
) -> Union{Nothing, Function}

Return the dual function associated with the upper bounds of control constraints.

Arguments

  • model::DualModel: A model including dual variables for control upper bounds.

Returns

A function mapping time t to a vector of dual values, or nothing if not set.

control_constraints_ub_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, Function}

Return the upper bound dual of the control constraints.

control_dimension

CTModels.OCP.control_dimensionFunction
control_dimension(ocp::CTModels.OCP.Model) -> Int64

Return the control dimension.

control_dimension(sol::CTModels.OCP.Solution) -> Int64

Return the dimension of the control.

control_name

CTModels.OCP.control_nameFunction
control_name(ocp::CTModels.OCP.Model) -> String

Return the name of the control.

control_name(sol::CTModels.OCP.Solution) -> String

Return the name of the control.

costate

CTModels.OCP.costateFunction
costate(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:CTModels.OCP.AbstractModel, Co<:Function}
) -> Function

Return the costate as a function of time.

julia> p  = costate(sol)
julia> t0 = time_grid(sol)[1]
julia> p0 = p(t0) # costate at the initial time

definition

CTModels.OCP.definitionFunction
definition(ocp::CTModels.OCP.Model) -> Expr

Return the model definition of the optimal control problem.

Arguments

  • ocp::Model: The built optimal control problem model.

Returns

  • Expr: The symbolic expression defining the problem.
definition(
    ocp::CTModels.OCP.PreModel
) -> Union{Nothing, Expr}

Return the model definition of the optimal control problem or nothing.

Arguments

  • ocp::PreModel: The pre-model (may not have a definition set).

Returns

  • Union{Expr, Nothing}: The symbolic expression or nothing if not set.

definition!

CTModels.OCP.definition!Function
definition!(ocp::CTModels.OCP.PreModel, definition::Expr)

Set the model definition of the optimal control problem.

Arguments

  • ocp::PreModel: The pre-model to modify.
  • definition::Expr: The symbolic expression defining the problem.

Returns

  • Nothing

dim_boundary_constraints_nl

CTModels.OCP.dim_boundary_constraints_nlFunction
dim_boundary_constraints_nl(
    model::CTModels.OCP.ConstraintsModel
) -> Int64

Return the dimension of nonlinear boundary constraints.

Arguments

  • model: The constraints model from which to retrieve the dimension of boundary constraints.

Returns

  • Dimension: The dimension of the nonlinear boundary constraints.

Example

# Example of getting the dimension of nonlinear boundary constraints
julia> model = ConstraintsModel(...)
julia> dim_boundary = dim_boundary_constraints_nl(model)
dim_boundary_constraints_nl(
    ocp::CTModels.OCP.Model
) -> Int64

Return the dimension of the boundary constraints.

dim_boundary_constraints_nl(
    sol::CTModels.OCP.Solution
) -> Int64

Return the dimension of the boundary constraints.

dim_control_constraints_box

CTModels.OCP.dim_control_constraints_boxFunction
dim_control_constraints_box(
    model::CTModels.OCP.ConstraintsModel
) -> Int64

Return the dimension of control box constraints.

Arguments

  • model: The constraints model from which to retrieve the dimension of control box constraints.

Returns

  • Dimension: The dimension of the control box constraints.

Example

julia> # Example of getting the dimension of control box constraints
julia> model = ConstraintsModel(...)
julia> dim_control = dim_control_constraints_box(model)
dim_control_constraints_box(
    ocp::CTModels.OCP.Model
) -> Int64

Return the dimension of box constraints on control.

dim_control_constraints_box(
    sol::CTModels.OCP.Solution
) -> Int64

Return the dimension of box constraints on control.

dim_path_constraints_nl

CTModels.OCP.dim_path_constraints_nlFunction
dim_path_constraints_nl(
    model::CTModels.OCP.ConstraintsModel
) -> Int64

Return the dimension of nonlinear path constraints.

Arguments

  • model: The constraints model from which to retrieve the dimension of path constraints.

Returns

  • Dimension: The dimension of the nonlinear path constraints.

Example

# Example of getting the dimension of nonlinear path constraints
julia> model = ConstraintsModel(...)
julia> dim_path = dim_path_constraints_nl(model)
dim_path_constraints_nl(ocp::CTModels.OCP.Model) -> Int64

Return the dimension of nonlinear path constraints.

dim_path_constraints_nl(sol::CTModels.OCP.Solution) -> Int64

Return the dimension of the path constraints.

dim_state_constraints_box

CTModels.OCP.dim_state_constraints_boxFunction
dim_state_constraints_box(
    model::CTModels.OCP.ConstraintsModel
) -> Int64

Return the dimension of state box constraints.

Arguments

  • model: The constraints model from which to retrieve the dimension of state box constraints.

Returns

  • Dimension: The dimension of the state box constraints.

Example

julia> # Example of getting the dimension of state box constraints
julia> model = ConstraintsModel(...)
julia> dim_state = dim_state_constraints_box(model)
dim_state_constraints_box(ocp::CTModels.OCP.Model) -> Int64

Return the dimension of box constraints on state.

dim_state_constraints_box(
    sol::CTModels.OCP.Solution
) -> Int64

Return the dimension of box constraints on state.

dim_variable_constraints_box

CTModels.OCP.dim_variable_constraints_boxFunction
dim_variable_constraints_box(
    model::CTModels.OCP.ConstraintsModel
) -> Int64

Return the dimension of variable box constraints.

Arguments

  • model: The constraints model from which to retrieve the dimension of variable box constraints.

Returns

  • Dimension: The dimension of the variable box constraints.

Example

julia> # Example of getting the dimension of variable box constraints
julia> model = ConstraintsModel(...)
julia> dim_variable = dim_variable_constraints_box(model)
dim_variable_constraints_box(
    ocp::CTModels.OCP.Model
) -> Int64

Return the dimension of box constraints on variable.

dim_variable_constraints_box(
    sol::CTModels.OCP.Solution
) -> Int64

Return the dimension of the variable box constraints.

dual

CTModels.OCP.dualFunction
dual(
    sol::CTModels.OCP.Solution,
    model::CTModels.OCP.Model,
    label::Symbol
) -> Any

Return the dual variable associated with a constraint identified by its label.

Searches through all constraint types (path, boundary, state, control, and variable constraints) defined in the model and returns the corresponding dual value from the solution.

Arguments

  • sol::Solution: Solution object containing dual variables.
  • model::Model: Model containing constraint definitions.
  • label::Symbol: Symbol corresponding to a constraint label.

Returns

A function of time t for time-dependent constraints, or a scalar/vector for time-invariant duals. If the label is not found, throws an IncorrectArgument exception.

dynamics

CTModels.OCP.dynamicsFunction
dynamics(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, D<:Function}
) -> Function

Return the dynamics.

final_time

CTModels.OCP.final_timeFunction
final_time(
    model::CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, <:CTModels.OCP.FixedTimeModel{T<:Real}}
) -> Real

Get the final time from the times model, from a fixed final time model.

final_time(
    model::CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, CTModels.OCP.FreeTimeModel},
    variable::AbstractArray{T<:Real, 1}
) -> Any

Get the final time from the times model, from a free final time model.

final_time(ocp::CTModels.OCP.AbstractModel) -> Any

Throw an error for unsupported final time access.

final_time(
    ocp::CTModels.OCP.AbstractModel,
    variable::AbstractVector
) -> Any

Throw an error for unsupported final time access with variable.

final_time(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, CTModels.OCP.FixedTimeModel{T<:Real}}}
) -> Any

Return the final time, for a fixed final time.

final_time(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, CTModels.OCP.FreeTimeModel}},
    variable::AbstractArray{T<:Real, 1}
) -> Any

Return the final time, for a free final time.

final_time(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, CTModels.OCP.FreeTimeModel}},
    variable::Real
) -> Real

Return the final time, for a free final time.

final_time(sol::CTModels.OCP.Solution) -> Real

Return the final time of the solution.

final_time_name

CTModels.OCP.final_time_nameFunction
final_time_name(model::CTModels.OCP.TimesModel) -> String

Get the name of the final time from the times model.

final_time_name(ocp::CTModels.OCP.Model) -> String

Return the name of the final time.

final_time_name(sol::CTModels.OCP.Solution) -> String

Return the name of the final time.

get_build_examodel

CTModels.OCP.get_build_examodelFunction
get_build_examodel(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.AbstractObjectiveModel, <:CTModels.OCP.AbstractConstraintsModel, BE<:Function}
) -> Function

Return the build_examodel.

get_build_examodel(
    _::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.AbstractObjectiveModel, <:CTModels.OCP.AbstractConstraintsModel, <:Nothing}
)

Return an error (PreconditionError) since the model is not built with the :exa backend.

has_fixed_final_time

CTModels.OCP.has_fixed_final_timeFunction
has_fixed_final_time(
    times::CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, <:CTModels.OCP.FixedTimeModel{T<:Real}}
) -> Bool

Check if the final time is fixed. Return true.

has_fixed_final_time(
    times::CTModels.OCP.TimesModel{<:CTModels.OCP.AbstractTimeModel, CTModels.OCP.FreeTimeModel}
) -> Bool

Check if the final time is free. Return false.

has_fixed_final_time(ocp::CTModels.OCP.Model) -> Bool

Check if the final time is fixed.

has_fixed_final_time(sol::CTModels.OCP.Solution) -> Bool

Check if the final time is fixed.

has_free_final_time

CTModels.OCP.has_free_final_timeFunction
has_free_final_time(times::CTModels.OCP.TimesModel) -> Bool

Check if the final time is free.

has_free_final_time(ocp::CTModels.OCP.Model) -> Bool

Check if the final time is free.

has_free_final_time(sol::CTModels.OCP.Solution) -> Bool

Check if the final time is free.

has_free_initial_time

CTModels.OCP.has_free_initial_timeFunction
has_free_initial_time(
    times::CTModels.OCP.TimesModel
) -> Bool

Check if the final time is free.

has_free_initial_time(ocp::CTModels.OCP.Model) -> Bool

Check if the initial time is free.

has_free_initial_time(sol::CTModels.OCP.Solution) -> Bool

Check if the initial time is free.

infos

CTModels.OCP.infosFunction
infos(sol::CTModels.OCP.Solution) -> Dict{Symbol, Any}

Return a dictionary of additional infos depending on the solver or nothing.

initial_time_name

CTModels.OCP.initial_time_nameFunction
initial_time_name(model::CTModels.OCP.TimesModel) -> String

Get the name of the initial time from the times model.

initial_time_name(ocp::CTModels.OCP.Model) -> String

Return the name of the initial time.

initial_time_name(sol::CTModels.OCP.Solution) -> String

Return the name of the initial time.

iterations

CTModels.OCP.iterationsFunction
iterations(sol::CTModels.OCP.Solution) -> Int64

Return the number of iterations (if solved by an iterative method).

lagrange

CTModels.OCP.lagrangeFunction
lagrange(
    model::CTModels.OCP.LagrangeObjectiveModel{L<:Function}
) -> Function

Return the Lagrange function.

lagrange(
    model::CTModels.OCP.BolzaObjectiveModel{<:Function, L<:Function}
) -> Function

Return the Lagrange function.

lagrange(ocp::CTModels.OCP.AbstractModel) -> Function

Throw an error when accessing Lagrange cost on a model without one.

lagrange(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, CTModels.OCP.LagrangeObjectiveModel{L<:Function}}
) -> Function

Return the Lagrange cost.

lagrange(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.BolzaObjectiveModel{<:Function, L<:Function}}
) -> Any

Return the Lagrange cost.

mayer

CTModels.OCP.mayerFunction
mayer(
    model::CTModels.OCP.MayerObjectiveModel{M<:Function}
) -> Function

Return the Mayer function.

mayer(
    model::CTModels.OCP.BolzaObjectiveModel{M<:Function}
) -> Function

Return the Mayer function.

mayer(ocp::CTModels.OCP.AbstractModel) -> Any

Throw an error when accessing Mayer cost on a model without one.

mayer(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.MayerObjectiveModel{M<:Function}}
) -> Any

Return the Mayer cost.

mayer(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.BolzaObjectiveModel{M<:Function}}
) -> Any

Return the Mayer cost.

message

CTModels.OCP.messageFunction
message(sol::CTModels.OCP.Solution) -> String

Return the message associated to the status criterion.

model

CTModels.OCP.modelFunction
model(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, M<:CTModels.OCP.AbstractModel}
) -> CTModels.OCP.AbstractModel

Return the model of the optimal control problem.

objective

CTModels.OCP.objectiveFunction
objective(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, O<:CTModels.OCP.AbstractObjectiveModel}
) -> CTModels.OCP.AbstractObjectiveModel

Return the objective struct.

objective(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:CTModels.OCP.AbstractModel, <:Function, O<:Real}
) -> Real

Return the objective value.

path_constraints_dual

CTModels.OCP.path_constraints_dualFunction
path_constraints_dual(
    model::CTModels.OCP.DualModel{PC_Dual<:Union{Nothing, Function}}
) -> Union{Nothing, Function}

Return the dual function associated with the nonlinear path constraints.

Arguments

  • model::DualModel: A model including dual variables for path constraints.

Returns

A function mapping time t to the vector of dual values, or nothing if not set.

path_constraints_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, Function}

Return the dual of the path constraints.

state

CTModels.OCP.stateFunction
state(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel, T<:CTModels.OCP.AbstractStateModel}
) -> CTModels.OCP.AbstractStateModel

Return the state struct.

state(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.StateModelSolution{TS<:Function}}
) -> Function

Return the state as a function of time.

julia> x  = state(sol)
julia> t0 = time_grid(sol)[1]
julia> x0 = x(t0) # state at the initial time
state(init::CTModels.Init.AbstractInitialGuess) -> Any

Return the state trajectory from an initial guess.

state(sol::CTModels.OCP.AbstractSolution) -> Function

Return the state trajectory from a solution.

state_components

CTModels.OCP.state_componentsFunction
state_components(ocp::CTModels.OCP.Model) -> Vector{String}

Return the names of the components of the state.

state_components(
    sol::CTModels.OCP.Solution
) -> Vector{String}

Return the names of the components of the state.

state_constraints_box

CTModels.OCP.state_constraints_boxFunction
state_constraints_box(
    model::CTModels.OCP.ConstraintsModel{<:Tuple, <:Tuple, TS}
) -> Any

Get the state box constraints from the model.

Arguments

  • model: The constraints model from which to retrieve the state box constraints.

Returns

  • The state box constraints.

Example

# Example of retrieving state box constraints
julia> model = ConstraintsModel(...)
julia> state_constraints = state_constraints_box(model)
state_constraints_box(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.AbstractVariableModel, <:Function, <:CTModels.OCP.AbstractObjectiveModel, <:CTModels.OCP.ConstraintsModel{<:Tuple, <:Tuple, TS<:Tuple}}
) -> Any

Return the box constraints on state.

state_constraints_lb_dual

CTModels.OCP.state_constraints_lb_dualFunction
state_constraints_lb_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, SC_LB_Dual<:Union{Nothing, Function}}
) -> Union{Nothing, Function}

Return the dual function associated with the lower bounds of state constraints.

Arguments

  • model::DualModel: A model including dual variables for state lower bounds.

Returns

A function mapping time t to a vector of dual values, or nothing if not set.

state_constraints_lb_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, Function}

Return the lower bound dual of the state constraints.

state_constraints_ub_dual

CTModels.OCP.state_constraints_ub_dualFunction
state_constraints_ub_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, <:Union{Nothing, Function}, SC_UB_Dual<:Union{Nothing, Function}}
) -> Union{Nothing, Function}

Return the dual function associated with the upper bounds of state constraints.

Arguments

  • model::DualModel: A model including dual variables for state upper bounds.

Returns

A function mapping time t to a vector of dual values, or nothing if not set.

state_constraints_ub_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, Function}

Return the upper bound dual of the state constraints.

state_name

CTModels.OCP.state_nameFunction
state_name(ocp::CTModels.OCP.Model) -> String

Return the name of the state.

state_name(sol::CTModels.OCP.Solution) -> String

Return the name of the state.

status

CTModels.OCP.statusFunction
status(sol::CTModels.OCP.Solution) -> Symbol

Return the status criterion (a Symbol).

successful

time_grid

CTModels.OCP.time_gridFunction
time_grid(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.UnifiedTimeGridModel{T<:Union{StepRangeLen, AbstractVector{<:Real}}}}
) -> Union{StepRangeLen, AbstractVector{<:Real}}

Return the time grid for solutions with unified time grid.

time_grid(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.UnifiedTimeGridModel{T<:Union{StepRangeLen, AbstractVector{<:Real}}}},
    component::Symbol
) -> Union{StepRangeLen, AbstractVector{<:Real}}

Return the time grid for a specific component.

Arguments

  • sol::Solution: The solution (unified or multiple time grids)
  • component::Symbol: The component (:state, :control, :path) Also accepted: :costate/:costates (→ :state), :dual/:duals (→ :path), :stateboxconstraint(s) (→ :state), :controlboxconstraint(s) (→ :control), plural forms (:states, :controls)

Returns

  • TimesDisc: The time grid for the specified component

Behavior

  • For UnifiedTimeGridModel: Returns the unique time grid for any component
  • For MultipleTimeGridModel: Returns the specific time grid for the component

Throws

  • IncorrectArgument: If component is not one of the valid symbols

Examples

julia> time_grid(sol, :state)   # Works for both unified and multiple grids
julia> time_grid(sol, :control) # Works for both unified and multiple grids
julia> time_grid(sol, :costate) # Maps to :state grid
julia> time_grid(sol, :dual)    # Maps to :path grid
time_grid(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.MultipleTimeGridModel},
    component::Symbol
) -> Union{StepRangeLen, AbstractVector{<:Real}}

Return the time grid for a specific component in solutions with multiple time grids.

Arguments

  • sol::Solution: The solution with multiple time grids
  • component::Symbol: The component (:state, :control, :path) Also accepted: :costate/:costates (→ :state), :dual/:duals (→ :path), :stateboxconstraint(s) (→ :state), :controlboxconstraint(s) (→ :control), plural forms (:states, :controls)

Returns

  • TimesDisc: The time grid for the specified component

Throws

  • IncorrectArgument: If component is not one of the valid symbols

Examples

julia> time_grid(sol, :state)   # Get state time grid
julia> time_grid(sol, :control) # Get control time grid
julia> time_grid(sol, :costate) # Maps to state time grid
julia> time_grid(sol, :dual)    # Maps to path time grid
time_grid(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.MultipleTimeGridModel}
)

Return the time grid for solutions with multiple time grids (component must be specified).

Throws

  • IncorrectArgument: Always thrown for MultipleTimeGridModel without component specification

Notes

This method enforces explicit component specification for solutions with multiple time grids to avoid ambiguity about which grid is being accessed.

Examples

julia> time_grid(sol)  # ❌ Error for MultipleTimeGridModel
julia> time_grid(sol, :state)  # ✅ Correct usage

time_name

CTModels.OCP.time_nameFunction
time_name(model::CTModels.OCP.TimesModel) -> String

Get the name of the time variable from the times model.

time_name(ocp::CTModels.OCP.Model) -> String

Return the name of the time.

time_name(sol::CTModels.OCP.Solution) -> String

Return the name of the time component.

times

CTModels.OCP.timesFunction
times(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, T<:CTModels.OCP.TimesModel}
) -> CTModels.OCP.TimesModel

Return the times struct.

times(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, TM<:CTModels.OCP.AbstractTimesModel}
) -> CTModels.OCP.AbstractTimesModel

Return the times model.

variable

CTModels.OCP.variableFunction
variable(
    ocp::CTModels.OCP.Model{<:CTModels.OCP.TimeDependence, <:CTModels.OCP.TimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, T<:CTModels.OCP.AbstractVariableModel}
) -> CTModels.OCP.AbstractVariableModel

Return the variable struct.

variable(
    sol::CTModels.OCP.Solution{<:CTModels.OCP.AbstractTimeGridModel, <:CTModels.OCP.AbstractTimesModel, <:CTModels.OCP.AbstractStateModel, <:CTModels.OCP.AbstractControlModel, <:CTModels.OCP.VariableModelSolution{TS<:Union{Real, AbstractVector{<:Real}}}}
) -> Union{Real, AbstractVector{<:Real}}

Return the variable or nothing.

julia> v  = variable(sol)
variable(init::CTModels.Init.AbstractInitialGuess) -> Any

Return the variable value from an initial guess.

variable_components

CTModels.OCP.variable_componentsFunction
variable_components(
    ocp::CTModels.OCP.Model
) -> Vector{String}

Return the names of the components of the variable.

variable_components(
    sol::CTModels.OCP.Solution
) -> Vector{String}

Return the names of the components of the variable.

variable_constraints_lb_dual

CTModels.OCP.variable_constraints_lb_dualFunction
variable_constraints_lb_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, VC_LB_Dual<:Union{Nothing, AbstractVector{<:Real}}}
) -> Union{Nothing, AbstractVector{<:Real}}

Return the dual vector associated with the lower bounds of variable constraints.

Arguments

  • model::DualModel: A model including dual variables for variable lower bounds.

Returns

A vector of dual values, or nothing if not set.

variable_constraints_lb_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, AbstractVector{<:Real}}

Return the lower bound dual of the variable constraints.

variable_constraints_ub_dual

CTModels.OCP.variable_constraints_ub_dualFunction
variable_constraints_ub_dual(
    model::CTModels.OCP.DualModel{<:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, <:Union{Nothing, Function}, <:Union{Nothing, AbstractVector{<:Real}}, VC_UB_Dual<:Union{Nothing, AbstractVector{<:Real}}}
) -> Union{Nothing, AbstractVector{<:Real}}

Return the dual vector associated with the upper bounds of variable constraints.

Arguments

  • model::DualModel: A model including dual variables for variable upper bounds.

Returns

A vector of dual values, or nothing if not set.

variable_constraints_ub_dual(
    sol::CTModels.OCP.Solution
) -> Union{Nothing, AbstractVector{<:Real}}

Return the upper bound dual of the variable constraints.

variable_dimension

CTModels.OCP.variable_dimensionFunction
variable_dimension(ocp::CTModels.OCP.Model) -> Int64

Return the variable dimension.

variable_dimension(sol::CTModels.OCP.Solution) -> Int64

Return the dimension of the variable.

variable_name

CTModels.OCP.variable_nameFunction
variable_name(ocp::CTModels.OCP.Model) -> String

Return the name of the variable.

variable_name(sol::CTModels.OCP.Solution) -> String

Return the name of the variable.