Skip to content

Commit e8fc84b

Browse files
PierreMartinonPierreMartinon
committed
todo: fix shuttle error (test_comparison:89, * with nh vector)
1 parent d886522 commit e8fc84b

File tree

9 files changed

+123
-110
lines changed

9 files changed

+123
-110
lines changed

‎ext/JuMPModels/moonlander.jl‎

Lines changed: 9 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -23,7 +23,7 @@ function OptimalControlProblems.moonlander(
2323
@variables(
2424
model,
2525
begin
26-
0.0 <= step, (start = 0.1)
26+
0.01 <= tf, (start = 0.1)
2727
# state variables
2828
p1[k=0:nh], (start = 0.1)
2929
p2[k=0:nh], (start = 0.1)
@@ -55,6 +55,13 @@ function OptimalControlProblems.moonlander(
5555
)
5656

5757
#dynamics
58+
@expressions(
59+
model,
60+
begin
61+
step, tf / nh
62+
end
63+
)
64+
5865
@expressions(
5966
model,
6067
begin
@@ -95,7 +102,7 @@ function OptimalControlProblems.moonlander(
95102
end
96103
)
97104

98-
@objective(model, Min, step * nh)
105+
@objective(model, Min, tf)
99106

100107
return model
101108
end

‎ext/MetaData/cart_pendulum.jl‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,6 @@ cart_pendulum_meta = Dict(
66
:minimize => true,
77
:state_name => ["x", "dx", "theta", "omega"],
88
:costate_name => ["d_x", "d_dx", "d_theta", "d_omega"],
9-
:control_name => "Fex",
9+
:control_name => ["Fex"],
1010
:time => ("final_time", "tf", nothing)
11-
)
11+
)

‎ext/MetaData/glider.jl‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,6 @@ glider_meta = Dict(
66
:minimize => false,
77
:state_name => ["x", "y", "vx", "vy"],
88
:costate_name => ["x_eqn", "y_eqn", "vx_eqn", "vy_eqn"],
9-
:control_name => "cL",
9+
:control_name => ["cL"],
1010
:time => ("final_time", "tf", nothing)
11-
)
11+
)

‎test/figures/cart_pendulum.png‎

147 KB
Loading

‎test/figures/electrical_vehicle.png‎

-85.6 KB
Binary file not shown.

‎test/figures/glider.png‎

125 KB
Loading

‎test/figures/moonlander.png‎

182 KB
Loading

‎test/runtests.jl‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -47,8 +47,8 @@ const verbose = true # print or not details during tests
4747
#:aqua,
4848
#:JuMP, # convergence tests for JuMP models
4949
#:OptimalControl, # convergence tests for OptimalControl models
50-
#:Comparison, # comparison between OptimalControl and JuMP
51-
:quick, # quick comparison: objective rel error only
50+
:Comparison, # comparison between OptimalControl and JuMP
51+
#:quick, # quick comparison: objective rel error only
5252
)
5353
@testset "$(name)" verbose=verbose begin
5454
test_name = Symbol(:test_, name)

‎test/test_Comparison.jl‎

Lines changed: 108 additions & 102 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,7 @@
11
function test_Comparison()
22

3+
test_init = false
4+
35
# Comparison Parameters
46
ε = 1e-2
57
p = 2
@@ -35,125 +37,127 @@ function test_Comparison()
3537
println("############ TEST $f #############")
3638
println()
3739

38-
#================== INIT ======================#
40+
if test_init
41+
#================== INIT ======================#
3942

40-
########## OptimalControl ##########
43+
########## OptimalControl ##########
4144

42-
# Set up the OptimalControl model
43-
docp_init, OC_model_init = OptimalControlProblems.eval(f)(OptimalControlBackend()) # +++ UPDATE
45+
# Set up the OptimalControl model
46+
docp_init, OC_model_init = OptimalControlProblems.eval(f)(OptimalControlBackend()) # +++ UPDATE
4447

45-
# Solve the problem
46-
nlp_sol_init = NLPModelsIpopt.ipopt(OC_model_init; kwargs_init...)
48+
# Solve the problem
49+
nlp_sol_init = NLPModelsIpopt.ipopt(OC_model_init; kwargs_init...)
4750

