*irreducible infeasible subsystems*(IISes) and feasOpt() to find a “cheapest” relaxation of constraints and bounds that makes the model feasible. They also provide access to what is known as a

*Farkas certificate*.

Thanks to Bo Jensen, who reminded me about Farkas certificates in a comment yesterday. He also provided a link to MOSEK's description of them, which saves me quite a bit of typing. Farkas certificates draw their name from Farkas's lemma, a well-known result in linear programming theory. While it is possible for both a primal linear program and its dual to be infeasible, that condition is rare in practice. The (much) more common case is that an infeasible LP has an unbounded dual. In essence, a Farkas certificate is a recession direction for the dual problem. It's existence proves the infeasibility of the primal.

For infeasible LPs, CPLEX provides the IloCplex.dualFarkas() method to get a Farkas certificate. (As always, I'll stick to the Java API, but all the APIs should provide some version of this.) The method has two arguments used for returning results (and

*not*for providing inputs). One argument returns a list of primal constraints with non-zero dual values; the other vector holds the corresponding dual values.

Continuing yesterday's Java code (without repeating it), I'll start by constructing the new primal problem:

```
IloCplex primal = new IloCplex();
IloNumVar u = primal.numVar(0, Double.MAX_VALUE, "u");
IloNumVar v = primal.numVar(0, Double.MAX_VALUE, "v");
IloNumVar w = primal.numVar(0, Double.MAX_VALUE, "w");
IloNumVar[] vars = new IloNumVar[] {u, v, w};
primal.addMinimize(
primal.scalProd(new double[] {1.0, 1.0, -2.0}, vars));
IloConstraint c1 =
primal.addGe(
primal.scalProd(new double[] {1.0, -2.0, -1.0}, vars),
3.0);
IloConstraint c2 =
primal.addGe(
primal.scalProd(new double[] {-2.0, 1.0, -1.0}, vars),
2.0);
```

Since I cannot be sure which constraints dualFarkas will return, nor in what order, it will be handy to construct a map associating each constraint with the corresponding dual variable.

```
HashMap<IloConstraint, IloNumVar> map =
new HashMap<IloConstraint, IloNumVar>();
map.put(c1, x);
map.put(c2, y);
```

Now I turn off the presolver, solve the model, and verify that it is infeasible. (As was the case yesterday with an unbounded primal, fail to turn off the presolver will cause CPLEX to throw an exception when I try to access the Farkas certificate.) I also need to use the dual simplex algorithm; any other algorithm (primal simplex, barrier, ...) will cause CPLEX to throw an exception. There is no explicit statement in the code since CPLEX's default choice for algorithm

*on this problem*is dual simplex, but in general it should probably be specified just in case CPLEX prefers some other algorithm.

```
primal.setParam(IloCplex.BooleanParam.PreInd, false);
primal.solve();
System.out.println("CPLEX status: "
+ primal.getCplexStatus());
```

The next step is to access the Farkas certificate:

```
IloRange[] farkasConstraints = new IloRange[2];
double[] farkasValues = new double[2];
primal.dualFarkas(farkasConstraints, farkasValues);
System.out.print("Farkas certificate:");
for (int i = 0; i < 2; i++) {
System.out.print("\t" + map.get(farkasConstraints[i])
+ " -> " + farkasValues[i]);
}
System.out.println("");
```

To find the vertex at which it is anchored, I can access the last feasible dual solution obtained before CPLEX stopped.

```
double[] q = primal.getDuals(farkasConstraints);
System.out.print("Dual vertex:");
for (int i = 0; i < 2; i++) {
System.out.print("\t" + map.get(farkasConstraints[i])
+ " -> " + q[i]);
}
System.out.println("");
```

I still have one concern. The definition of a Farkas certificate only guarantees a recession direction for the dual; it does not need to be an extreme dual direction to qualify as a Farkas certificate. I suspect, but do not know with certainty, that the mechanism by which solvers select a Farkas certificate will consistently lead to an extreme direction.

