I have implemented a set of "or" constraints in my ILP using binary decision variables (as in this method). It works fine for smaller problems, but when I try to increase the number of variables it gets very slow very fast, such that it is not feasible for the size of problem I need to solve. Is there a way to more cleverly implement this method, so that I can get a solution in hours instead of months?

Without the “or” constraint, the solution is found in a matter of seconds, even for a problem twice the size that I need, so it’s not just a matter of the number of variables (because the decision variables make the simplex no longer convex). I have assigned the big constant to be as small as possible while still satisfying the constraints, but I am not sure what else I can do. All my variables are binary, I am using a GLPK solver from R (Rglpk), running on a professional-grade laptop. I can successfully solve the problem for $1,000$ variables in a couple seconds, while $10,000$ takes a couple hours. My application calls for a maximum of $1,000,000$ variables.

Thanks for your advice.

EDITED to add information about the model:

My variables are the entries of a binary matching matrix $B$, weighted by a similarity matrix $S$, the constraints are just row sums and column sums:

$$\max_{B} \sum_{i,j} S_{i,j} B_{i,j}$$


  1. $B_{i,j}$ binary
  2. $\forall i, \sum_j B_{i,j} = a$
  3. $\forall j, \sum_i B_{i,j} \le b$
  4. $\forall j, \sum_i B_{i,j} \ge c$ or $\sum_i B_{i,j} \le 0$
  • $\begingroup$ Does this answer your question? Linear Programming with additional "if-then"/"Default to zero" constraints? $\endgroup$ – RobPratt Aug 5 at 15:55
  • $\begingroup$ Thanks for your comment RobPratt. That is indeed the same method I am using. My question is about improving the speed: is there anything I can do to be clever to solve the problem faster? $\endgroup$ – Animik Aug 5 at 16:01
  • 2
    $\begingroup$ You can reduce the number of constraint coefficients by introducing a column sum variable $C_j$ with bounds $[0,b]$ and constraint $C_j = \sum_i B_{i,j}$, and then omit your constraint 3 and use $C_j$ in your constraint 4. This could improve the underlying LP solves. $\endgroup$ – RobPratt Aug 5 at 17:41
  • $\begingroup$ Can you get an academic license for Gurobi? There are significant performance differences between open source and industrial solvers for ILP. $\endgroup$ – Richard Aug 6 at 15:39
  • $\begingroup$ Or just use Gurobi over NEOS... $\endgroup$ – kurtosis Aug 17 at 19:04

Answers to the linked question mention both big-M constraints and semicontinuous variables. To speed up the big-M approach, you might consider introducing the constraints dynamically only as they are violated ("row generation" or "cut generation"). Explicitly:

  1. Omit all big-M constraints and the associated binary variables.
  2. Solve the problem with the current set of big-M constraints and associated binary variables.
  3. For any violations (variables $x_i$ with $0 < x_i < c$), introduce big-M constraints and the associated binary variables.
  4. If there are no violations, stop. Else go to step 2.

Depending on the rest of your model, there might be additional ways to improve the solve time.

| improve this answer | |
  • $\begingroup$ I edited my question to add information about the model. If I solve without the last constraint, usually more than half the columns violate that constraint, so that’s still too big of a problem to run all together. Do you suggest to add each constraint one at a time? $\endgroup$ – Animik Aug 5 at 16:48
  • $\begingroup$ I had assumed that any implementation of semi-continuous variables would be equally slow to what I am doing, but maybe that is an incorrect assumption! $\endgroup$ – Animik Aug 5 at 16:51
  • $\begingroup$ If the solver just converts to big-M constraints, the performance should be the same, but I believe that some solvers omit the big-M constraints and instead explicitly branch on the semicontinuous variables. That is an approach you might consider coding up yourself. $\endgroup$ – RobPratt Aug 5 at 17:35
  • $\begingroup$ I'm not sure if it's possible to feed an initialization into any of the open source ILP solvers in R, which is too bad, but I suppose that would be the topic for a different question. I am wondering if there might be something clever to do with the dual that could help decrease solve time? $\endgroup$ – Animik Aug 5 at 20:33

GLPK is not the best performing MILP solver. Instead, you could give one of the leading commercial MILP solvers a try (e.g. Gurobi). You can also try open-source solvers like SCIP a try. Those solvers should be faster out of the box.

You can quickly evaluate different solvers with your model by writing it out as .MPS file. Every MILP solver I know of can read these files.

Some solvers, like Gurobi, have special APIs to formulate OR constraints and indicator constraints. By using these APIs you can sidestep finding a suitable Big-M yourself and let the solver find appropriate values (which may even tightened during the solving process).

Still, a million variables is probably challenging for these solvers as well.

| improve this answer | |
  • $\begingroup$ Thanks for your suggestions Simon, unfortunately for this project it is important to remain open source with an R wrapper. I had understood that rGLPK was the fastest of the open source solvers available wrapped at a higher level in R (clpAPI doesn't support IP, and lpSolve is slower than rGLPK, etc), is this no longer correct? $\endgroup$ – Animik Aug 13 at 15:07
  • $\begingroup$ @Animik: Sorry, didn't get that. Have you tried github.com/dirkschumacher/rcbc ? Also, are you aware that Gurobi is giving away free licenses for academic research? $\endgroup$ – Simon Aug 14 at 8:54

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.