48-
# Build the solution
49-
sol_init = build_OCP_solution(docp_init; primal=nlp_sol_init.solution, dual=nlp_sol_init.multipliers)
51+
# Build the solution
52+
sol_init = build_OCP_solution(docp_init; primal=nlp_sol_init.solution, dual=nlp_sol_init.multipliers)
5053

51-
# Retrieves values of variables
52-
x_init_oc = state(sol_init)
53-
p_init_oc = costate(sol_init)
54-
u_init_oc = control(sol_init)
54+
# Retrieves values of variables
55+
x_init_oc = state(sol_init)
56+
p_init_oc = costate(sol_init)
57+
u_init_oc = control(sol_init)
5558

56-
############### JuMP ###############
59+
############### JuMP ###############
5760

58-
# Set up the JuMP model
59-
JuMP_init_model = OptimalControlProblems.eval(f)(JuMPBackend())
60-
set_optimizer(JuMP_init_model, Ipopt.Optimizer)
61-
set_silent(JuMP_init_model)
62-
set_optimizer_attribute(JuMP_init_model, "tol", tol)
63-
set_optimizer_attribute(JuMP_init_model, "max_iter", 0)
64-
set_optimizer_attribute(JuMP_init_model, "mu_strategy", mu_strategy)
65-
set_optimizer_attribute(JuMP_init_model, "linear_solver", "mumps")
66-
set_optimizer_attribute(JuMP_init_model, "max_wall_time", max_wall_time)
67-
set_optimizer_attribute(JuMP_init_model, "sb", sb)
61+
# Set up the JuMP model
62+
JuMP_init_model = OptimalControlProblems.eval(f)(JuMPBackend())
63+
set_optimizer(JuMP_init_model, Ipopt.Optimizer)
64+
set_silent(JuMP_init_model)
65+
set_optimizer_attribute(JuMP_init_model, "tol", tol)
66+
set_optimizer_attribute(JuMP_init_model, "max_iter", 0)
67+
set_optimizer_attribute(JuMP_init_model, "mu_strategy", mu_strategy)
68+
set_optimizer_attribute(JuMP_init_model, "linear_solver", "mumps")
69+
set_optimizer_attribute(JuMP_init_model, "max_wall_time", max_wall_time)
70+
set_optimizer_attribute(JuMP_init_model, "sb", sb)
6871

69-
# Solve the model
70-
optimize!(JuMP_init_model)
72+
# Solve the model
73+
optimize!(JuMP_init_model)
7174

