Support Gurobi nonlinear expressions in the direct/persistent solvers
Summary
Gurobi v12 supports nonlinear constraints, it would be good if the direct/persistent solvers supported adding these directly to the gurobipy.Model object.
Rationale
Right now, to use this feature through Pyomo, one must use the gurobinl solver which uses the AMPL interface via an NL file.
Description
I would like to leverage the Pyomo expressions to either create an expression tree to pass this to gurobipy via Model.addGenConstrNLAdv or even create a NLExpr whatever is easiest. This requires a bit of knowledge of the Pyomo expressions and I haven't had time to dig into it yet.
Additional information
This can be used to translate nonlinear functions such as pyomo.kernel.sin as well as polynomials, hence it would solve #3249
I am very interested in this functionality and would be happy to contribute in some way. I can certainly help with testing and maybe with some coding, if someone with more knowledge of the internals of Pyomo and the gurobipy API could sketch a design.
I'm working on a direct Gurobi interface that supports MINLP as a contrib package on this branch: https://github.com/emma58/pyomo/tree/gurobi-minlp
It's not quite ready for people to beat on yet, but I can update this issue when it is.
With Gurobi v12 supporting expression trees, I believe a simpler approach is possible. Since interest in MINLP using Pyomo+Gurobi is increasing over the past weeks, I've made a start in https://github.com/ronaldvdv-gurobi/pyomo (sin/cos/exp/sum/product) and added two unit tests as well. Will continue to add the missing expression types and functions and add more tests. Please feel free to comment on this approach!
@ronaldvdv-gurobi, thanks for starting this off. I am very interested in testing this on our MINLP problems. I had previously tried using the AMPL interface, but I could not be sure that was always using the most efficient methods of interacting with Gurobi and since I have my problem coded in Pyomo I am keen to try this out. I believe that there will be significant improvements in MINLP support when Gurobi 13 is release in the next couple of months.
Yes, MINLP is very much a focus area for v13. What non-linear functions are you most interested in testing? I could definitely use a few (small) examples as testcases.
Here is an example of a very small problem: https://github.com/coin-or/SHOT/issues/166 if it is of any use. Generally, we have objective and constraint functions that have pow() or log()/exp().
@grahamsparrow-8451 I've just added that issue as a unit test and fixed a small bug.
@ronaldvdv-gurobi, thank you for looking into this! While your approach is indeed a far simpler implementation, mainly for performance reasons, this ultimately needs to be designed as a writer that walks the Pyomo expression tree and builds up the expression in gurobipy. In addition (and apologies that this isn't well advertised at the moment), the dev team has plans to eventually drop support for the solver interfaces currently in pyomo.solvers: We are working on a rewrite of all our interfaces in pyomo.contrib.solvers. The plan for the moment is to implement a direct interface to gurobi that supports MINLP in contrib.solvers, but keep it separate from the current gurobi direct rewrite (again because of performance). We will probably integrate them eventually, but are planning that as a separate step. I'm getting very close for a first draft of that interface on the branch I mentioned above. I'll merge in your tests from your branch--they're a great addition. Hoping to have a PR open this week!
Hey @emma58 , really glad to hear your branch is still WIP and could be available soon! Combining your Pyomo expertise with the testcases from my branch seems like the best of both worlds, especially if you believe a PR could be open (and merged) short-term. The reason for doing a simplified attempt was the request from several clients for a working MINLP interface to Gurobi - any approach that gets them started quickly would be great. If you want to discuss anything on this work, don't hesitate to reach out to me directly. As you've probably guessed, I work at Gurobi and would be happy to collaborate on this.
Just curious, what's the reason for implementing a new Gurobi direct interface in the contrib.solvers folder instead of extending the existing one with your walker-approach? I know from experience people are struggling with picking the right Pyomo-Gurobi interface and maintaining one more interface (until the full rewrite is completed and merged) can be tricky. How would the timeline for the full rewrites compare to the availability of your PR? Is it worth having a separate MINLP interface?
@ronaldvdv-gurobi, the reason to do a new one is that the current gurobi_direct_v2 only supports linear problems and can handle very large linear problems since it is aware of templatized Pyomo expressions. So I'm hesitant to pollute it with catching the case where things go nonlinear at the moment since we want to do that in such a way that we don't lose these gains. We're also getting close to expanding all that to include quadratic, so I think it makes sense to get that in first. However, once we have, all the pieces will exist, so perhaps it's not such a big lift to integrate them all together? We'd be very happy to discuss that with you at some point--and we do love help! And in the long term, I agree that a separate MINLP interface is clunky, but I do think that's the path to getting something usable in the short (next few week) time horizon.
@ronaldvdv-gurobi, if you have a few minutes, could I have your help? I merged in your tests (on this branch), and 4 of them are failing with my current implementation: In particular, I dug into test_gurobi_minlp_sincosexp and I can't find what about the resulting model is incorrect. The Pyomo model is:
m = ConcreteModel(name="test")
m.x = Var(bounds=(-1, 4))
m.o = Objective(expr=sin(m.x) + cos(2 * m.x) + 1)
m.c = Constraint(expr=0.25 * exp(m.x) - m.x <= 0)
and if I write the lp file for the resulting Gurobi model (using this code):
from pyomo.contrib.solver.solvers.gurobi_direct_minlp import GurobiDirectMINLP
from pyomo.opt import WriterFactory
grb_model, var_map, obj, grb_cons, pyo_cons = WriterFactory(
'gurobi_minlp'
).write(m, symbolic_solver_labels=True)
grb_model.write('debug.lp')
I get:
\ LP format - for model browsing. Use MPS format to capture full model detail.
Minimize
0 x + C1
Subject To
R0: C2 <= 0
Bounds
-1 <= x <= 4
General Constraints
\ C1 = (sin(x) + cos(2 * x)) + 1
GC0: C1 = NL : ( PLUS , -1 , -1 ) ( PLUS , -1 , 0 ) ( SIN , -1 , 1 )
( VARIABLE , x , 2 ) ( COS , -1 , 1 ) ( MULTIPLY , -1 , 4 )
( CONSTANT , 2 , 5 ) ( VARIABLE , x , 5 ) ( CONSTANT , 1 , 0 )
\ C2 = (0.25 * exp(x)) + (-1 * x)
GC1: C2 = NL : ( PLUS , -1 , -1 ) ( MULTIPLY , -1 , 0 )
( CONSTANT , 0.25 , 1 ) ( EXP , -1 , 1 ) ( VARIABLE , x , 3 )
( MULTIPLY , -1 , 0 ) ( CONSTANT , -1 , 5 ) ( VARIABLE , x , 5 )
End
And indeed if I call optimize on that gurobi model, I get the same wrong answer as the Pyomo interface is currently getting. But I'm not spotting the error?
(Also, if you'd like to take this offline, please email [email protected], and I can connect with you directly.)
Yep, I've had the same error: it's (likely) because the auxiliary variables generated by the solver interface, have their default lowerbound of 0 instead of -GRB.INFINITY . To confirm I made that change in my own interface and get the same (wrong) answer for the model above. And likewise, with infinite bounds, inspecting the solution shows that it requires the auxiliary (first) variable to have a negative value. Could you try if introducing infinite bounds solves the issue with your interface?