Documentation


Constraints

A constraint in Gurobi captures a restriction on the values that a set of variables may take. The simplest example is a linear constraint, which states that a linear expression on a set of variables take a value that is either less-than-or-equal, greater-than-or-equal, or equal to another linear expression. Recall that Gurobi works in finite-precision arithmetic, so constraints are only satisfied to tolerances. Tolerances can be tightened to reduce such violations, but there are limits to how small the violations can be -- errors are inherent in floating-point arithmetic.

The available constraint types are linear, SOS, quadratic (both convex and non-convex), and general.

Linear Constraints

A linear constraint allows you to restrict the value of a linear expression. For example, you may require that any feasible solution satisfies the constraint <span>$</span>3 x + 4 y \leq 5z<span>$</span>. Note that the matrix-oriented Gurobi API's (C, MATLAB, and R) require the right-hand side of a linear constraint to be a constant, while the object-oriented APIs (C++, Java, .NET, and Python) allow arbitrary linear expressions on both sides of the comparator.

The computed solution should satisfy the stated constraint to within FeasibilityTol (although it may not in cases of numerical ill-conditioning -- we'll discuss this shortly).

Gurobi supports a limited set of comparators. Specifically, you can constrain an expression to be less-than-or-equal, greater-than-or-equal, or equal another. We do not support strict less-than, strict greater-than, or not-equal comparators. While these other comparators may seem appropriate for mathematical programming, we exclude them to avoid potential confusion related to numerical tolerances. Consider a simple example of a strict inequality constraint on a pair of continuous variables: <span>$</span>x > y<span>$</span>. How large would <span>$</span>x-y<span>$</span> need to be in order to satisfy the constraint? Rather than trying to embed a subtle and potentially confusing strategy for handling such constraints into the solver, we've chosen not to support them instead.

SOS Constraints

A Special-Ordered Set, or SOS constraint, is a highly specialized constraint that places restrictions on the values that variables in a given list can take. There are two types of SOS constraints. In an SOS constraint of type 1 (an SOS1 constraint), at most one variable in the specified list is allowed to take a non-zero value. In an SOS constraint of type 2 (an SOS2 constraint), at most two variables in the specified, ordered list are allowed to take a non-zero value, and those non-zero variables must be contiguous in the list. The variables in an SOS constraint can be continuous, integer, or binary.

Again, tolerances play an important role in SOS constraints. Specifically, variables that take values less than IntFeasTol (in absolute value) are considered to be zero for the purposes of determining whether an SOS constraint is satisfied.

An SOS constraint is described using a list of variables and a list of corresponding weights. While the weights have historically had intuitive meanings associated with them, we simply use them to order the list of variables. The weights should be unique. This is especially important for an SOS2 constraint, which relies on the notion of contiguous variables. Since the variables in the SOS are ordered by weight, contiguity becomes ambiguous when multiple variables have the same weight.

It is often more efficient to capture SOS structure using linear constraints rather than SOS constraints. The optimizer will often perform this conversion automatically. This is controlled with two parameters: PreSOS1BigM and PreSOS2BigM. The conversion is done by adding constraints of the form <span>$</span>x <= M b<span>$</span>, where <span>$</span>x<span>$</span> is the variable that participates in the SOS constraint, <span>$</span>b<span>$</span> is a binary variable, and <span>$</span>M<span>$</span> is an upper bound on the value of variable <span>$</span>x<span>$</span>. Large values of <span>$</span>M<span>$</span> can lead to numerical issues, so these parameters control the maximum value of <span>$</span>M<span>$</span> that can be introduced by this conversion. SOS constraints that would require a larger value aren't converted.

Quadratic Constraints

A quadratic constraint allows you to restrict the value of a quadratic expression. For example, you may require that any feasible solution satisfy the constraint <span>$</span>3 x^2 + 4 y^2 + 5 z \leq 10<span>$</span>. Note that the matrix-oriented Gurobi API's (C, MATLAB, and R) require the right-hand side of a quadratic constraint to be a constant, while the object-oriented APIs (C++, Java, .NET, and Python) allow arbitrary quadratic expressions on both sides of the comparator.

The computed solution should satisfy the stated constraint to within FeasibilityTol. Quadratic constraints are often much more challenging to satisfy than linear constraints, so tightening the parameter may increase runtimes dramatically.

Gurobi can handle both convex and non-convex quadratic constraints. However, there are some subtle and important differences in how the different constraint types are handled. The default algorithms in Gurobi only accept a few forms of quadratic constraints that are known to have convex feasible regions. Constraints of the following forms are always accepted:

  • <span>$</span>x^TQx + q^Tx \le b<span>$</span>, where <span>$</span>Q<span>$</span> is Positive Semi-Definite (PSD)
  • <span>$</span>x^Tx \le y^{2}<span>$</span>, where <span>$</span>x<span>$</span> is a vector of variables, and <span>$</span>y<span>$</span> is a non-negative variable (a Second-Order Cone constraint)
  • <span>$</span>x^Tx \le y z<span>$</span>, where <span>$</span>x<span>$</span> is a vector of variables, and <span>$</span>y<span>$</span> and <span>$</span>z<span>$</span> are non-negative variables (a rotated Second-Order Cone constraint)
To be more precise, a constraint will be accepted if presolve is able to transform it into one of these forms. Note that if the quadratic terms each contain at least one binary variable, then presolve will always be able to transform it.

If you add a constraint that can't be transformed into one of these forms, then with default settings you will get an error (GRB_ERROR_Q_NOT_PSD) when you try to solve the model. Quadratic equality constraints are always non-convex; they will give a GRB_ERROR_QCP_EQUALITY_CONSTRAINT error with default settings.

Why distinguish between quadratic constraints in this form and other types of quadratic constraints? Solving models with non-convex quadratic constraints is typically much more expensive. To avoid accidentally solving a much harder problem than may have been intended, Gurobi rejects such constrains by default. If you set the NonConvex parameter to 2, however, then Gurobi will accept arbitrary quadratic constraints and attempt to solve the resulting model.

Note that other non-convex quadratic solvers often only find locally optimal solutions. The algorithms in Gurobi explore the entire search space, so they provide a globally valid lower bound on the optimal objective value, and given enough time they will find a globally optimal solution (subject to tolerances).

We would like to note a subtle point here regarding terminology. A quadratic constraint that involves only products of disjoint pairs of variables is often called a bilinear constraint, and a model that contains bilinear constraints is often called a bilinear program. Bilinear constraints are a special case of non-convex quadratic constraints, and the algorithms Gurobi uses to handle the latter are also well suited to solving bilinear programming problems.

General Constraints

The previously-described constraints are typically handled directly by the underlying optimization algorithms (but not always). Gurobi includes an additional set of constraints, which we collectively refer to as general constraints. General constraints are mostly a convenience feature, designed to allow you to define certain variable relationships easily without having to immerse yourself in the often esoteric details of how to model these relationships in terms of the more fundamental constraints of MIP. Capturing a single one of these general constraints can often require a large set of linear and SOS constraints, plus a number of auxiliary decision variables. By supporting them directly in the Gurobi API, we simplify the modeling process by performing the transformation to a corresponding MIP formulation automatically and transparently during the solution process.

What sorts of variable relationships can be captured with general constraints? We think of them as belonging to two types: function constraints and simple constraints. Function constraints allow you to state a relationship <span>$</span>y =
f(x)<span>$</span>, where <span>$</span>x<span>$</span> and <span>$</span>y<span>$</span> are Gurobi decision variables and <span>$</span>f()<span>$</span> is chosen from a predefined list of functions. Gurobi performs a piecewise-linear approximation of that function within the domain of <span>$</span>x<span>$</span>. Simple general constraints allow you to state common but more direct relationships between decision variables. The translation that goes on under the hood for these is much simpler, and the result is an exact representation of the original constraint (not an approximation).

Simple General Constraints

Gurobi supports the following simple general constraints, each with its own syntax and semantics:

  • MAX constraint: The constraint <span>$</span>r = \max\{x_1,\ldots,x_k,c\}<span>$</span> states that the resultant variable <span>$</span>r<span>$</span> should be equal to the maximum of the operand variables <span>$</span>x_1,\ldots,x_k<span>$</span> and the constant <span>$</span>c<span>$</span>. For example, a solution <span>$</span>(r=3, x_1=2, x_2=3, x_3=0)<span>$</span> would be feasible for the constraint <span>$</span>r = \max\{x_1,x_2,x_3,1.7\}<span>$</span> because <span>$</span>3<span>$</span> is indeed the maximum of <span>$</span>2<span>$</span>, <span>$</span>3<span>$</span>, <span>$</span>0<span>$</span>, and <span>$</span>1.7<span>$</span>.
  • MIN constraint: Similar to a MAX constraint, the constraint <span>$</span>r = \min\{x_1,\ldots,x_k,c\}<span>$</span> states that the resultant variable <span>$</span>r<span>$</span> should be equal to the minimum of the operand variables <span>$</span>x_1,\ldots,x_k<span>$</span> and the constant <span>$</span>c<span>$</span>.
  • ABS constraint: The constraint <span>$</span>r = \mbox{abs}\{x\}<span>$</span> states that the resultant variable <span>$</span>r<span>$</span> should be equal to the absolute value of the operand variable <span>$</span>x<span>$</span>. For example, a solution <span>$</span>(r=3, x=-3)<span>$</span> would be feasible for the constraint <span>$</span>r = \mbox{abs}\{x\}<span>$</span>.
  • AND constraint: The constraint <span>$</span>r = \mbox{and}\{x_1,\ldots,x_k\}<span>$</span> states that the binary resultant variable <span>$</span>r<span>$</span> should be equal <span>$</span>1<span>$</span> if and only if all of the binary operand variables <span>$</span>x_1,\ldots,x_k<span>$</span> are equal to <span>$</span>1<span>$</span>. For example, a solution <span>$</span>(r=1, x_1=1, x_2=1, x_3=1)<span>$</span> would be feasible for the constraint <span>$</span>r = \mbox{and}\{x_1,x_2,x_3\}<span>$</span>. Note that declaring an AND constraint implicitly declares all involved variables to be of binary type.
  • OR constraint: Similar to an AND constraint, the constraint <span>$</span>r = \mbox{or}\{x_1,\ldots,x_k\}<span>$</span> states that the binary resultant variable <span>$</span>r<span>$</span> should be <span>$</span>1<span>$</span> if and only if at least one of the binary operand variables <span>$</span>x_1,\ldots,x_k<span>$</span> is equal to <span>$</span>1<span>$</span>. Note that declaring an OR constraint implicitly declares all involved variables to be of binary type.
  • INDICATOR constraints: An indicator constraint <span>$</span>y = f \rightarrow a^Tx \leq b<span>$</span> states that if the binary indicator variable <span>$</span>y<span>$</span> has the value <span>$</span>f \in \{0,1\}<span>$</span> in a given solution, then the linear constraint <span>$</span>a^Tx \leq b<span>$</span> has to be satisfied. On the other hand, if <span>$</span>y \neq f<span>$</span> (i.e., <span>$</span>y = 1-f<span>$</span>) then the linear constraint may be violated. Note that the sense of the linear constraint can also be <span>$</span>=<span>$</span> or <span>$</span>\geq<span>$</span>; refer to this earlier section for a more detailed description of linear constraints. Note also that declaring an INDICATOR constraint implicitly declares the indicator variable to be of binary type.
  • Piecewise-linear constraints: A piecewise-linear constraint <span>$</span>y =
f(x)<span>$</span> states that the point <span>$</span>(x, y)<span>$</span> must lie on the piecewise-linear function <span>$</span>f()<span>$</span> defined by a set of points <span>$</span>(x_1, y_1), (x_2, y_2), ..., (x_n, y_n)<span>$</span>. Refer to the description of piecewise-linear objectives for details of how piecewise-linear functions are defined.
Note that adding any of these constraints to an otherwise continuous model will transform it into a MIP

As stated above, each general constraint has an equivalent MIP formulation that consists of linear and SOS constraints, and possibly auxiliary variables. Thus, you could always model such constraints yourself without using a Gurobi general constraint. For example, the MAX constraint <span>$</span>r = \max\{x_1,\ldots,x_k,c\}<span>$</span> can be modeled as follows:

\begin{displaymath}
\begin{array}{rcll}
r & = & x_j + s_j & \mbox{ for all } j =...
...& \in & \{0,1\} & \mbox{ for all } j = 1,\ldots,k+1
\end{array}\end{displaymath}

The first two constraints state that <span>$</span>r \geq \max\{x_1,\ldots,x_k,c\}<span>$</span>, i.e., that the resultant variable <span>$</span>r<span>$</span> has to be at least as large as each of the operand variables <span>$</span>x_j<span>$</span> and the constant <span>$</span>c<span>$</span>. This can be modeled using inequalities, but we turned them into equations by introducing explicit continuous slack variables <span>$</span>s_j \geq 0<span>$</span>, which we will reuse below.

Those slack variables and the remaining constraints model <span>$</span>r \leq \max\{x_1,\ldots,x_k,c\}<span>$</span>, which is more complicated. In addition to the explicit slacks, this requires the introduction of binary auxiliary variables <span>$</span>z_j \in \{0,1\}<span>$</span>. The SOS1 constraints state that at most one of the two variables <span>$</span>s_j<span>$</span> and <span>$</span>z_j<span>$</span> can be non-zero, which models the implication <span>$</span>z_j = 1 \rightarrow s_j = 0<span>$</span>. Due to the third constraint, one <span>$</span>z_j<span>$</span> will be equal to <span>$</span>1<span>$</span> and thus at least one <span>$</span>s_j<span>$</span> will be zero. Hence, <span>$</span>r = x_j<span>$</span> for at least one <span>$</span>j<span>$</span> due to the first constraint, or <span>$</span>r = c<span>$</span> due to the second constraint.

Tolerances play a role in general constraints, although as you might expect, the exact role depends on the constraint type. Generally, violations in the resultant will be smaller than the feasibility tolerance, and integrality violations in integer resultants will also satisfy the integrality tolerance.

By most measures, general constraints are just a means of concisely capturing relationships between variables while removing the burden of creating an equivalent MIP formulation. However, general constraints have another potential advantage: Gurobi might be able to simplify the MIP formulation if it can prove during presolve that the simplified version suffices for the correctness of the model. For this reason, Gurobi might be able to produce a smaller or tighter representation of the general constraint than you would get from the most general formulation. For example, it might be the case that <span>$</span>r \leq \max\{x_1,\ldots,x_k,c\}<span>$</span> is already implied by the other constraints in the model, so that a simple set of inequalities

\begin{displaymath}
\begin{array}{rcll}
r & \geq & x_j \;\;\mbox{ for all } j = 1,\ldots,k \\
r & \geq & c
\end{array}\end{displaymath}

to describe <span>$</span>r \geq \max\{x_1,\ldots,x_k,c\}<span>$</span> suffices to model the relevant part of the MAX constraint.

Function Constraints

Gurobi supports the following function constraints, each with somewhat different syntax and semantics (<span>$</span>x<span>$</span> and <span>$</span>y<span>$</span> below are Gurobi decision variables, and other terms are constants provided as input when the constraint is added to the model):

  • Polynomial: <span>$</span>y = p_0 x^n + p_1 x^{n-1} + ... + p_n x + p_{n+1}<span>$</span>
  • Natural exponential: <span>$</span>y = exp(x)<span>$</span> or <span>$</span>y = e^x<span>$</span>
  • Exponential: <span>$</span>y = a^x<span>$</span>, where <span>$</span>a > 0<span>$</span> is the base for the exponential function
  • Natural logarithm: <span>$</span>y = log_e(x)<span>$</span> or <span>$</span>y = ln(x)<span>$</span>
  • Logarithm: <span>$</span>y = log_a(x)<span>$</span>, where <span>$</span>a > 0<span>$</span> is the base for the logarithmic function
  • Power: <span>$</span>y = x^a<span>$</span>, where <span>$</span>a >= 0<span>$</span>
  • Sine: <span>$</span>y = sin(x)<span>$</span>
  • Cosine: <span>$</span>y = cos(x)<span>$</span>
  • Tangent: <span>$</span>y = tan(x)<span>$</span>

As noted earlier, Gurobi will automatically add a piecewise-linear approximation of the function to the model. You face a fundamental cost-versus-accuracy tradeoff when performing such an approximation, though: adding more pieces produces smaller approximation errors, but also increases the cost of solving the problem. The tradeoff can be complex. Gurobi provides a set of three attributes that help to navigate this tradeoff: FuncPieces, FuncPieceLength, FuncPieceError. They are used as follows:

  • If you would like to choose the number of pieces to use for the approximation, set the FuncPieces attribute to the desired value. All pieces will have equal width. This approach allows you to control the size of the approximation.
  • If you would like to choose the width of each piece, set the FuncPieces attribute to a special value of 1 and set the FuncPieceLength attribute equal to the desired width of each piece. This approach provides some control over both the size and the error of the approximation. While this may appear to be a minor variation of the first option, note that presolve may tighten the domain of <span>$</span>x<span>$</span>, often substantially, which can make it difficult to predict the relationship between the width of each piece and the number of pieces.
  • If you would like to set the maximum error you are willing to tolerate in the approximation, set the FuncPieces attribute to a special value of -1 and set the FuncPieceError attribute equal to the maximum absolute approximation you are willing to tolerate. Gurobi will choose pieces, typically of different sizes, to achieve that error bound. Note that the number of pieces required may be quite large if you set a tight error tolerance. You can control the maximum relative error rather than the absolute error by setting the FuncPieces attribute to -2 instead of -1.
These are attributes on the general constraints, so you can choose different values for each individual constraint.

The other relevant attribute is FuncPieceRatio, which controls whether the approximation is an underestimate of the function (0.0), an overestimate (1.0), or somewhere in between (any value strictly between 0.0 and 1.0). You can also choose the special value of -1, which will choose points that are on the original function.

Consider the following simple example:

\scalebox{1.0}{\includegraphics[width=2.5in]{graphics/func}}

The goal is to find an approximation of the polynomial <span>$</span>y=x^2<span>$</span>. We've set FuncPieces to <span>$</span>1<span>$</span> and FuncPieceLength to <span>$</span>1.0<span>$</span>, so we're performing an approximation with fixed-width pieces of width 1.0. The domain of <span>$</span>x<span>$</span> is <span>$</span>[0,2]<span>$</span>, so the approximation has two pieces. The figure shows 6 points: <span>$</span>P_{u1}(0,0), P_{u2}(1,1), P_{u3}(2,4)<span>$</span>, and <span>$</span>P_{l1}(0,-0.25), P_{l2}(1, 0.75), P_{l3}(2,3.75)<span>$</span>. If FuncPieceRatio is set to <span>$</span>0.0<span>$</span>, the approximation would be built from the points below the function ( <span>$</span>P_{l1}, P_{l2}<span>$</span>, and <span>$</span>P_{l3}<span>$</span>). Similarly, if it is set to <span>$</span>1.0<span>$</span>, the approximation would be built from the points above the function ( <span>$</span>P_{u1}, P_{u2}<span>$</span>, and <span>$</span>P_{u3}<span>$</span>). A value of <span>$</span>0.6<span>$</span> would use weighted combinations of the points: <span>$</span>0.6<span>$</span> times <span>$</span>P_{ui}<span>$</span> plus <span>$</span>0.4<span>$</span> times <span>$</span>P_{li}<span>$</span>. In this case, the line segments would be built from the points <span>$</span>(0, -0.1), (1,
0.9)<span>$</span>, and <span>$</span>(2, 3.9)<span>$</span>. If FuncPieceRatio is set to <span>$</span>-1<span>$</span>, meaning that the approximation would be built from points that are on the original function, in this case the upper points ( <span>$</span>P_{u1}, P_{u2}<span>$</span>, and <span>$</span>P_{u3}<span>$</span>) fit the bill. This will always be the case for a convex function.

Recall that you can set FuncPieces to <span>$</span>-1<span>$</span> to control the maximum absolute error. In this case, choosing a FuncPieceError value of <span>$</span>0.25<span>$</span> would give the piecewise approximation shown in the figure, since the distance between the upper and lower curves is always <span>$</span>0.25<span>$</span>. A smaller error value would of course lead to more pieces. We should add that piece widths will typically be non-uniform when limiting the maximum approximation error. The approximation algorithms we use try to limit the number of pieces needed to meet the error targets, which often requires more refinement in some portions of the domain than in others.

Note that the approximations are guaranteed to be under- and over-estimates in all cases except for polynomials of degree greater than 5. Finding the roots of higher-degree polynomials, which would be required to guarantee this property, is quite difficult.

If you wish to experiment with different approaches to approximating a set of functions, it is often convenient to be able to change the approach for all functions at once. We provide a set of parameters with the same names as the attributes to make it easier to do this: FuncPieces, FuncPieceLength, FuncPieceError, and FuncPieceRatio. If you set the FuncPieces attribute on a function constraint to <span>$</span>0<span>$</span>, then the approximation approach for that constraint will be determined by the parameter settings instead.

For some of the supported functions, modest <span>$</span>x<span>$</span> values can lead to enormous <span>$</span>y<span>$</span> values (and vice-versa). This can cause numerical issues when solving the resulting piecewise-linear MIP model. To avoid such issues, we limit the range of any <span>$</span>x<span>$</span> or <span>$</span>y<span>$</span> that participates in a function constraint to [-1e+6, 1e+6]. The parameter FuncMaxVal allows you to change these limits, but we recommend that you proceed with caution.

We should point out that PWL approximations can sometimes cause unexpected results, including sub-optimal solutions or even infeasible conclusions on feasible models. Consider a simple example with two constraints: <span>$</span>y = 2x-1<span>$</span> and <span>$</span>y=x^2<span>$</span>. Clearly <span>$</span>(x, y) = (1,
1)<span>$</span> is a feasible solution, but a piecewise-linear approximation could introduce breakpoints at <span>$</span>x=0.9<span>$</span> and <span>$</span>x=1.1<span>$</span>. The resulting approximation gives a <span>$</span>y<span>$</span> value of <span>$</span>1.01<span>$</span> at <span>$</span>x=1<span>$</span>, which is sufficiently far from the actual function value that Gurobi will not consider that a valid solution and declare the model infeasible, since there are no other solutions to the constraints. Reducing the maximum approximation error (by setting FuncPieces to -1 and FuncPieceError to a much smaller value) would help, but this isn't always the best way to address the problem, since tighter error tolerances can substantially increase the number of pieces in the approximation and thus the cost. We recommend the following approach when you encounter unexpected results. For inequalities, you should ask for an approximation that always overestimates or underestimates the function (depending on the sense of the constraint), to ensure that your approximation will always satisfy the constraint. The FuncPieceRatio parameter allows you to do this. For equalities, if you have a sense of where your solution is likely to lie, one option for managing the size of the approximation is to introduce additional variables to capture your function in different ranges, and then perform approximations with different levels of accuracy on these different pieces.

While users could perform piecewise-linear approximations themselves, there are several advantages to asking Gurobi to do it instead. First, Gurobi can often reduce the domains of variables, by using bound strengthening in presolve, or by exploiting repetition in periodic functions like sine or cosine. Smaller domains means fewer pieces to achieve the same accuracy. Gurobi also provides many options to make experimentation easier (for error control, piece length, etc.). These options can be quite difficult to implement and maintain.