Some questions about an MOI wrapper

Specific Domains Optimization (Mathematical)
1

Hi everyone. I’m currently working on an SQP solver and now trying to use that to solve optimal control problem. The optimal control problems are trancripted to NLPs in the framework of JuMP, which means if I want to use my solver for the NLP, I better develop a wrapper so that my SQP solver can be used under JuMP.

The thing is my current solver needs constriants to be formulated in the form as
c(x) = 0 or
c(x) <= 0
However, I don’t really know how to transform a constraint of set saying Interval{Float64} to Lessthan(0.0) with MOI.
I tried to use MOI.transform, but it seems not able to real with nonlinear models that defined by MOI.Nonlinear.Model().

Inside my solver, constraints are defined in the form of equalities and then inequalites.
For example,

function eq_con(x) # if no equality constraint
return Float64
end

function ineq_con(x)
con1 = 3x[1]^2 - 2x[1]*x[2] + x[2]^2 - 1
return [con1]
end

add_con!(prob, eq_con, ineq_con)

It will be helpful if I can have some advice about this contraint handling. Whether it is better to do that inside of MOI, or making my solver itself able to detail with a general formated constraint, instead of strictly <=0 or =0.

1 Like
4

Thanks your for your reply. But when I use add_constrint to formulate this problem, I didn’t find it that straightforward when it comes to implementation.
I probably should clarify before, the problem is on implementation instead of mathmatical wise.

It is a easy fix for constraints by postprocessing evaluation, but it could be diffcult when it comes to evaluate and postprocess other related things, like jacobian. That’s why I’m think if there is a systematic way of doing that similar to MOI.transorm

5

Hi @Lester, welcome to the forum :smile:

Do you have a link to the solver you want to wrap?

See Implementing a solver interface · MathOptInterface

But not it is a nontrivial undertaking. Take a look at the source code of Ipopt.jl as an example.

You want to implement support for ScalarNonlinearFunction in LessThan and ScalarNonlinearFunction in EqualTo. JuMP will do the rest of the transformations automatically.

I can provide more advice once i know the solver.

ps you can ignore the other reply. (If you’re a real user hello and my apologies! But i think they might be a bot.)

6

@odow
I did a quick test. That’s true, I only need to support EqualTo and LessThan and it does the tranformation for me.

I’m using the wrapper of Ipopt and Optim when I was writing mine, really thank you for your work.

For the code, the solver is stil on its very early stage of development, so I haven’t uploaded the solver to Github yet. Probably I should?

I’ll put the code related to constraints here for now. It is a filterSQP like solver, and it now needs to define equality and inequalies seperately, in order to form subQP and evaluate contraint violations and etc.

function add_con!(prob, eq_con, ineq_con)
    prob.eq_con = eq_con
    prob.ineq_con = ineq_con
end

where both eq_con and ineq_con are functions, which given the current point, it returns the constraint evaluation as an vector.

function eval_cons(prob::prob, x::Vector{Float64})
    c_eq = prob.m_eq == 0 ? Float64[] : prob.eq_con(x)
    c_ineq = prob.m_ineq == 0 ? Float64[] : prob.ineq_con(x)
    return c_eq, c_ineq
end

where m_eq and m_ineq are the number of eq_con and ineq_con, defined by

prob.m_eq = length(prob.eq_con(ini))
prob.m_ineq = length(prob.ineq_con(ini))

for constraint violation h(x),

function h(prob::prob, x::Vector{Float64})
    c_eq, c_ineq = eval_cons(prob, x)

    if prob.settings.h_norm == 1
        h1 = c_eq == [] ? 0. : sum(abs(ceqi) for ceqi in c_eq) 
        h2 = c_ineq == [] ? 0. : sum(max(0, cineqi) for cineqi in c_ineq)
        return h1 + h2
    else
        error("Unsupported norm for constraint violation evaluation.")
    end
end

for jacobian, it is defined by

A = x -> begin
            A_eq = ForwardDiff.jacobian(prob.eq_con, x)
            A_ineq = ForwardDiff.jacobian(prob.ineq_con, x)
        return [A_eq; A_ineq]
    end

then I can pick selected lines based on the number of eq_con and ineq_con

7

JuMP can provide the derivatives. What oracles do you need? Is the Hessian dense or sparse?

The general approach would be:

  • Declare support for a set of constraints
  • Add them to a MOI.Nonlinear.Model
  • Turn the model an MOI.Nonlinear.Evaluator object
  • Query the oracles

The closest wrapper to look at is probably NLopt, which also supports only == 0 and <= 0:

NLopt.jl/ext/NLoptMathOptInterfaceExt.jl at master · jump-dev/NLopt.jl · GitHub.

8

Hi @Lester,
Solver developer here :slight_smile: Out of curiosity, what are the features of your SQP solver? Can’t you use an existing solver?

1 Like
9

Thanks. The NLopt wrapper is indeed very close, I almost only need some small changes and it can then work for my solver.

1 Like
10

Hi Charlie.

Honestly the solver is more about learning how a solver works, since I’ve just started my PhD and there are loads of things to explore.

Why SQP? I basically want a solver that espically favours warm-starting and use that for MPC, so an SQP-based solver can be considered. However, I reckon currently there’s no SQP solver written or interfaced to Julia for me to start with.

I think you are the one who is developing Uno? Uno should be very good I think, if one day it is in its full state. I did give it a go several months before, but the performance was not extremely good when I compare it with other solvers like MadNLP (benchmarked using OptimizationProblems.jl) and I think it was because Uno is using SLP for its Julia version?

1 Like
11

This interface to SNOPT is not maintained, but perhaps you can get it to work? GitHub - snopt/SNOPT7.jl: Julia interface for SNOPT7

Absolutely, Uno_jll.jl is only interfaced to the HiGHS LP solver at the moment, so you probably tested the SLP method (which indeed is not particularly impressive). I need to write an interface to the HiGHS QP (active-set) solver. @amontoison is also preparing precompiled binaries of the QP solver BQPD, which is our go-to QP solver at the moment.

2 Likes