The properties of the farkas ray returned, depends on the solver used. If you use an interior point optimizer (which allows you access to the interior point solution before basis identification), then the ray will not be an extreme ray if multiple solution exists. If the optimizer used is simplex, then the ray is an extreme ray. All the information needed is right there in the simplex method when it is either infeasible or unbounded, you just have to get it back through the postsolve. I think it only gets more confusing for the user if they have to switch off presolve and use certain algorithms etc. The optimizer should return farkas rays pr. default and report the quality of them, then there can be no discussion if the problem is infeasible or not, it is a proof. Farkas rays can be used in many cases such as DW, where you practically can implement phase one as phase two.

ReplyDeleteIf you're using primal simplex on an unbounded problem or dual simplex on an infeasible problem, the Farkas certificate found should be an extreme direction. What if you're using dual simplex on an unbounded problem or primal simplex (phase I) on an infeasible one? Do you still have an obvious Farkas certificate in the last iteration, or do you need to do additional computation? If it takes more work to get the certificate in those cases, then someone solving a large problem who only cares that the model is unbounded/infeasible (not looking for a proof) might not want to spend the extra CPU cycles.

ReplyDeleteThat said, I agree with you that forcing the user to jump through hoops to get the certificate is not a good idea. It seems particularly problematic to me that in CPLEX the user has to turn off presolving. Assuming the user is not anticipating an infeasible/unbounded result, the user will want presolving on for the gain in computational efficiency; so CPLEX is forcing her/him to solve the problem a second time to get the certificate (which, again, is a PITA if the model is large).

Do you have a blog? An entry about the use of Farkas certificates in DW would be make for interesting reading.

If you are using dual simplex and the problem is dual infeasible, then you still get a primal farkas ray for free, so thats no excuse for being sloppy :-)

ReplyDeleteI have a empty blog, I never found time to write and I think trafic is too limited. Actually I think OR blogs should be gathered in one place such as or-exchange, which would make it a lot more active. Blogs is meant to have a lot of comments and share ideas, but most blogs even Dr. Trick's have few or none comments, which is a shame.

The information you have written about infeasibility/unboundedness is very useful for many people, the right topics for blogs is everything :-) I love blogs with a practical view.

I should add, if you use some other phase one strategy than minimizing sum of infeasibilities (such as penalty heuristics), then you don't have a certificate, but on the other hand you can't either for sure declare the problem infeasible. Most commercial codes use some variant of minimizing sum of infeasibilities.

ReplyDeleteI have answered in blog

ReplyDeletehttp://erlingdandersen.blogspot.com/2010/08/exploiting-farkas-certificate-in.html

how to exploit the Farkas certificate in DW decomposition.

@Erling: Excellent! I think the post is very informative, and I hope that people working with D-W also find it. It does seem to be one of those useful things that somehow slips through the cracks.

ReplyDeleteHi Paul.Thanks for the post.Please do you know if I can get the extreme ray/certificate of infeasibility information from CPLEX API for MATLAB?If yes, what code do I use?I have been looking for that information for sometime now. Thanks.

ReplyDeleteEmmanuel

Emmanuel: According to this FAQ, it is not possible from MATLAB.

DeleteDear Prof. Rubin,

ReplyDeleteA friend of mine used dualFarkas method without any concern on the ordering of returned ray values. (he assumes the returned constraints are in order) and he finds the correct result. Is this a coincidence? Should I be careful about mapping? Thanks in advance.

I'm not sure. The corresponding call in the C API does *not* pass a vector of constraints, which implies that the returned vector of dual variables corresponds to the constraints in the order they appear in the model. So that would suggest your friend's assumption was sound. On the other hand, the Java API *does* pass a vector of constraints, which suggests that your friend just got lucky (else why bother with that argument?). The C++ API also has that argument.

DeleteThe higher level APIs (C++, Java, ...) don't always behave the same way that the C API does, so my guess is that CPLEX does not guarantee a particular order for the certificate obtained via dualFarkas(). It's just a guess, but I would be careful about the ordering if I were doing it.

Dear Prof. Rubin,

ReplyDeletethank you for your post but I have one question, what if I'm getting a feasible solution returned by the solver for the master problem and I need to get the dual prices of the constraints in order to pass them to the pricing problem !? the dualFarkas() method requires an infeasible solution .. what could be the way to get the dual prices in this case ?

Thank you very much in advance,

Invoke the getDual or getDuals method on the master problem. You can find a call to getDuals in my sample code above.

Delete