Add AbstractOracle to make a type-stable interface

[odow]

In addition, add a CI job to test JuliaC with trimming to ensure that we maintain support for trimming going forward.

This is breaking in the sense that we are adding a new type parameter. It could also be argued that this is a non-breaking change because the only case where it would break is someone dispatching on ::Type{IpoptProblem} which I can't see any reason for doing.

Closes #532
Closes #534

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 100.00%. Comparing base (1f460a9) to head (34ff621).

Additional details and impacted files

@@            Coverage Diff            @@
##            master      #535   +/-   ##
=========================================
  Coverage   100.00%   100.00%           
=========================================
  Files            5         5           
  Lines         1192      1219   +27     
=========================================
+ Hits          1192      1219   +27     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[amontoison]

@odow Can you explain why it is more stable with a FunctionOracle ?
I think we we should add some tests related to type inference.

[odow]

The FunctionOracle is not type stable. Using JuliaC requires a separate type-stable <:AbstractOracle

[odow]

I also don't know that we should do this. I just wanted to see if it works

[odow]

If we're making changes, I'd rather do this than #532.

[odow]

Okay, I think I've managed to make this a strictly non-breaking feature addition.

It would be breaking only if someone was:

• accessing the prob.eval_f, prob.eval_g, etc. fields directly. This is non-breaking because accessing undocumented fields of a struct is not public behavior
• Constructing IpoptProblem(args...) directly instead of CreateIpoptProblem, which is also not public behaviour.

Decision needed:

If we're making changes like this, there's an opportunity for us to expose all of the underlying C arguments.

This would resolve #313 and #322

[odow]

I should also try updating the MOI interface

[odow]

@hfytr do you want to try this branch with your JuliaC code?

[odow]

I should also try updating the MOI interface

I looked into this. There's not much point because there are a few type instabilities we can't get rid of:

• The evaluator in MOI.NLPBlockData is abstract: https://github.com/jump-dev/MathOptInterface.jl/blob/65e4a61fda04845bbdfc646bb980a2cc6491a22e/src/nlp.jl#L176
• The AD backend is abstract:

  Ipopt.jl/ext/IpoptMathOptInterfaceExt/MOI_wrapper.jl

  Line 59 in 1f460a9

  ad_backend::MOI.Nonlinear.AbstractAutomaticDifferentiation

So supporting JuliaC in MOI is going to require some mechanism for freezing methods at compile time.

[odow]

We could do something like this, but now the inner is not type stable.

