jump-dev / Ipopt.jl

A Julia interface to the Ipopt nonlinear solver
https://github.com/coin-or/ipopt
Other
152 stars 58 forks source link

BoundsError calling MOI.set #428

Closed jalving closed 1 month ago

jalving commented 1 month ago

I am running into an error with the EpsilonConstraint algorithm for a nonlinear unconstrained optimization problem.

My problem consists of two nonlinear external functions that i register as nonlinear operators (also with external gradients).

model = build_coupled_model()
mean_displacement = model[:displacement]  # calculated from registered external function with gradient
mean_stress = model[:stress]              # calculated from a different external function with gradient
@objective(model, Min, [mean_stress/1e5,  -mean_displacement])

ipopt_optimizer = optimizer_with_attributes(
    Ipopt.Optimizer,
    "hessian_approximation" => "limited-memory",
    "print_level" => 5,
    "tol" => 1e-2,
    "max_iter" => 100
)
set_optimizer(model, () -> MOA.Optimizer(ipopt_optimizer))
set_attribute(model, MOA.Algorithm(), MOA.EpsilonConstraint())
optimize!(model)

The algorithm runs for a few iterations before the error. I included the last Ipopt solve before the error.

This is Ipopt version 3.14.16, running with linear solver MUMPS 5.7.3.

Number of nonzeros in equality constraint Jacobian...:        0
Number of nonzeros in inequality constraint Jacobian.:        4
Number of nonzeros in Lagrangian Hessian.............:        0

Total number of variables............................:        4
                     variables with only lower bounds:        0
                variables with lower and upper bounds:        4
                     variables with only upper bounds:        0
Total number of equality constraints.................:        0
Total number of inequality constraints...............:        1
        inequality constraints with only lower bounds:        0
   inequality constraints with lower and upper bounds:        0
        inequality constraints with only upper bounds:        1

iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
   0 -1.3451500e-01 0.00e+00 1.39e-01   0.0 0.00e+00    -  0.00e+00 0.00e+00   0
   1 -1.2258170e-01 0.00e+00 2.41e-01  -0.8 8.79e-02    -  9.91e-01 1.00e+00f  1
   2 -1.2542720e-01 0.00e+00 3.80e-01  -0.6 1.34e-02    -  1.00e+00 1.00e+00h  1
   3 -1.2664752e-01 0.00e+00 1.47e-01  -1.3 1.42e-02    -  1.00e+00 1.00e+00H  1
   4 -1.3681531e-01 0.00e+00 7.57e-01  -2.0 5.04e-02    -  1.00e+00 1.00e+00f  1
   5 -1.7346826e-01 0.00e+00 1.60e-01  -2.5 2.38e-01    -  9.99e-01 5.46e-01f  1
   6 -1.7426194e-01 0.00e+00 1.02e-01  -3.4 2.19e-02    -  1.00e+00 1.00e+00h  1
   7 -1.7457616e-01 0.00e+00 6.37e-01  -3.8 4.53e-02    -  9.34e-01 1.00e+00H  1
   8 -1.7566504e-01 0.00e+00 5.34e-01  -4.3 1.09e-02    -  1.00e+00 1.00e+00f  1
   9 -1.7594168e-01 0.00e+00 1.56e-01  -4.7 1.57e-02    -  1.00e+00 1.00e+00h  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  10 -1.7690110e-01 0.00e+00 4.21e-01  -3.1 1.31e-01    -  1.00e+00 2.29e-01f  1
  11 -1.7712431e-01 0.00e+00 1.35e-01  -3.7 5.27e-03    -  1.00e+00 7.52e-01h  1
  12 -1.7719035e-01 0.00e+00 9.27e-02  -5.3 1.52e-03    -  1.00e+00 1.00e+00h  1
  13 -1.7725310e-01 0.00e+00 3.02e-01  -5.7 9.26e-02    -  1.00e+00 2.15e-02h  6
  14 -1.7725439e-01 0.00e+00 1.86e-01  -6.4 1.50e-03    -  1.00e+00 1.00e+00h  1
  15 -1.7736329e-01 0.00e+00 1.33e-01  -6.3 5.30e-04    -  1.00e+00 1.00e+00h  1
  16 -1.7738020e-01 0.00e+00 1.70e-01  -7.9 2.95e-03    -  1.00e+00 2.40e-01h  3
  17 -1.7742823e-01 0.00e+00 1.00e-01  -6.6 2.93e-04    -  1.00e+00 1.00e+00h  1
  18 -1.7746374e-01 0.00e+00 1.83e-01  -7.1 1.89e-03    -  1.00e+00 8.44e-01h  1
  19 -1.7750655e-01 0.00e+00 1.56e-01  -6.3 2.65e-04    -  1.00e+00 1.00e+00h  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  20 -1.7724258e-01 0.00e+00 2.19e-01  -8.1 1.43e-03    -  1.00e+00 1.00e+00H  1
  21 -1.7726637e-01 0.00e+00 2.03e-01  -9.3 8.23e-04    -  1.00e+00 1.00e-01f  1
  22 -1.7753009e-01 0.00e+00 1.31e-01  -7.9 8.86e-04    -  1.00e+00 9.99e-01h  1
  23 -1.7746446e-01 0.00e+00 2.33e-01  -9.5 7.25e-04    -  1.00e+00 1.00e+00h  1
  24 -1.7753243e-01 0.00e+00 8.51e-02  -9.4 4.43e-04    -  1.00e+00 1.00e+00h  1