72-
# Retrieves values of variables
73-
time_data, time_var_name, time_value = OptimalControlProblems.metadata[f][:time]
74-
if time_data == "final_time"
75-
if time_value !== nothing
76-
tf_init = time_value
77-
else
78-
tf_init = value.(JuMP_init_model[Symbol(time_var_name)])
79-
end
80-
h_init = tf_init / nh
81-
elseif time_data == "step"
82-
if time_value !== nothing
83-
h_init = time_value
84-
else
85-
h_init = value.(JuMP_init_model[Symbol(time_var_name)])
75+
# Retrieves values of variables
76+
time_data, time_var_name, time_value = OptimalControlProblems.metadata[f][:time]
77+
if time_data == "final_time"
78+
if time_value !== nothing
79+
tf_init = time_value
80+
else
81+
tf_init = value.(JuMP_init_model[Symbol(time_var_name)])
82+
end
83+
h_init = tf_init / nh
84+
elseif time_data == "step"
85+
if time_value !== nothing
86+
h_init = time_value
87+
else
88+
h_init = value.(JuMP_init_model[Symbol(time_var_name)])
89+
end
90+
tf_init = h_init * nh
8691
end
87-
tf_init = h_init * nh
88-
end
89-
t_init = Vector((0:nh) * h_init)
90-
91-
x_vars = OptimalControlProblems.metadata[f][:state_name]
92-
p_vars = OptimalControlProblems.metadata[f][:costate_name]
93-
u_vars = OptimalControlProblems.metadata[f][:control_name]
94-
95-
x_jmp_vars_init = [JuMP.value.(JuMP_init_model[Symbol(xv)]) for xv in x_vars]
96-
inds_x_init = axes(x_jmp_vars_init[1], 1)
97-
x_jmp_init = [[x_jmp_vars_init[j][i] for j in 1:length(x_vars)] for i in inds_x_init]
98-
99-
u_jmp_vars_init = [JuMP.value.(JuMP_init_model[Symbol(uv)]) for uv in u_vars]
100-
inds_u_init = axes(u_jmp_vars_init[1], 1)
101-
u_jmp_init = [[u_jmp_vars_init[j][i] for j in 1:length(u_vars)] for i in inds_u_init]
102-
103-
p_jmp_vars_init = [JuMP.dual.(JuMP_init_model[Symbol(pv)]) for pv in p_vars]
104-
inds_p_init = axes(p_jmp_vars_init[1], 1)
105-
p_jmp_init = -[[p_jmp_vars_init[j][i] for j in 1:length(p_vars)] for i in inds_p_init]
106-
p_jmp_init = costateInterpolation(p_jmp_init, t_init)
107-
108-
############ TEST ############
109-
@testset "init" verbose=verbose begin
110-
print("Init:\n")
111-
for k in 1:length(x_jmp_init[1])
112-
dist_x_init = abs(x_init_oc(0)[k] - x_jmp_init[1][k])
113-
print(" Test x$k : ")
114-
@testset "x$k" verbose=verbose begin
115-
if !(dist_x_init < ε)
116-
print("$dist_x_init < $ε \033[1;33mTest Broken\033[0m\n")
117-
@test dist_x_init < ε broken=true
118-
global list_of_problems_final
119-
list_of_problems_final = setdiff(list_of_problems_final, [f])
120-
else
121-
print("$dist_x_init < $ε \033[1;32mTest Passed\033[0m\n")
122-
@test dist_x_init < ε
92+
t_init = Vector((0:nh) * h_init)
93+
94+
x_vars = OptimalControlProblems.metadata[f][:state_name]
95+
p_vars = OptimalControlProblems.metadata[f][:costate_name]
96+
u_vars = OptimalControlProblems.metadata[f][:control_name]
97+
98+
x_jmp_vars_init = [JuMP.value.(JuMP_init_model[Symbol(xv)]) for xv in x_vars]
99+
inds_x_init = axes(x_jmp_vars_init[1], 1)
100+
x_jmp_init = [[x_jmp_vars_init[j][i] for j in 1:length(x_vars)] for i in inds_x_init]
101+
102+
u_jmp_vars_init = [JuMP.value.(JuMP_init_model[Symbol(uv)]) for uv in u_vars]
103+
inds_u_init = axes(u_jmp_vars_init[1], 1)
104+
u_jmp_init = [[u_jmp_vars_init[j][i] for j in 1:length(u_vars)] for i in inds_u_init]
105+
106+
p_jmp_vars_init = [JuMP.dual.(JuMP_init_model[Symbol(pv)]) for pv in p_vars]
107+
inds_p_init = axes(p_jmp_vars_init[1], 1)
108+
p_jmp_init = -[[p_jmp_vars_init[j][i] for j in 1:length(p_vars)] for i in inds_p_init]
109+
p_jmp_init = costateInterpolation(p_jmp_init, t_init)
110+
111+
############ TEST ############
112+
@testset "init" verbose=verbose begin
113+
print("Init:\n")
114+
for k in 1:length(x_jmp_init[1])
115+
dist_x_init = abs(x_init_oc(0)[k] - x_jmp_init[1][k])
116+
print(" Test x$k : ")
117+
@testset "x$k" verbose=verbose begin
118+
if !(dist_x_init < ε)
119+
print("$dist_x_init < $ε \033[1;33mTest Broken\033[0m\n")
120+
@test dist_x_init < ε broken=true
121+
global list_of_problems_final
122+
list_of_problems_final = setdiff(list_of_problems_final, [f])
123+
else
124+
print("$dist_x_init < $ε \033[1;32mTest Passed\033[0m\n")
125+
@test dist_x_init < ε
126+
end
123127
end
124128
end
125-
end
126129

