Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
92 changes: 36 additions & 56 deletions docs/problems.jl
Original file line number Diff line number Diff line change
Expand Up @@ -28,15 +28,6 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
DRAFT = draft_meta(draft)
LEFT_MARGIN = get_left_margin(Symbol(PROBLEM))

VARIABLE_COMPONENTS = isnothing(OptimalControlProblems.metadata(Symbol(PROBLEM))[:variable_components]) ? "" :
"""
The variable components are named:

```@example main
metadata(:$PROBLEM)[:variable_components]
```
"""

documentation=DRAFT * """
# $TITLE

Expand Down Expand Up @@ -78,20 +69,12 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni

## Metadata

The state components are named:

```@example main
metadata(:$PROBLEM)[:state_components]
```

The control components are named:
The default number of time steps is:

```@example main
metadata(:$PROBLEM)[:control_components]
metadata(:$PROBLEM)[:grid_size]
```

$VARIABLE_COMPONENTS

The default values of the parameters are:

```@example main
Expand Down Expand Up @@ -120,25 +103,22 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
```@example main
function plot_initial_guess(problem)

# -----------------------------
# Extract dimensions from metadata
# -----------------------------
x_vars = metadata(problem)[:state_components]
u_vars = metadata(problem)[:control_components]
n_states = length(x_vars)
n_controls = length(u_vars)

# -----------------------------
# Build OptimalControl problem
# -----------------------------
ocp_model = eval(problem)(OptimalControlBackend())
nlp_oc = nlp_model(ocp_model)
docp = eval(problem)(OptimalControlBackend())
nlp_oc = nlp_model(docp)
ocp_oc = ocp_model(docp)

# 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)
ocp_sol = build_ocp_solution(docp, nlp_oc_sol)

# get dimensions
n = state_dimension(ocp_oc)
m = control_dimension(ocp_oc)

# -----------------------------
# Plot OptimalControl solution
Expand All @@ -150,13 +130,13 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
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)),
size=(816, 220*(n+m)),
label="OptimalControl",
leftmargin=$LEFT_MARGIN,
)

# Hide legend for additional state plots
for i in 2:n_states
for i in 2:n
plot!(plt[i]; legend=:none)
end

Expand All @@ -171,28 +151,28 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
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)
t_grid = time_grid(nlp_jp)
x_fun = state(nlp_jp)
u_fun = control(nlp_jp)
p_fun = costate(nlp_jp)

# -----------------------------
# Plot JuMP solution on top
# -----------------------------
# States
for i in 1:n_states
for i in 1:n
label = i == 1 ? "JuMP" : :none
plot!(plt[i], t_grid, t -> x_fun(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)
for i in 1:n
plot!(plt[n+i], t_grid, t -> -p_fun(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)
for i in 1:m
plot!(plt[2*n+i], t_grid, t -> u_fun(t)[i]; color=2, linestyle=:dash, label=:none)
end

return plt
Expand Down Expand Up @@ -319,16 +299,16 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
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)
t_jp = time_grid(nlp_jp)
x_jp = state(nlp_jp).(t_jp)
u_jp = control(nlp_jp).(t_jp)
o_jp = objective(nlp_jp)
v_jp = variable(nlp_jp)
i_jp = iterations(nlp_jp)

x_vars = metadata(problem)[:state_components]
u_vars = metadata(problem)[:control_components]
v_vars = metadata(problem)[:variable_components]
x_vars = state_components(nlp_jp)
u_vars = control_components(nlp_jp)
v_vars = variable_components(nlp_jp)

println("┌─ ", string(problem))
println("│")
Expand Down Expand Up @@ -397,8 +377,8 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
ocp_sol = build_ocp_solution(docp, nlp_oc_sol)

# dimensions
n = state_dimension(ocp_sol) # or length(metadata(:$PROBLEM)[:state_components])
m = control_dimension(ocp_sol) # or length(metadata(:$PROBLEM)[:control_components])
n = state_dimension(ocp_sol)
m = control_dimension(ocp_sol)

# from OptimalControl solution
plt = plot(
Expand All @@ -413,10 +393,10 @@ function generate_documentation(PROBLEM::String, DESCRIPTION::String; draft::Uni
end

# from JuMP solution
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
t = time_grid(nlp_jp) # t0, ..., tN = tf
x = state(nlp_jp) # function of time
u = control(nlp_jp) # function of time
p = costate(nlp_jp) # function of time

for i in 1:n # state
label = i == 1 ? "JuMP" : :none
Expand Down
47 changes: 21 additions & 26 deletions docs/src/dev-add.md
Original file line number Diff line number Diff line change
Expand Up @@ -7,15 +7,6 @@ To add a new problem to **OptimalControlProblems**, you must follow these steps:
```julia
new_problem_meta = OrderedDict(
:grid_size => 100, # Number of steps
:state_components => ["x₁", "x₂"], # Names of the state components
:costate_components => ["∂x₁", "∂x₂"], # Names of the dynamics constraints (for the costate)
:control_components => ["u"], # Names of the control components
:variable_components => ["v"], # Names of the optimisation variables
:time_grid_names => Dict(
:initial_time => "t0", # Name of the initial time
:final_time => "tf", # Name of the final time
:grid_size => "N", # Name of the grid size
),
:parameters => (
t0 = 0, # Value of the initial time
tf = 1, # Value of the final time
Expand Down Expand Up @@ -109,20 +100,27 @@ function OptimalControlProblems.new_problem(
# model
model = JuMP.Model(args...; kwargs...)

# ------------------------------------------------
# expressions to get grid time infos
@expressions(
# metadata: required
model[:time_grid] = () -> range(t0, tf, grid_size+1) # tf is a fixed
model[:state_components] = ["x₁", "x₂"]
model[:costate_components] = ["∂x₁", "∂x₂"]
model[:control_components] = ["u"]
model[:variable_components] = String[] # no variable

# N = grid_size
@expression(model, N, grid_size)

# variables and initial guess
@variables(
model,
begin
t0, t0 # (required if the initial time is fixed)
tf, tf # (required if the final time is fixed)
N, grid_size # (required)
x₁[0:N], (start = 0.5) # consistent with model[:state_components]
x₂[0:N], (start = 0.1) # consistent with model[:state_components]
u[0:N], (start = 0.1) # consistent with model[:control_components]
end
)
# ------------------------------------------------

# define the problem
# @variables, @constraints, @objective...

# @constraints, @objective...

return model

Expand All @@ -131,14 +129,11 @@ end

!!! warning

The names `t0`, `tf` and `N` in `@expressions` command are the same as in the metadata:

- The metadata in JuMP are required and must be consistent with the other models.
- If `tf` is free, then define:

```julia
:time_grid_names => Dict(
:initial_time => "t0",
:final_time => "tf",
:grid_size => "N",
),
model[:time_grid] = () -> range(t0, value(model[:tf]), grid_size+1) # tf is a free
```

!!! tip
Expand Down
20 changes: 10 additions & 10 deletions docs/src/tutorial-solve.md
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ From `ocp_sol` you can also plot the state, control, and costate trajectories. F

```@example main
using Plots
plt = plot(ocp_sol; color=1, size=(800, 700), control_style=(label="OptimalControl", ))
plt = plot(ocp_sol; color=1, size=(700, 600), control_style=(label="OptimalControl", ))
```

## Solving from a JuMP model
Expand Down Expand Up @@ -122,12 +122,12 @@ objective_value(nlp)
To get the time grid, state, control, and costate, but also to retrieve the number of iterations and the objective value, OptimalControlProblems provides the following getters:

```@example main
t = time_grid(:beam, nlp) # t0, ..., tN = tf
x = state(:beam, nlp) # function of time
u = control(:beam, nlp) # function of time
p = costate(:beam, nlp) # function of time
o = objective(:beam, nlp) # scalar objective value
i = iterations(:beam, nlp) # number of iteration
t = time_grid(nlp) # t0, ..., tN = tf
x = state(nlp) # function of time
u = control(nlp) # function of time
p = costate(nlp) # function of time
o = objective(nlp) # scalar objective value
i = iterations(nlp) # number of iteration

tf = t[end]
println("tf = ", tf)
Expand All @@ -142,14 +142,14 @@ println("iterations: ", i)
If the `problem` includes additional optimisation variables, such as the final time, you can retrieve them with:

```julia
v = variable(problem, nlp)
v = variable(nlp)
```

Now, we can add the state, costate, and control to the plot:

```@example main
n = length(metadata(:beam)[:state_components]) # dimension of the state
m = length(metadata(:beam)[:control_components]) # dimension of the control
n = state_dimension(nlp)
m = control_dimension(nlp)

for i in 1:n # state
plot!(plt[i], t, t -> x(t)[i]; color=2, linestyle=:dash, label=:none)
Expand Down
Loading
Loading