Hello,
I am trying to solve an NLP modeled with InfiniteOpt/JuMP using the SLSQP solver in NLopt.jl
However, the solver fails and I am struggling to investigate the root cause of this failure.
Is there a way to have a more verbose output from NLopt, similar to that of Ipopt with print_level set to 5 or greater? It can be either printed to the REPL/stdout, written on a log file, or retrieved using a custom callback.
Thanks!
There’s no easy solution for this. NLopt does not have a verbose output setting. The usual answer is to add a print to your oracle functions, but you don’t have access to these via JuMP.
Just use Ipopt instead?
Can you post your JuMP model?
Ok, thanks.
Yes, I can use Ipopt, but I was interested in testing other solvers since Ipopt also fails in some of my problems.
This is an example which is easily solved by Ipopt, but on which NLopt with SLSQP fails.
It is the optimization of a launch trajectory from the Moon surface, where the rocket starts at rest and needs to reach a given altitude with a given tangential velocity.
using InfiniteOpt, Ipopt, HSL_jll, NLopt
using Interpolations
using Printf
# optimizer
opt = NLopt.Optimizer; # change to Ipopt to make it work
# physical constants
mu = 4.9028e12; # Moon standard gravitational parameter [m^3/s^2]
rm = 1737.4e3; # Moon radius [m]
g0 = 9.80665; # sea-level gravity [m/s^2]
Isp = 450.0; # specific impulse [s]
twr = 2.1; # thrust/weight ratio [-]
# scaling parameters
lc = rm; # characteristic length [m]
tc = sqrt(rm^3 / mu); # characteristic time [s]
vc = lc / tc; # characteristic velocity [m/s]
# boundary conditions
r0 = 1.0; # initial radius [-]
θ0 = 0.0; # initial angle [-]
u0 = 0.0; # initial radial velocity [-]
v0 = 0.0; # initial tangential velocity [-]
m0 = 1.0; # initial mass [-]
rf = 1.05; # target radius [-]
uf = 0.0; # target radial velocity [-]
vf = sqrt(1.0 / rf); # target tangential velocity [-]
θfg = π/18.0; # initial guess for final angle [-]
mfg = 0.7; # initial guess for final mass [-]
αsg = [π/3.0, -π/6.0]; # initial guess for thrust direction [-]
# time of flight
tof0 = 0.5; # initial guess for time of flight [s]
N = 30; # number of discrete points [-]
# constants
ṁ = -twr / (Isp * g0 / vc); # mass flow rate [-]
# model
mdl = InfiniteModel(opt);
# optimizer options
if mdl.optimizer_constructor == Ipopt.Optimizer
set_optimizer_attribute(mdl, "hsllib", HSL_jll.libhsl_path);
set_optimizer_attribute(mdl, "tol", 1e-12);
set_optimizer_attribute(mdl, "linear_solver", "ma57");
set_optimizer_attribute(mdl, "print_level", 0);
elseif mdl.optimizer_constructor == NLopt.Optimizer
set_optimizer_attribute(mdl, "algorithm", :LD_SLSQP);
end
# time of flight discretization
s = [0.0, 1.0]; # pseudo-time interval
@infinite_parameter(mdl, t in s, num_supports=N);
# state and control variables
@variable(mdl, 0.25 <= tof <= 0.75, start=tof0); # time of flight
@variables(mdl, begin
r >= 1.0, Infinite(t) # radial position
θ >= 0.0, Infinite(t) # angular position
u, Infinite(t) # radial velocity
v, Infinite(t) # tangential velocity
m >= 0.0, Infinite(t) # mass
-π/2.0 <= α <= π/2.0, Infinite(t) # thrust direction
end);
# initial guess
set_start_value_function(r, t -> linear_interpolation(s, [r0, rf])(t));
set_start_value_function(θ, t -> linear_interpolation(s, [θ0, θfg])(t));
set_start_value_function(u, t -> linear_interpolation(s, [u0, uf])(t));
set_start_value_function(v, t -> linear_interpolation(s, [v0, vf])(t));
set_start_value_function(m, t -> linear_interpolation(s, [m0, mfg])(t));
set_start_value_function(α, t -> linear_interpolation(s, αsg)(t));
# dynamics
@constraint(mdl, ∂(r, t) == tof * u);
@constraint(mdl, ∂(θ, t) == tof * (v / r));
@constraint(mdl, ∂(u, t) == tof * (-1.0 / r^2 + v^2 / r + twr / m * sin(α)));
@constraint(mdl, ∂(v, t) == tof * (-(u * v) / r + twr / m * cos(α)));
@constraint(mdl, ∂(m, t) == tof * ṁ);
# boundary constraints
@constraint(mdl, r(s[1]) == r0);
@constraint(mdl, θ(s[1]) == θ0);
@constraint(mdl, u(s[1]) == u0);
@constraint(mdl, v(s[1]) == v0);
@constraint(mdl, m(s[1]) == m0);
@constraint(mdl, r(s[2]) == rf);
@constraint(mdl, u(s[2]) == uf);
@constraint(mdl, v(s[2]) == vf);
# objective
@objective(mdl, Min, tof);
build_optimizer_model!(mdl);
# solution
InfiniteOpt.optimize!(mdl);
display(solution_summary(mdl.optimizer_model));
tv = value(t) * value(tof);
rv = value(r);
θv = value(θ);
uv = value(u);
vv = value(v);
mv = value(m);
αv = value(α);
αv[1] = αv[2];
# output
@printf "Time of flight: %.3f s\n" tc * tv[end];
@printf "Final mass: %.3f kg\n" mv[end];
Output with Ipopt
* Solver : Ipopt
* Status
Result count : 1
Termination status : LOCALLY_SOLVED
Message from the solver:
"Solve_Succeeded"
* Candidate solution (result #1)
Primal status : FEASIBLE_POINT
Dual status : FEASIBLE_POINT
Objective value : 4.57374e-01
Dual objective value : 6.46409e-01
* Work counters
Solve time (sec) : 1.23205e+00
Barrier iterations : 11
Time of flight: 473.041 s
Final mass: 0.634 kg
Output with NLopt (it returns the intial guess)
* Solver : NLopt
* Status
Result count : 1
Termination status : OTHER_ERROR
Message from the solver:
"FAILURE"
* Candidate solution (result #1)
Primal status : UNKNOWN_RESULT_STATUS
Dual status : NO_SOLUTION
Objective value : 5.00000e-01
* Work counters
Solve time (sec) : 5.93421e-01
Time of flight: 517.128 s
Final mass: 0.700 kg