127-
for k in 1:length(p_jmp_init[1])
128-
dist_p_init = abs(p_init_oc(0)[k] - p_jmp_init[1][k])
129-
print(" Test p$k : ")
130-
@testset "p$k" verbose=verbose begin
131-
if !(dist_p_init < ε)
132-
print("$dist_p_init < $ε \033[1;33mTest Broken\033[0m\n")
133-
@test dist_p_init < ε broken=true
134-
else
135-
print("$dist_p_init < $ε \033[1;32mTest Passed\033[0m\n")
136-
@test dist_p_init < ε
130+
for k in 1:length(p_jmp_init[1])
131+
dist_p_init = abs(p_init_oc(0)[k] - p_jmp_init[1][k])
132+
print(" Test p$k : ")
133+
@testset "p$k" verbose=verbose begin
134+
if !(dist_p_init < ε)
135+
print("$dist_p_init < $ε \033[1;33mTest Broken\033[0m\n")
136+
@test dist_p_init < ε broken=true
137+
else
138+
print("$dist_p_init < $ε \033[1;32mTest Passed\033[0m\n")
139+
@test dist_p_init < ε
140+
end
137141
end
138142
end
139-
end
140143

141-
for k in 1:length(u_jmp_init[1])
142-
dist_u_init = abs(u_init_oc(0)[k] - u_jmp_init[1][k])
143-
print(" Test u$k : ")
144-
@testset "u$k" verbose=verbose begin
145-
if !(dist_u_init < ε)
146-
print("$dist_u_init < $ε \033[1;33mTest Broken\033[0m\n")
147-
@test dist_u_init < ε broken=true
148-
global list_of_problems_final
149-
list_of_problems_final = setdiff(list_of_problems_final, [f])
150-
else
151-
print("$dist_u_init < $ε \033[1;32mTest Passed\033[0m\n")
152-
@test dist_u_init < ε
144+
for k in 1:length(u_jmp_init[1])
145+
dist_u_init = abs(u_init_oc(0)[k] - u_jmp_init[1][k])
146+
print(" Test u$k : ")
147+
@testset "u$k" verbose=verbose begin
148+
if !(dist_u_init < ε)
149+
print("$dist_u_init < $ε \033[1;33mTest Broken\033[0m\n")
150+
@test dist_u_init < ε broken=true
151+
global list_of_problems_final
152+
list_of_problems_final = setdiff(list_of_problems_final, [f])
153+
else
154+
print("$dist_u_init < $ε \033[1;32mTest Passed\033[0m\n")
155+
@test dist_u_init < ε
156+
end
153157
end
154158
end
155-
end
156-
end
159+
end
160+
end
157161

158162
#======================= END INIT ========================#
159163

@@ -225,12 +229,13 @@ function test_Comparison()
225229

226230
obj_jmp = objective_value(JuMP_model)
227231

228-
dist_obj = abs(obj_oc - obj_jmp)
232+
dist_obj = abs(obj_oc - obj_jmp) / (obj_oc + obj_jmp) / 2.
229233
@testset "objective" verbose=verbose begin
230234
print("Objective:\n")
231235
print(" Test objective : ")
232236
if !(dist_obj < ε)
233237
print("$dist_obj < $ε \033[1;33mTest Broken\033[0m\n")
238+
println("Jump ", obj_jmp, " vs OC ", obj_oc)
234239
@test dist_obj < ε broken=true
235240
global list_of_problems_final
236241
list_of_problems_final = setdiff(list_of_problems_final, [f])
@@ -244,6 +249,7 @@ function test_Comparison()
244249
plots_p = Vector{Any}()
245250
plots_u = Vector{Any}()
246251

252+
# +++ use relative error
247253
@testset "norm_L$p" verbose=verbose begin
248254
print("Norm_L$p:\n")
249255
for k in 1:length(x_jmp[1])

0 commit comments

Comments
 (0)