I am currently working with some complicated problems where I want to make a small change to my objective function (just increasing a constant) if some statement is not satisfied in an iteration. I want to check this while IPOPT is running (before optimize! (model)is finished) Is this possible? Can I use the callback function to change the objective function from one iteration to the next?
IPOPT has a callback implemented, but it gives you access only to a handful information, see for instance ( https://github.com/jump-dev/Ipopt.jl ):
using JuMP, Ipopt, Test
model = Model(Ipopt.Optimizer)
set_silent(model)
@variable(model, x >= 1)
@objective(model, Min, x + 0.5)
x_vals = Float64[]
function my_callback(
alg_mod::Cint,
iter_count::Cint,
obj_value::Float64,
inf_pr::Float64,
inf_du::Float64,
mu::Float64,
d_norm::Float64,
regularization_size::Float64,
alpha_du::Float64,
alpha_pr::Float64,
ls_trials::Cint,
)
push!(x_vals, callback_value(model, x))
@test isapprox(obj_value, 1.0 * x_vals[end] + 0.5, atol = 1e-1)
# return `true` to keep going, or `false` to terminate the optimization.
return iter_count < 1
end
MOI.set(model, Ipopt.CallbackFunction(), my_callback)
optimize!(model)
@test MOI.get(model, MOI.TerminationStatus()) == MOI.INTERRUPTED
@test length(x_vals) == 2
However, you cannot modify the problem. One of the main reasons is that IPOPT is not coded in Julia, so a conversion in C is done along the way.
The only solution I see is to use a flexible modeling structure (one that would let you modify the objective) and a solver coded in pure Julia. Still, it is a bit of work.
For instance, if you deal with unconstrained optimization problems you could use lbfgs as follows:
using ADNLPModels, JSOSolvers, ForwardDiff, NLPModels
# Using ADNLPModel for its flexibility.
# But we need to deactivate the pre-computation done on the objective function
# https://juliasmoothoptimizers.github.io/tutorials/generic-adnlpmodels/
struct GenericGradientBackend <: ADNLPModels.ADBackend end
GenericGradientBackend(args...; kwargs...) = GenericGradientBackend()
function ADNLPModels.gradient!(::GenericGradientBackend, g, f, x)
return ForwardDiff.gradient!(g, f, x)
end
nlp = ADNLPModel(
x -> 100 * (x[2] - 2)^2 + 2 * (x[1] - 1)^2 + (x[1] - x[2]^2)^2,
[-1.2; 1.0],
gradient_backend = GenericGradientBackend,
)
# Then, we define the callback
function cb(nlp, solver, stats)
@show stats.iter
@show solver.x
if stats.iter > 1
nlp.f = x -> 2 * (x[2] - 2)^2 + 100 * (x[1] - 1)^2 + (x[1] - x[2]^2)^2
end
@show obj(nlp, solver.x)
end
# call JSOSolvers.lbfgs https://juliasmoothoptimizers.github.io/tutorials/introduction-to-jsosolvers/
stats = lbfgs(nlp, callback = cb, verbose = 1)
the output will be
julia> stats = lbfgs(nlp, callback = cb, verbose = 1)
[ Info: iter f(x) ‖∇f‖ ∇fᵀd bk
stats.iter = 0
solver.x = [-1.2, 1.0]
obj(nlp, solver.x) = 114.52000000000001
[ Info: 0 1.1e+02 1.9e+02 -3.7e+04 6
stats.iter = 1
solver.x = [-1.1459328, 1.7831552000000002]
obj(nlp, solver.x) = 32.62282328590102
[ Info: 1 3.3e+01 2.1e+01 -2.2e+00 0
stats.iter = 2
solver.x = [-1.0645542968025667, 1.8450655818506445]
obj(nlp, solver.x) = 446.25681758709055
[ Info: 2 3.1e+01 1.7e+01 -5.0e+00 25
stats.iter = 3
solver.x = [-1.0645542967686918, 1.8450655818604536]
obj(nlp, solver.x) = 446.25681757311787
[ Info: 3 4.5e+02 4.2e+02 -3.1e+03 1
stats.iter = 4
solver.x = [1.9246380432427963, 2.7281329116323207]
obj(nlp, solver.x) = 117.00501528447751
[ Info: 4 1.2e+02 1.8e+02 -1.9e+02 0
stats.iter = 5
solver.x = [0.9801608161355717, 2.3095374587434887]
obj(nlp, solver.x) = 19.186582034724903
[ Info: 5 1.9e+01 4.3e+01 -4.8e+00 0
stats.iter = 6
solver.x = [0.9652544626973018, 2.189462921736651]
obj(nlp, solver.x) = 14.849879523336808
[ Info: 6 1.5e+01 3.7e+01 -1.9e+01 0
stats.iter = 7
solver.x = [1.0002184791527038, 1.6363671414275927]
obj(nlp, solver.x) = 3.078398087147021
[ Info: 7 3.1e+00 1.0e+01 -2.1e+00 0
stats.iter = 8
solver.x = [1.0089309580533612, 1.4229017930719876]
obj(nlp, solver.x) = 1.7057450643634193
[ Info: 8 1.7e+00 3.5e+00 -4.3e-01 0
stats.iter = 9
solver.x = [1.0053404604019518, 1.2982558231123305]
obj(nlp, solver.x) = 1.4503155493415034
[ Info: 9 1.5e+00 7.8e-01 -2.4e-02 0
stats.iter = 10
solver.x = [1.0064110790068321, 1.2657336397712304]
obj(nlp, solver.x) = 1.4372277943232419
[ Info: 10 1.4e+00 1.2e-01 -8.8e-04 1
stats.iter = 11
solver.x = [1.0047305130737216, 1.2632117110472276]
obj(nlp, solver.x) = 1.4372011985565385
[ Info: 11 1.4e+00 2.4e-01 -3.3e-04 0
stats.iter = 12
solver.x = [1.00584984478302, 1.2615152914236158]
obj(nlp, solver.x) = 1.4370347773092513
[ Info: 12 1.4e+00 1.5e-03 -4.5e-08 0
stats.iter = 13
solver.x = [1.0058545642521763, 1.2614703004604053]
obj(nlp, solver.x) = 1.4370347544382267
[ Info: 13 1.4e+00 9.6e-06 -2.5e-12 0
stats.iter = 14
solver.x = [1.005854516585087, 1.26146995704064]
obj(nlp, solver.x) = 1.4370347544370015
[ Info: 14 1.4e+00 5.6e-07
"Execution stats: first-order stationary"