diff --git a/docs/source/index.rst b/docs/source/index.rst index 52c02b4..3ab2c6c 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -78,6 +78,7 @@ If you encounter issues with Gurobi or ``gurobipy`` please contact performance naming + nonlinear advanced typing diff --git a/docs/source/nonlinear.rst b/docs/source/nonlinear.rst new file mode 100644 index 0000000..d314f19 --- /dev/null +++ b/docs/source/nonlinear.rst @@ -0,0 +1,119 @@ +Adding Nonlinear Constraints +============================ + +Gurobi 12 supports adding nonlinear constraints, using the ``NLExpr`` object to +capture nonlinear expressions. ``gurobipy-pandas`` supports adding a ``Series`` +of nonlinear constraints to a model via ``add_constrs``. Note that ``gurobipy`` +only supports nonlinear constraints of the form :math:`y = f(\bar{x})` and +``gurobipy-pandas`` applies the same restriction. + +.. note:: + + To add nonlinear constraints, you must have at least ``gurobipy`` version + 12.0.0 installed. + +Nonlinear Equality Constraints +------------------------------ + +This example builds the constraint set :math:`y_i = \log(\frac{1}{x_i})`, for +each :math:`i` in the index. + +.. doctest:: [nonlinear] + + >>> import pandas as pd + >>> import gurobipy as gp + >>> from gurobipy import GRB + >>> from gurobipy import nlfunc + >>> import gurobipy_pandas as gppd + >>> gppd.set_interactive() + + >>> model = gp.Model() + >>> index = pd.RangeIndex(5) + >>> x = gppd.add_vars(model, index, lb=1.0, name="x") + >>> y = gppd.add_vars(model, index, lb=-GRB.INFINITY, name="y") + +You can create nonlinear expressions using standard Python operators. A +nonlinear expression is any expression which is not a polynomial of degree at +most 2. For example, dividing a constant by a series of ``Var`` objects produces +a series of nonlinear expressions. + +.. doctest:: [nonlinear] + + >>> 1 / x + 0 1.0 / x[0] + 1 1.0 / x[1] + 2 1.0 / x[2] + 3 1.0 / x[3] + 4 1.0 / x[4] + Name: x, dtype: object + +The ``nlfunc`` module of gurobipy is used to create nonlinear expressions +involving mathematical functions. You can use ``apply`` to construct a series +representing the logarithm of the above result. + +.. doctest:: [nonlinear] + + >>> (1 / x).apply(nlfunc.log) + 0 log(1.0 / x[0]) + 1 log(1.0 / x[1]) + 2 log(1.0 / x[2]) + 3 log(1.0 / x[3]) + 4 log(1.0 / x[4]) + Name: x, dtype: object + +This series of expressions can be added as constraints using ``add_constrs``. +Note that you can only provide nonlinear constraints in the form :math:`y = +f(\bar{x})` where :math:`f` may be a univariate or multivariate compound +function. Therefore the ``lhs`` argument must be a series of ``Var`` objects, +while the ``sense`` argument must be ``GRB.EQUAL``. + +.. doctest:: [nonlinear] + + >>> gppd.add_constrs(model, y, GRB.EQUAL, (1 / x).apply(nlfunc.log), name="log_x") + 0 + 1 + 2 + 3 + 4 + Name: log_x, dtype: object + +Nonlinear Inequality Constraints +-------------------------------- + +As noted above, nonlinear constraints can only be provided in the form :math:`y= +f(\bar{x})`. To formulate inequality constraints, you must introduce bounded +intermediate variables. The following example formulates the set of constraints +:math:`\log(x_i^2 + 1) \le 2` by introducing intermediate variables :math:`z_i` +with no lower bound and an upper bound of :math:`2.0`. + +.. doctest:: [nonlinear] + + >>> z = gppd.add_vars(model, index, lb=-GRB.INFINITY, ub=2.0, name="z") + >>> constrs = gppd.add_constrs(model, z, GRB.EQUAL, (x**2 + 1).apply(nlfunc.log)) + +To see the effect of this constraint, you can set a maximization objective +:math:`\sum_{i=0}^{4} x_i` and solve the resulting model. You will find that the +original variables :math:`x_i` are maximized to :math:`\sqrt{e^2 - 1}` in +the optimal solution. The intermediate variables :math:`z_i` are at their upper +bounds, indicating that the constraint is satisfied with equality. + +.. doctest:: [nonlinear] + + >>> model.setObjective(x.sum(), sense=GRB.MAXIMIZE) + >>> model.optimize() # doctest: +ELLIPSIS + Gurobi Optimizer ... + Optimal solution found ... + >>> x.gppd.X.round(3) + 0 2.528 + 1 2.528 + 2 2.528 + 3 2.528 + 4 2.528 + Name: x, dtype: float64 + >>> z.gppd.X.round(3) + 0 2.0 + 1 2.0 + 2 2.0 + 3 2.0 + 4 2.0 + Name: z, dtype: float64 diff --git a/src/gurobipy_pandas/constraints.py b/src/gurobipy_pandas/constraints.py index b2a0a12..d7ce0f3 100644 --- a/src/gurobipy_pandas/constraints.py +++ b/src/gurobipy_pandas/constraints.py @@ -13,6 +13,7 @@ from gurobipy_pandas.util import create_names, gppd_global_options CONSTRAINT_SENSES = frozenset([GRB.LESS_EQUAL, GRB.EQUAL, GRB.GREATER_EQUAL]) +GUROBIPY_MAJOR_VERSION, *_ = gp.gurobi.version() def add_constrs_from_dataframe( @@ -138,6 +139,15 @@ def _add_constr(model, lhs, sense, rhs, name): name = "" if not isinstance(sense, str) or sense[0] not in CONSTRAINT_SENSES: raise ValueError(f"'{sense}' is not a valid constraint sense") + if GUROBIPY_MAJOR_VERSION >= 12: + # Check for nonlinear constraints; accept only y = f(x) + if isinstance(lhs, gp.NLExpr): + raise ValueError("Nonlinear constraints must be in the form y = f(x)") + if isinstance(rhs, gp.NLExpr): + if isinstance(lhs, gp.Var) and sense == GRB.EQUAL: + return model.addGenConstrNL(lhs, rhs, name=name) + else: + raise ValueError("Nonlinear constraints must be in the form y = f(x)") if isinstance(lhs, gp.QuadExpr) or isinstance(rhs, gp.QuadExpr): return model.addQConstr(lhs, sense, rhs, name=name) return model.addLConstr(lhs, sense, rhs, name=name) diff --git a/tests/test_api.py b/tests/test_api.py index 02a8f3f..f6a43c1 100644 --- a/tests/test_api.py +++ b/tests/test_api.py @@ -3,6 +3,10 @@ tests of data types, errors, etc, are done on the lower-level functions. """ +import math +import unittest + +import gurobipy as gp import pandas as pd from gurobipy import GRB from pandas.testing import assert_index_equal, assert_series_equal @@ -10,6 +14,8 @@ import gurobipy_pandas as gppd from tests.utils import GurobiModelTestCase +GUROBIPY_MAJOR_VERSION, *_ = gp.gurobi.version() + class TestAddVars(GurobiModelTestCase): def test_from_dataframe(self): @@ -269,6 +275,136 @@ def test_sense_series(self): self.assertEqual(constr.RHS, -1.0) +@unittest.skipIf( + GUROBIPY_MAJOR_VERSION < 12, + "Nonlinear constraints are only supported for Gurobi 12 and later", +) +class TestNonlinear(GurobiModelTestCase): + def assert_approx_equal(self, value, expected, tolerance=1e-6): + difference = abs(value - expected) + self.assertLessEqual(difference, tolerance) + + def test_log(self): + # max sum y_i + # s.t. y_i = log(x_i) + # 1.0 <= x <= 2.0 + from gurobipy import nlfunc + + index = pd.RangeIndex(3) + x = gppd.add_vars(self.model, index, lb=1.0, ub=2.0, name="x") + y = gppd.add_vars(self.model, index, obj=1.0, name="y") + self.model.ModelSense = GRB.MAXIMIZE + + gppd.add_constrs(self.model, y, GRB.EQUAL, x.apply(nlfunc.log), name="log_x") + + self.model.optimize() + self.assert_approx_equal(self.model.ObjVal, 3 * math.log(2.0)) + + def test_inequality(self): + # max sum x_i + # s.t. log2(x_i^2 + 1) <= 2.0 + # 0.0 <= x <= 1.0 + # + # Formulated as + # + # max sum x_i + # s.t. log2(x_i^2 + 1) == z_i + # 0.0 <= x <= 1.0 + # -GRB.INFINITY <= z_i <= 2 + from gurobipy import nlfunc + + index = pd.RangeIndex(3) + x = gppd.add_vars(self.model, index, name="x") + z = gppd.add_vars(self.model, index, lb=-GRB.INFINITY, ub=2.0, name="z") + self.model.setObjective(x.sum(), sense=GRB.MAXIMIZE) + + gppd.add_constrs(self.model, z, GRB.EQUAL, (x**2 + 1).apply(nlfunc.log)) + + self.model.optimize() + self.model.write("model.lp") + self.model.write("model.sol") + + x_sol = x.gppd.X + z_sol = z.gppd.X + for i in range(3): + self.assert_approx_equal(x_sol[i], math.sqrt(math.exp(2.0) - 1)) + self.assert_approx_equal(z_sol[i], 2.0) + + def test_wrong_usage(self): + index = pd.RangeIndex(3) + x = gppd.add_vars(self.model, index, name="x") + y = gppd.add_vars(self.model, index, name="y") + + with self.assertRaisesRegex( + gp.GurobiError, "Objective must be linear or quadratic" + ): + self.model.setObjective((x / y).sum()) + + with self.assertRaisesRegex( + ValueError, "Nonlinear constraints must be in the form" + ): + gppd.add_constrs(self.model, y, GRB.LESS_EQUAL, x**4) + + with self.assertRaisesRegex( + ValueError, "Nonlinear constraints must be in the form" + ): + gppd.add_constrs(self.model, y + x**4, GRB.EQUAL, 1) + + with self.assertRaisesRegex( + ValueError, "Nonlinear constraints must be in the form" + ): + gppd.add_constrs(self.model, y**4, GRB.EQUAL, x) + + with self.assertRaisesRegex( + ValueError, "Nonlinear constraints must be in the form" + ): + x.to_frame().gppd.add_constrs(self.model, "x**3 == 1", name="bad") + + def test_eval(self): + index = pd.RangeIndex(3) + df = ( + gppd.add_vars(self.model, index, name="x") + .to_frame() + .gppd.add_vars(self.model, name="y") + .gppd.add_constrs(self.model, "y == x**3", name="nlconstr") + ) + + self.model.update() + for row in df.itertuples(index=False): + self.assert_nlconstr_equal( + row.nlconstr, + row.y, + [ + (GRB.OPCODE_POW, -1.0, -1), + (GRB.OPCODE_VARIABLE, row.x, 0), + (GRB.OPCODE_CONSTANT, 3.0, 0), + ], + ) + + def test_frame(self): + from gurobipy import nlfunc + + index = pd.RangeIndex(3) + df = ( + gppd.add_vars(self.model, index, name="x") + .to_frame() + .gppd.add_vars(self.model, name="y") + .assign(exp_x=lambda df: df["x"].apply(nlfunc.exp)) + .gppd.add_constrs(self.model, "y", GRB.EQUAL, "exp_x", name="nlconstr") + ) + + self.model.update() + for row in df.itertuples(index=False): + self.assert_nlconstr_equal( + row.nlconstr, + row.y, + [ + (GRB.OPCODE_EXP, -1.0, -1), + (GRB.OPCODE_VARIABLE, row.x, 0), + ], + ) + + class TestDataValidation(GurobiModelTestCase): # Test that we throw some more informative errors, instead of obscure # ones from the underlying gurobipy library diff --git a/tests/test_constraints.py b/tests/test_constraints.py index 46fd2f5..92db2cb 100644 --- a/tests/test_constraints.py +++ b/tests/test_constraints.py @@ -15,6 +15,8 @@ from .utils import GurobiModelTestCase +GUROBIPY_MAJOR_VERSION, *_ = gp.gurobi.version() + class TestAddConstrsFromDataFrame(GurobiModelTestCase): def test_scalar_rhs(self): @@ -317,6 +319,95 @@ def test_nonpython_columnnames(self): self.assertEqual(row.getCoeff(0), 1.0) +@unittest.skipIf( + GUROBIPY_MAJOR_VERSION < 12, + "Nonlinear constraints are only supported for Gurobi 12 and later", +) +class TestNonlinearFromDataFrame(GurobiModelTestCase): + def test_varseries_nlseries(self): + # y_i == exp(x_i) + from gurobipy import nlfunc + + index = pd.RangeIndex(5) + df = pd.DataFrame( + dict( + x=add_vars_from_index(self.model, index, name="x"), + y=add_vars_from_index(self.model, index, name="y"), + ) + ) + + df["expr"] = df["x"].apply(nlfunc.exp) + + constrs = add_constrs_from_dataframe( + self.model, df, "y", GRB.EQUAL, "expr", name="nonlinear" + ) + + # Constraints added to model + self.model.update() + self.assertEqual(self.model.NumGenConstrs, 5) + + # Returned series has the right structure + self.assertIsInstance(constrs, pd.Series) + self.assertEqual(constrs.name, "nonlinear") + assert_index_equal(constrs.index, index) + + # Check data in model + for ind, constr in constrs.items(): + self.assertIsInstance(constr, gp.GenConstr) + self.assertEqual(constr.GenConstrName, f"nonlinear[{ind}]") + self.assertEqual(constr.GenConstrType, GRB.GENCONSTR_NL) + self.assert_nlconstr_equal( + constr, + df.loc[ind, "y"], + [ + (GRB.OPCODE_EXP, -1.0, -1), + (GRB.OPCODE_VARIABLE, df.loc[ind, "x"], 0), + ], + ) + + def test_left_nlexpr(self): + # exp(x_i) == y_i + from gurobipy import nlfunc + + index = pd.RangeIndex(5) + df = pd.DataFrame( + dict( + x=add_vars_from_index(self.model, index, name="x"), + y=add_vars_from_index(self.model, index, name="y"), + ) + ) + + df["expr"] = df["x"].apply(nlfunc.exp) + + with self.assertRaisesRegex( + ValueError, "Nonlinear constraints must be in the form" + ): + add_constrs_from_dataframe( + self.model, df, "expr", GRB.EQUAL, "y", name="nonlinear" + ) + + def test_inequality(self): + # y_i <= exp(x_i) + from gurobipy import nlfunc + + index = pd.RangeIndex(5) + df = pd.DataFrame( + dict( + x=add_vars_from_index(self.model, index, name="x"), + y=add_vars_from_index(self.model, index, name="y"), + ) + ) + + df["expr"] = df["x"].apply(nlfunc.exp) + + with self.assertRaisesRegex( + ValueError, "Nonlinear constraints must be in the form" + ): + add_constrs_from_dataframe( + self.model, df, "y", GRB.LESS_EQUAL, "expr", name="nonlinear" + ) + + class TestAddConstrsFromSeries(GurobiModelTestCase): def test_bothseries(self): # linear series <= linear series diff --git a/tests/test_docs.py b/tests/test_docs.py index bfabe70..af9e0b4 100644 --- a/tests/test_docs.py +++ b/tests/test_docs.py @@ -31,5 +31,7 @@ def load_tests(loader, tests, ignore): here = pathlib.Path(__file__).parent docs = here.joinpath("../docs/source") for docfile in docs.rglob("*.rst"): + if docfile.stem == "nonlinear" and GUROBIPY_MAJOR_VERSION < 12: + continue tests.addTests(doctest.DocFileSuite(str(docfile.relative_to(here)))) return tests diff --git a/tests/utils.py b/tests/utils.py index 70aeac3..2643fc6 100644 --- a/tests/utils.py +++ b/tests/utils.py @@ -1,6 +1,7 @@ import unittest import gurobipy as gp +from gurobipy import GRB class GurobiModelTestCase(unittest.TestCase): @@ -36,3 +37,18 @@ def assert_quadexpr_equal(self, expr1, expr2): self.assertTrue(expr1.getVar1(i).sameAs(expr2.getVar1(i))) self.assertTrue(expr1.getVar2(i).sameAs(expr2.getVar2(i))) self.assertEqual(expr1.getCoeff(i), expr2.getCoeff(i)) + + def assert_nlconstr_equal(self, genconstr, resvar_true, tree): + resvar, opcode_array, data_array, parent_array = self.model.getGenConstrNLAdv( + genconstr + ) + self.assertIs(resvar, resvar_true) + for opcode, data, parent, (opcode_true, data_true, parent_true) in zip( + opcode_array, data_array, parent_array, tree + ): + self.assertEqual(opcode, opcode_true) + if opcode == GRB.OPCODE_VARIABLE: + self.assertIs(data, data_true) + else: + self.assertEqual(data, data_true) + self.assertEqual(parent, parent_true)