Constraint Relaxation
OngoingHi all,
I am currently building the MILP model and would like to know how to make relaxations (RHS of constraints) if the model becomes infeasible.
For example, let's say the model has two different constraints : Con1 and Con2.
If the model turns out infeasible then I want Gurobi to change the RHS of Con1 for feasibility.
However RHS of Con2 can't be relaxed (i.e. the relaxation can't be applied to Con2).
Would you please let me know how to setup relaxation condition?
Thank you for your help in advance.
Sung
-
Official comment
Hi Sung -
Please have a look at our Knowledge Base article "How do I determine why my model is infeasible?". I think you'll find that Model.feasRelax() will help with what you would like to do.
Gwyneth
-
Thank you Gwyneth,
Yes, feasRelax ( relaxobjtype, minrelax, vars, lbpen, ubpen, constrs, rhspen ) is the function I need to use.
According to the documentation, constrs argument tells Gurobi which constraints could be relxaed.
Could you please tell me how I can use constrs more in detail?
For example, if I have 10 constraints named as con1, con2,,,,,, con9, and con10.
(Constraints con1, con2, and con3 could be relaxed.)
Would constrs argument be constrs= [con1, con2, con3] ?
Thank you for your help in advance.
Sung
0 -
Yes, that's correct. For example, if you want to construct a feasibility relaxation to minimize the sum of the absolute values of constraint violations for only those three constraints (\( \texttt{relaxobjtype=0} \), \( \texttt{rhspen=[1, 1, 1]} \)) , you could use the following:
torelax = [con1, con2, con3]
conpens = [1]*len(torelax)
model.feasRelax(relaxobjtype=0, minrelax=False, vars=None, lbpen=None, ubpen=None, constrs=torelax, rhspen=conpens)
model.optimize()0 -
Thank you Eli for your comment!
0 -
Hi Eli,
I have one more question for you. How do I get con1, con2, and con3 for the argument constrs?
For example, if the model constraints look as below, then how can I define torelax list?
for j in Job:
model.addConstr( quicksum(schedule[i] >= k[j] for i in Index if i.student_id > 1000), "con1" )
model.addConstr( quicksum(schedule[i] <= p[j] for i in Index if i.offering_id > 0), "con2" )
Do I have to use getConstrs() or getConstrByName() to get the constraint information for the constrs?
(getConstrByName only returns randomly selected constraint if multiple linear constraints have the same name, which is my case, so I guess this is not ).
Thank you for your help.
- Sung
0 -
The constraint sense (e.g., \( \geq \)) should never show up inside of the quicksum() function. Also, like you noticed, you use the same constraint name for each individual \( \texttt{con1} \) and \( \texttt{con2} \) constraint. Constraint and variable names should be unique, otherwise you can encounter problems when retrieving constraints/variables by name or writing/reading a model file of the problem. When adding constraints inside of a \( \texttt{for} \) loop, it's good practice to add some indexing information to ensure constraint names are unique. In this case, you could use \( \texttt{f"con1_{j}"} \) instead of \( \texttt{con1} \).
To collect the \( \texttt{con1} \) and \( \texttt{con2} \) constraints into an argument for Model.feasRelax(), I would suggest using Model.addConstrs(). This method returns a tupledict that maps the constraint indices to their respective Constr objects. The method also adds indices to the end of the constraint names so the names are unique. Afterwards, you can collect the Constr objects into a list to use as the \( \texttt{constrs} \) keyword argument of Model.feasRelax(). For example:
con1 = model.addConstrs((quicksum(schedule[i] for i in Index if i.student_id > 1000) >= k[j] for j in Job), "con1")
con2 = model.addConstrs((quicksum(schedule[i] for i in Index if i.offering_id > 0) <= p[j] for j in Job), "con2")
# Create list of con1/con2 Constr objects
torelax = con1.values() + con2.values()
conpens = [1]*len(torelax)
model.feasRelax(relaxobjtype=0, minrelax=False, vars=None, lbpen=None, ubpen=None, constrs=torelax, rhspen=conpens)
model.optimize()If you want to add constraints in a \( \texttt{for} \) loop using Model.addConstr(), you need to create a list beforehand and iteratively add each Constr object to that list:
con1, con2 = [], []
for j in Job:
con1.append(model.addConstr(quicksum(schedule[i] for i in Index if i.student_id > 1000) >= k[j], f"con1_{j}"))
con2.append(model.addConstr(quicksum(schedule[i] for i in Index if i.offering_id > 0) <= p[j], f"con2_{j}"))
# Combine con1 and con2 Constr object lists
torelax = con1 + con2
conpens = [1]*len(torelax)
model.feasRelax(relaxobjtype=0, minrelax=False, vars=None, lbpen=None, ubpen=None, constrs=torelax, rhspen=conpens)
model.optimize()The \( \texttt{rhspen} \) argument should be a list of numeric values of equal length to the list provided as the \( \texttt{constrs} \) argument. Each element represents the penalty incurred for violating the corresponding constraint in the \( \texttt{constrs} \) list. If all constraints should incur the same penalty when violated, this can simply be a list of ones, as done in the code snippets above.
The feasibility relaxation adds nonnegative variables to the left-hand sides of the relaxed constraints. This process depends on the sense and name of each constraint. Consider a constraint named \( \texttt{mycon} \) that will be relaxed in the feasibility relaxation.
- If \( \texttt{mycon} \) is a \( \leq \) constraint, Gurobi introduces a nonnegative variable named \( \texttt{ArtN_mycon} \) to the model, then adds the term \( \texttt{-ArtN_mycon} \) to the constraint's left-hand side.
- If \( \texttt{mycon} \) is a \( \geq \) constraint, Gurobi introduces a nonnegative variable named \( \texttt{ArtP_mycon} \) to the model, then adds the term \( \texttt{ArtP_mycon} \) to the constraint's left-hand side.
- If \( \texttt{mycon} \) is an \( = \) constraint, Gurobi introduces nonnegative variables named \( \texttt{ArtN_mycon} \) and \( \texttt{ArtP_mycon} \) to the model, then adds \( \texttt{ArtP_mycon - ArtN_mycon} \) to the constraint's left-hand side.
After optimizing the feasibility relaxation, these new variables take on nonzero values if the corresponding constraint is relaxed in the solution. So, to determine the amount by which each constraint was relaxed, check the solution values of variables starting with \( \texttt{ArtP_} \) or \( \texttt{ArtN_} \). The constraint name can be extracted from the variable name.
for v in model.getVars():
if (v.VarName.startswith('ArtN_') or v.VarName.startswith('ArtP_')) and v.X > 1e-8:
print(f"Constraint {v.VarName[5:]} relaxed by {v.X}")0 -
Hi Eli,
Instead of passing con1+con2, is there a way to set constraint priorities (like Mandatory, High, Medium & Low) so that the relaxer automatically chooses constraints to relax based on the priority? Any inputs are appreciated!
0 -
This can be accomplished with Model.feasRelax() by (i) assigning different penalty values for constraint violations via the \( \texttt{constrs} \) and \( \texttt{rhspen} \) keyword arguments and (ii) defining how these penalties should be applied via the \( \texttt{relaxobjtype} \) keyword argument.
Assigning penalty values to constraints
In the code above, the list passed as the \( \texttt{constrs} \) keyword argument contains the constraints that are allowed to be violated. The \( \texttt{rhspen} \) keyword argument is a list of equal length representing the associated penalties. In the example, the \( \texttt{rhspen} \) list was all \( \texttt{1} \)'s, essentially meaning each constraint was of equal priority.
For the priorities you describe, you could assign \( \texttt{rhspen} \) values as follows:
- "Mandatory" constraints are not added to the list of constraints that are allowed to be violated.
- "High"-priority constraints can be violated but are assigned a very high \( \texttt{rhspen} \) value like \(100\).
- "Medium"-priority constraints can be violated and are assigned a medium \( \texttt{rhspen} \) value like \(10\).
- "Low"-priority constraints can be violated and are assigned a very low \( \texttt{rhspen} \) value like \(1\).
The exact numbers are completely dependent on the application and how you quantify these priorities.
Choosing a penalty metric
You also need to decide how these \( \texttt{rhspen} \) values should be applied when a constraint is violated. As an example, consider the constraint \( x + y = 5 \). Let's say we allow this constraint to be relaxed and assign it a \( \texttt{rhspen} \) value of \( 10 \). Now, consider a solution to the feasibility relaxation model with \( (x, y) = (4, 4) \). This solution violates the right-hand side of our constraint by \( 3 \). The exact penalty this violation induces in the objective function of the feasibility relaxation model is governed by the \( \texttt{relaxobjtype} \) keyword argument. Specifically:
- If you set \( \texttt{relaxoobjtype} \) to \( 0 \), the contribution of this constraint violation to the objective function of the feasibility relaxation model is the weighted magnitude of the constraint violation: \( 10 \cdot 3 = 30 \).
- If you set \( \texttt{relaxoobjtype} \) to \( 1 \), the objective function contribution is the squared weighted constraint violation: \( 10 \cdot 3^2 = 90 \).
- If you set \( \texttt{relaxobjtype} \) to \( 2 \), the objective function contribution is simply the \( \texttt{rhspen} \) value we assigned to this constraint, regardless of the amount by which the constraint is violated: \( 10 \).
For more information, see the documentation for Model.feasRelax().
0 -
Hi there!
(How to visualise what constraints are violated with feasRelax)
There is lots of helpful comments on the different ways to use feasRelax() and feasRelaxS().
When I call it, it does indeed find an optimal solution. However, I am really struggling with getting it to identify which of the constraints it has relaxed or broken?
Apologies if I am meant to be on a different thread
0 -
The article How do I change variable and/or constraint bounds to make an infeasible model feasible using feasRelax? gives an example of how to identify which variable constraints/bounds were relaxed. In particular, Gurobi's feasibility relaxation adds nonnegative variables to the model that correspond to "slack" on each constraint. So, you can:
- Before building the feasibility relaxation, query the number of variables in your model with the NumVars model attribute.
- After solving the feasibility relaxation, print the names of newly added variables with nonzero values. The names of these variables correspond to the names of constraints and/or variable bounds that were relaxed.
og_numvars = m.NumVars
# Relax variable bound and constraints
m.feasRelaxS(0, False, True, True)
m.optimize()
if m.Status == GRB.OPTIMAL:
# print the values of the artificial variables of the relaxation
print("\nSlack values:")
slacks = m.getVars()[og_numvars:]
for sv in slacks:
if sv.X > 1e-9:
print(f"{sv.VarName} = {sv.X}")The output will look something like this:
Slack values: ArtU_varname1 = 1 ArtL_varname2 = 2.79911e-06 ArtP_constrname1 = 1.57294e-04 ArtN_constrname2 = 2
0
Please sign in to leave a comment.
Comments
10 comments