diff --git a/ext/IpoptMathOptInterfaceExt/MOI_wrapper.jl b/ext/IpoptMathOptInterfaceExt/MOI_wrapper.jl
index 819027f..9e82019 100644
--- a/ext/IpoptMathOptInterfaceExt/MOI_wrapper.jl
+++ b/ext/IpoptMathOptInterfaceExt/MOI_wrapper.jl
@@ -36,7 +36,7 @@ end
 Create a new Ipopt optimizer.
 """
 mutable struct Optimizer <: MOI.AbstractOptimizer
-    inner::Union{Nothing,Ipopt.IpoptProblem}
+    inner::Union{Nothing,Ipopt.Problem}
     name::String
     invalid_model::Bool
     silent::Bool
@@ -773,13 +773,14 @@ function MOI.get(
 ) where {F<:MOI.VectorOfVariables,S<:MOI.VectorNonlinearOracle{Float64}}
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
+    inner = _inner(model)
     sign = -_dual_multiplier(model)
     f, s = model.vector_nonlinear_oracle_constraints[ci.value]
-    λ = model.inner.mult_g[row(model, ci)]
+    λ = inner.mult_g[row(model, ci)]
     J = Tuple{Int,Int}[]
     _jacobian_structure(J, 0, f, s)
     J_val = zeros(length(J))
-    _eval_constraint_jacobian(J_val, 0, model.inner.x, f, s)
+    _eval_constraint_jacobian(J_val, 0, inner.x, f, s)
     dual = zeros(MOI.dimension(s.set))
     # dual = λ' * J(x)
     col_to_index = Dict(x.value => j for (j, x) in enumerate(f.variables))
@@ -796,7 +797,8 @@ function MOI.get(
 ) where {F<:MOI.VectorOfVariables,S<:MOI.VectorNonlinearOracle{Float64}}
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    return -_dual_multiplier(model) * model.inner.mult_g[row(model, ci)]
+    inner = _inner(model)
+    return -_dual_multiplier(model) * inner.mult_g[row(model, ci)]
 end
 
 function MOI.supports(
@@ -1270,7 +1272,32 @@ end
 
 ### MOI.optimize!
 
-function _eval_jac_g_cb(model, x, rows, cols, values)
+struct Oracle{T} <: Ipopt.AbstractOracle
+    model::Optimizer
+    evaluator::T
+
+    function Oracle(model::Optimizer)
+        evaluator = model.nlp_data.evaluator
+        return new{typeof(evaluator)}(model, evaluator)
+    end
+end
+
+function Ipopt.eval_f(oracle::Oracle, x)
+    return MOI.eval_objective(oracle.model, x)
+end
+
+function Ipopt.eval_g(oracle::Oracle, x, g)
+    MOI.eval_constraint(oracle.model, g, x)
+    return
+end
+
+function Ipopt.eval_grad_f(oracle::Oracle, x, grad_f)
+    MOI.eval_objective_gradient(oracle.model, grad_f, x)
+    return
+end
+
+function Ipopt.eval_jac_g(oracle::Oracle, x, rows, cols, values)
+    model = oracle.model
     if values === nothing
         for i in 1:length(model.jacobian_sparsity)
             rows[i], cols[i] = model.jacobian_sparsity[i]
@@ -1281,7 +1308,12 @@ function _eval_jac_g_cb(model, x, rows, cols, values)
     return
 end
 
-function _eval_h_cb(model, x, rows, cols, obj_factor, lambda, values)
+function Ipopt.has_eval_h(oracle::Oracle)
+    return oracle.model.hessian_sparsity !== nothing
+end
+
+function Ipopt.eval_h(oracle::Oracle, x, rows, cols, obj_factor, lambda, values)
+    model = oracle.model
     if values === nothing
         for (i, v) in enumerate(model.hessian_sparsity::Vector{Tuple{Int,Int}})
             rows[i], cols[i] = v
@@ -1292,9 +1324,25 @@ function _eval_h_cb(model, x, rows, cols, obj_factor, lambda, values)
     return
 end
 
-function _setup_inner(model::Optimizer)::Ipopt.IpoptProblem
+function Ipopt.eval_intermediate(
+    ::Ipopt.Problem{<:Oracle},
+    oracle::Oracle,
+    args...,
+)
+    model = oracle.model
+    # iter_count is args[2]
+    model.barrier_iterations = args[2]
+    if model.callback !== nothing
+        return model.callback(args...)
+    end
+    return true
+end
+
+_inner(model::Optimizer) = model.inner::Ipopt.Problem{<:Oracle}
+
+function _setup_inner(model::Optimizer)
     if !model.needs_new_inner
-        return model.inner
+        return _inner(model)
     end
     g_L, g_U = copy(model.qp_data.g_L), copy(model.qp_data.g_U)
     for (_, s) in model.vector_nonlinear_oracle_constraints
@@ -1305,10 +1353,8 @@ function _setup_inner(model::Optimizer)::Ipopt.IpoptProblem
         push!(g_L, bound.lower)
         push!(g_U, bound.upper)
     end
-    function eval_h_cb(x, rows, cols, obj_factor, lambda, values)
-        return _eval_h_cb(model, x, rows, cols, obj_factor, lambda, values)
-    end
-    has_hessian = model.hessian_sparsity !== nothing
+    oracle = Oracle(model)
+    has_hessian = Ipopt.has_eval_h(oracle)
     model.inner = Ipopt.CreateIpoptProblem(
         length(model.variables.lower),
         model.variables.lower,
@@ -1318,14 +1364,9 @@ function _setup_inner(model::Optimizer)::Ipopt.IpoptProblem
         g_U,
         length(model.jacobian_sparsity),
         has_hessian ? length(model.hessian_sparsity) : 0,
-        (x) -> MOI.eval_objective(model, x),
-        (x, g) -> MOI.eval_constraint(model, g, x),
-        (x, grad_f) -> MOI.eval_objective_gradient(model, grad_f, x),
-        (x, rows, cols, values) ->
-            _eval_jac_g_cb(model, x, rows, cols, values),
-        has_hessian ? eval_h_cb : nothing,
+        oracle,
     )
-    inner = model.inner::Ipopt.IpoptProblem
+    inner = _inner(model)
     if model.sense == MOI.MIN_SENSE
         Ipopt.AddIpoptNumOption(inner, "obj_scaling_factor", 1.0)
     elseif model.sense == MOI.MAX_SENSE
@@ -1351,17 +1392,9 @@ function _setup_inner(model::Optimizer)::Ipopt.IpoptProblem
             Ipopt.AddIpoptStrOption(inner, "hessian_constant", "yes")
         end
     end
-    function _moi_callback(args...)
-        # iter_count is args[2]
-        model.barrier_iterations = args[2]
-        if model.callback !== nothing
-            return model.callback(args...)
-        end
-        return true
-    end
-    Ipopt.SetIntermediateCallback(inner, _moi_callback)
+    Ipopt.SetIntermediateCallback(inner)
     model.needs_new_inner = false
-    return model.inner
+    return _inner(model)
 end
 
 function _setup_model(model::Optimizer)
@@ -1531,7 +1564,8 @@ function MOI.get(model::Optimizer, ::MOI.TerminationStatus)
     elseif model.inner === nothing
         return MOI.OPTIMIZE_NOT_CALLED
     end
-    status, _ = _STATUS_CODES[Ipopt.ApplicationReturnStatus(model.inner.status)]
+    inner = _inner(model)
+    status, _ = _STATUS_CODES[Ipopt.ApplicationReturnStatus(inner.status)]
     return status
 end
 
@@ -1543,13 +1577,15 @@ function MOI.get(model::Optimizer, ::MOI.RawStatusString)
     elseif model.inner === nothing
         return "Optimize not called"
     end
-    return string(Ipopt.ApplicationReturnStatus(model.inner.status))
+    inner = _inner(model)
+    return string(Ipopt.ApplicationReturnStatus(inner.status))
 end
 
 ### MOI.PrimalStatus
 
 function _manually_evaluated_primal_status(model::Optimizer)
-    x, g = model.inner.x, model.inner.g
+    inner = _inner(model)
+    x, g = inner.x, inner.g
     x_L, x_U = model.variables.lower, model.variables.upper
     g_L, g_U = copy(model.qp_data.g_L), copy(model.qp_data.g_U)
     # Assuming constraints are guaranteed to be in the order:
@@ -1582,7 +1618,8 @@ function MOI.get(model::Optimizer, attr::MOI.PrimalStatus)
     if !(1 <= attr.result_index <= MOI.get(model, MOI.ResultCount()))
         return MOI.NO_SOLUTION
     end
-    _, status = _STATUS_CODES[Ipopt.ApplicationReturnStatus(model.inner.status)]
+    inner = _inner(model)
+    _, status = _STATUS_CODES[Ipopt.ApplicationReturnStatus(inner.status)]
     if status == MOI.UNKNOWN_RESULT_STATUS
         # Not sure. RestorationFailure can terminate at a feasible (but
         # non-stationary) point.
@@ -1597,7 +1634,8 @@ function MOI.get(model::Optimizer, attr::MOI.DualStatus)
     if !(1 <= attr.result_index <= MOI.get(model, MOI.ResultCount()))
         return MOI.NO_SOLUTION
     end
-    _, status = _STATUS_CODES[Ipopt.ApplicationReturnStatus(model.inner.status)]
+    inner = _inner(model)
+    _, status = _STATUS_CODES[Ipopt.ApplicationReturnStatus(inner.status)]
     return status
 end
 
@@ -1613,7 +1651,8 @@ MOI.get(model::Optimizer, ::MOI.BarrierIterations) = model.barrier_iterations
 
 function MOI.get(model::Optimizer, attr::MOI.ObjectiveValue)
     MOI.check_result_index_bounds(model, attr)
-    return model.inner.obj_val
+    inner = _inner(model)
+    return inner.obj_val
 end
 
 ### MOI.VariablePrimal
@@ -1629,7 +1668,8 @@ function MOI.get(
         p = model.parameters[vi]
         return model.nlp_model[p]
     end
-    return model.inner.x[Ipopt.column(vi)]
+    inner = _inner(model)
+    return inner.x[Ipopt.column(vi)]
 end
 
 ### MOI.ConstraintPrimal
@@ -1670,7 +1710,8 @@ function MOI.get(
 }
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    return model.inner.g[row(model, ci)]
+    inner = _inner(model)
+    return inner.g[row(model, ci)]
 end
 
 function MOI.get(
@@ -1680,7 +1721,8 @@ function MOI.get(
 )
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    return model.inner.x[ci.value]
+    inner = _inner(model)
+    return inner.x[ci.value]
 end
 
 ### MOI.ConstraintDual
@@ -1701,7 +1743,8 @@ function MOI.get(
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
     s = -_dual_multiplier(model)
-    return s * model.inner.mult_g[row(model, ci)]
+    inner = _inner(model)
+    return s * inner.mult_g[row(model, ci)]
 end
 
 function MOI.get(
@@ -1711,7 +1754,8 @@ function MOI.get(
 )
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    rc = model.inner.mult_x_L[ci.value] - model.inner.mult_x_U[ci.value]
+    inner = _inner(model)
+    rc = inner.mult_x_L[ci.value] - inner.mult_x_U[ci.value]
     return min(0.0, _dual_multiplier(model) * rc)
 end
 
@@ -1722,7 +1766,8 @@ function MOI.get(
 )
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    rc = model.inner.mult_x_L[ci.value] - model.inner.mult_x_U[ci.value]
+    inner = _inner(model)
+    rc = inner.mult_x_L[ci.value] - inner.mult_x_U[ci.value]
     return max(0.0, _dual_multiplier(model) * rc)
 end
 
@@ -1733,7 +1778,8 @@ function MOI.get(
 )
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    rc = model.inner.mult_x_L[ci.value] - model.inner.mult_x_U[ci.value]
+    inner = _inner(model)
+    rc = inner.mult_x_L[ci.value] - inner.mult_x_U[ci.value]
     return _dual_multiplier(model) * rc
 end
 
@@ -1744,7 +1790,8 @@ function MOI.get(
 )
     MOI.check_result_index_bounds(model, attr)
     MOI.throw_if_not_valid(model, ci)
-    rc = model.inner.mult_x_L[ci.value] - model.inner.mult_x_U[ci.value]
+    inner = _inner(model)
+    rc = inner.mult_x_L[ci.value] - inner.mult_x_U[ci.value]
     return _dual_multiplier(model) * rc
 end
 
@@ -1753,7 +1800,8 @@ end
 function MOI.get(model::Optimizer, attr::MOI.NLPBlockDual)
     MOI.check_result_index_bounds(model, attr)
     s = -_dual_multiplier(model)
-    return s .* model.inner.mult_g[(length(model.qp_data)+1):end]
+    inner = _inner(model)
+    return s .* inner.mult_g[(length(model.qp_data)+1):end]
 end
 
 ### Ipopt.CallbackFunction
@@ -1801,5 +1849,6 @@ function MOI.get(
     ::MOI.CallbackVariablePrimal,
     x::MOI.VariableIndex,
 )
-    return model.inner.x[Ipopt.column(x)]
+    inner = _inner(model)
+    return inner.x[Ipopt.column(x)]
 end

[odow]

@amontoison