Skip to content

Commit

Permalink
prepare CTDirect 0.7 for CTBase 0.10 (#116)
Browse files Browse the repository at this point in the history
* ok for free times, objective type and variable indicators

* constraints dims seem ok

* added test prof for ipopt stats

* test ok, seems ready to merge

* updated calls to former OCPInit

* tabs

* update CTBase to 0.10

* fixed docs
  • Loading branch information
PierreMartinon authored Jun 21, 2024
1 parent eb71e0e commit 6e3339e
Show file tree
Hide file tree
Showing 20 changed files with 453 additions and 739 deletions.
6 changes: 3 additions & 3 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "CTDirect"
uuid = "790bbbee-bee9-49ee-8912-a9de031322d5"
authors = ["Olivier Cots <[email protected]>"]
version = "0.6.1"
version = "0.7.0"

[deps]
ADNLPModels = "54578032-b7ea-4c30-94aa-7cbd1cce6c9a"
Expand All @@ -16,8 +16,8 @@ NLPModelsIpopt = "f4238b75-b362-5c4c-b852-0801c9a21d71"
StructTypes = "856f2bd8-1eba-4b0a-8007-ebc267875bd4"

[compat]
HSL = "0.4"
ADNLPModels = "0.8"
CTBase = "0.9"
CTBase = "0.10"
DocStringExtensions = "0.9"
HSL = "0.4"
julia = "1.10"
1 change: 1 addition & 0 deletions docs/make.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
using Documenter
using CTDirect
using CTBase
using Plots

DocMeta.setdocmeta!(CTBase, :DocTestSetup, :(using CTBase); recursive = true)
DocMeta.setdocmeta!(CTDirect, :DocTestSetup, :(using CTDirect); recursive = true)
Expand Down
70 changes: 13 additions & 57 deletions docs/src/continuation.md
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,6 @@ The most compact syntax to perform a discrete continuation is to use a function

```@setup main
using Printf
using Statistics
using Plots
```

Expand Down Expand Up @@ -42,7 +41,7 @@ nothing # hide
Then we perform the continuation with a simple *for* loop, using each solution to initialize the next problem.

```@example main
init1 = OCPInit()
init1 = OptimalControlInit()
iter_list = []
for T=1:5
ocp1 = ocp_T(T)
Expand Down Expand Up @@ -76,16 +75,21 @@ function F1(x)
end
ocp = Model(variable=true)
r0 = 1
v0 = 0
m0 = 1
mf = 0.6
x0=[r0,v0,m0]
vmax = 0.1
state!(ocp, 3)
control!(ocp, 1)
variable!(ocp, 1)
time!(ocp, 0, Index(1))
constraint!(ocp, :initial, [1,0,1], :initial_constraint)
constraint!(ocp, :final, Index(3), 0.6, :final_constraint)
constraint!(ocp, :state, 1:2:3, [1,0.6], [1.2,1], :state_box)
constraint!(ocp, :control, Index(1), 0, 1, :control_box)
constraint!(ocp, :variable, Index(1), 0.01, Inf, :variable_box)
constraint!(ocp, :state, Index(2), 0, Inf, :speed_limit)
time!(ocp, t0=0, indf=1)
constraint!(ocp, :initial, lb=x0, ub=x0)
constraint!(ocp, :final, rg=3, lb=mf, ub=Inf)
constraint!(ocp, :state, lb=[r0,v0,mf], ub=[r0+0.2,vmax,m0])
constraint!(ocp, :control, lb=0, ub=1)
constraint!(ocp, :variable, lb=0.01, ub=Inf)
objective!(ocp, :mayer, (x0, xf, v) -> xf[1], :max)
dynamics!(ocp, (x, u, v) -> F0(x) + u*F1(x) )
Expand Down Expand Up @@ -118,51 +122,3 @@ plot(sol0)
p = plot!(sol)
plot(pobj, p, layout=2)
```


## Manual constraint redefinition

Here we illustrate a slightly more involved way of modifying the OCP problem during the continuation.
Instead of just updating a global variable as before, we now remove and redefine one of the constraints (maximal speed).
```@example main
global Tmax = 3.5
vmax_list = []
obj_list = []
iter_list = []
print("vmax ")
for vmax=0.15:-0.01:0.05
print(vmax," ")
remove_constraint!(ocp, :speed_limit)
constraint!(ocp, :state, Index(2), 0, vmax, :speed_limit)
global sol = solve(ocp, print_level=0, init=sol)
push!(vmax_list, vmax)
push!(obj_list, sol.objective)
push!(iter_list, sol.iterations)
end
@printf("\nAverage iterations %d\n", mean(iter_list))
```

We now plot the objective with respect to the speed limit, as well as a comparison of the solutions for the unconstrained case and the *vmax*=0.05 case.

```@example main
pobj = plot(vmax_list, obj_list, label="r(tf)",seriestype=:scatter)
xlabel!("Speed limit (vmax)")
ylabel!("Maximal altitude r(tf)")
plot(sol0)
p = plot!(sol)
plot(pobj, p, layout=2)
```

We can compare with solving each problem with the default initial guess, which here gives the same solutions but takes more iterations overall.

```@example main
iter_list = []
for vmax=0.15:-0.01:0.05
print(vmax," ")
remove_constraint!(ocp, :speed_limit)
constraint!(ocp, :state, Index(2), 0, vmax, :speed_limit)
global sol = solve(ocp, print_level=0)
push!(iter_list, sol.iterations)
end
@printf("\nAverage iterations %d\n", mean(iter_list))
```
9 changes: 5 additions & 4 deletions docs/src/tutorial.md
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@ First import the CTDirect and CTBase modules
```@example main
using CTDirect
using CTBase
using Plots
```

Then define the OCP for the Goddard problem. Note that the free final time is modeled as an *optimization variable*, and has both a lower bound to prevent the optimization getting stuck at tf=0. In this particular case an upper bound is not needed for the final time since the final mass is prescribed.
Expand Down Expand Up @@ -78,22 +79,22 @@ Let us start with the simplest case, constant initialisation.
x_const = [1.05, 0.2, 0.8]
u_const = 0.5
v_const = 0.15
init1 = OCPInit(state=x_const, control=u_const, variable=v_const)
init1 = (state=x_const, control=u_const, variable=v_const)
sol2 = solve(ocp, print_level=0, init=init1)
println("Objective ", sol2.objective, " after ", sol2.iterations, " iterations")
```

Now we illustrate the functional initialisation, with some random functions. Note that we only consider the state and control variables, since the optimization variables are scalar and therefore a functional initialisation is not relevant. In the example notice that the call to **OCPInit** does not provide an argument for the optimization variables, therefore the default initial guess will be used.
Now we illustrate the functional initialisation, with some random functions. Note that we only consider the state and control variables, since the optimization variables are scalar and therefore a functional initialisation is not relevant. In the example notice that the initialization does not provide an argument for the optimization variables, therefore the default initial guess will be used.
```@example main
x_func = t->[1+t^2, sqrt(t), 1-t]
u_func = t->(cos(t)+1)*0.5
init2 = OCPInit(state=x_func, control=u_func)
init2 = (state=x_func, control=u_func)
sol3 = solve(ocp, print_level=0, init=init2)
println("Objective ", sol3.objective, " after ", sol3.iterations, " iterations")
```
More generally, the default, constant and functional initialisations can be mixed, as shown in the example below that uses a functional initial guess for the state, a constant initial guess for the control, and the default initial guess for the optimization variables.
```@example main
init3 = OCPInit(state=x_func, control=u_const)
init3 = (state=x_func, control=u_const)
sol4 = solve(ocp, print_level=0, init=init3)
println("Objective ", sol4.objective, " after ", sol4.iterations, " iterations")
```
Expand Down
Loading

0 comments on commit 6e3339e

Please sign in to comment.