Number of Iterations....: 24

                                   (scaled)                 (unscaled)
Objective...............:  -1.7753243446350098e-01   -1.7753243446350098e-01
Dual infeasibility......:   8.5145175588206404e-02    8.5145175588206404e-02
Constraint violation....:   0.0000000000000000e+00    0.0000000000000000e+00
Variable bound violation:   9.7088203032669185e-09    9.7088203032669185e-09
Complementarity.........:   4.0164787073937268e-10    4.0164787073937268e-10
Overall NLP error.......:   8.5145175588206404e-02    8.5145175588206404e-02

Number of objective function evaluations             = 39
Number of objective gradient evaluations             = 25
Number of equality constraint evaluations            = 0
Number of inequality constraint evaluations          = 39
Number of equality constraint Jacobian evaluations   = 0
Number of inequality constraint Jacobian evaluations = 25
Number of Lagrangian Hessian evaluations             = 0
Total seconds in IPOPT                               = 59.445

EXIT: Optimal Solution Found.
ERROR: BoundsError: attempt to access 0-element Vector{Float64} at index [1]
Stacktrace:
  [1] setindex!
    @ ./array.jl:1021 [inlined]
  [2] set
    @ ~/.julia/packages/Ipopt/3yCIJ/src/utils.jl:412 [inlined]
  [3] set
    @ ~/.julia/packages/Ipopt/3yCIJ/src/MOI_wrapper.jl:415 [inlined]
  [4] _replace_constraint_function_or_set(m::MathOptInterface.Utilities.CachingOptimizer{…}, attr::MathOptInterface.ConstraintSet, cindex::MathOptInterface.ConstraintIndex{…}, replacement::MathOptInterface.LessThan{…})
    @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/aJZbq/src/Utilities/cachingoptimizer.jl:600
  [5] set(m::MathOptInterface.Utilities.CachingOptimizer{…}, ::MathOptInterface.ConstraintSet, cindex::MathOptInterface.ConstraintIndex{…}, set::MathOptInterface.LessThan{…})
    @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/aJZbq/src/Utilities/cachingoptimizer.jl:632
  [6] set(::MultiObjectiveAlgorithms.Optimizer, ::MathOptInterface.ConstraintSet, ::MathOptInterface.ConstraintIndex{MathOptInterface.ScalarNonlinearFunction, MathOptInterface.LessThan{Float64}}, ::MathOptInterface.LessThan{Float64})
    @ MultiObjectiveAlgorithms ~/.julia/packages/MultiObjectiveAlgorithms/ojxR6/src/MultiObjectiveAlgorithms.jl:445
  [7] optimize_multiobjective!(algorithm::MultiObjectiveAlgorithms.EpsilonConstraint, model::MultiObjectiveAlgorithms.Optimizer)
    @ MultiObjectiveAlgorithms ~/.julia/packages/MultiObjectiveAlgorithms/ojxR6/src/algorithms/EpsilonConstraint.jl:121
  [8] optimize!(model::MultiObjectiveAlgorithms.Optimizer)
    @ MultiObjectiveAlgorithms ~/.julia/packages/MultiObjectiveAlgorithms/ojxR6/src/MultiObjectiveAlgorithms.jl:544
  [9] optimize!
    @ ~/.julia/packages/MathOptInterface/aJZbq/src/Bridges/bridge_optimizer.jl:367 [inlined]
 [10] optimize!
    @ ~/.julia/packages/MathOptInterface/aJZbq/src/MathOptInterface.jl:122 [inlined]
 [11] optimize!(m::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.Bridges.LazyBridgeOptimizer{MultiObjectiveAlgorithms.Optimizer}, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}})
    @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/aJZbq/src/Utilities/cachingoptimizer.jl:321
 [12] optimize!(model::Model; ignore_optimize_hook::Bool, _differentiation_backend::MathOptInterface.Nonlinear.SparseReverseMode, kwargs::@Kwargs{})
    @ JuMP ~/.julia/packages/JuMP/7rBNn/src/optimizer_interface.jl:595
 [13] optimize!(model::Model)
    @ JuMP ~/.julia/packages/JuMP/7rBNn/src/optimizer_interface.jl:546
 [14] top-level scope
    @ REPL[80]:1

I will try to create a full example, but the external functions i have are fairly complex and difficult to replicate in a minimal way.

I also noticed the algorithm (and maybe others?) terminates without a solution if the inner optimizer exceeds max iterations. Since my problem is unconstrained, I could still use a sub-optimal. I suspect the MOA algorithm simply does not return results for this termination status?

I'm happy to take a deeper look under the hood, but I was wondering if it's obvious what the error is on inspection.

odow commented 1 month ago

I think this is actually a bug in Ipopt because it reminds me of https://github.com/jump-dev/Ipopt.jl/issues/423. I'll take a look.

jalving commented 1 month ago

Great, thanks for checking into this!