In this post I present a Groovy framework for Linear Programming (hereafter LP) reductions. The term 'framework' denotes that its expressive power covers the reduction of any problem P into LP.
Given a Problem P with a domain D and range (or solution domain if you wish) R, we say that a problem P reduces to LP, if you can use an algorithm that solves LP to solve P. The reduction includes 3 steps (a nice visual representation can be found at page 9/46 of this):
1. given an instance x of P, map it to an instance y of LP;
2. solve y by an algorithm (in this case Simplex) for LP, to derive a solution, let's denote it with LP(y);
3. map LP(y) to the solution of x.
The mapping functions that are needed by steps (1) and (3) are passed as inputs in the form of Groovy closures, while step (2) is delegated to Apache Math3 Simplex implementation by default, although abstractions for replacing this with the solver of your wish are provided.
Of course, not any possible combination of map functions in step (1) and (3) lead to a reduction. You need to prove that a certain combination is proper. The framework presented here does not verify the combination is proper (there are certain limitations even if we wished to do so), but assumes it is proper.
Specifically for B=LP:
1. LP is polynomial. The ellipsoid algorithm solves this polynomially, although the simplex algorithm (which is exponential in the worst case) is much faster in practice.
2. The expressive power of LP and its massive application in the financial sector (e.g. see operations research) has led to very efficient, commercial simplex implementations offered by many business analyst vendors, like IBM.
Given a Problem P with a domain D and range (or solution domain if you wish) R, we say that a problem P reduces to LP, if you can use an algorithm that solves LP to solve P. The reduction includes 3 steps (a nice visual representation can be found at page 9/46 of this):
1. given an instance x of P, map it to an instance y of LP;
2. solve y by an algorithm (in this case Simplex) for LP, to derive a solution, let's denote it with LP(y);
3. map LP(y) to the solution of x.
The mapping functions that are needed by steps (1) and (3) are passed as inputs in the form of Groovy closures, while step (2) is delegated to Apache Math3 Simplex implementation by default, although abstractions for replacing this with the solver of your wish are provided.
Of course, not any possible combination of map functions in step (1) and (3) lead to a reduction. You need to prove that a certain combination is proper. The framework presented here does not verify the combination is proper (there are certain limitations even if we wished to do so), but assumes it is proper.
Background
Generally, algorithmic reductions are very useful for both theoretical and practical reasons. For instance, if A reduces to B polynomially (i.e. the mapping functions are polynomial) and B is a polynomial problem, then A is polynomial as well (for more details about the significance of reductions you can read this or this). In addition, if you have an efficient implementation for an algorithm B, but not for A, you can leverage this implementation.Specifically for B=LP:
1. LP is polynomial. The ellipsoid algorithm solves this polynomially, although the simplex algorithm (which is exponential in the worst case) is much faster in practice.
2. The expressive power of LP and its massive application in the financial sector (e.g. see operations research) has led to very efficient, commercial simplex implementations offered by many business analyst vendors, like IBM.
Source code, Dependencies, License
The lp-reduction framework is an open source project licensed under Apache License, Version 2.0.
The code is hosted by github under the lp-reduction repository.
Gradle is used for building and dependency management.
The only compile time dependency is the Apache Commons Mathematics Library, while testing is done with Spock.How it works
The implementation is straightforward as the main class (LPReduction) suggests:
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
class LPReduction<D, R> { | |
Closure vars | |
Closure objFunc | |
Closure constraints | |
Closure result | |
// default is maximisation problem | |
boolean max = true | |
// default is to allow negative solution values | |
boolean nonNegative = false | |
LPSolver solver | |
// caches the variables to avoid calling 'vars' closure more than once, | |
// in case 'constraints' closure needs this info. | |
def variables | |
/** | |
* Maps an instance x of P into an LP instance y | |
* | |
* @param x the domain object | |
* @return the reduced LP instance y | |
*/ | |
LPInstance mapToLP(D x) { | |
variables = vars(x) | |
new LPInstance(variables, | |
objFunc ? objFunc(x) : new LPFunction([:], 0), | |
constraints(x), | |
max, | |
nonNegative) | |
} | |
/** | |
* Solves the specified LP instance. | |
* If no solver is set, it uses ApacheSimplexSolver by default. | |
* | |
* @param lpi the LP instance to solve | |
* @return the solution of the LP instance | |
*/ | |
private Solution solve(LPInstance lpi) { | |
if (!solver) { | |
solver = new ApacheSimplexSolver() | |
} | |
solver.solve(lpi) | |
} | |
/** | |
* Computes P(x) by reducing P to LP | |
* | |
* @param x the domain object | |
* @return P(x) | |
*/ | |
R reduce(D x) { | |
LPInstance y = mapToLP(x) | |
Solution solution = solve(y) | |
result(solution) | |
} | |
} |
Step (1) is implemented by 'mapToLP' method. It takes as input an instance x of P and returns the LP instance y. To do so it calls the 'vars'/'objFunc'/'constraints' closures that are defined as class fields in order to create the LP variables/objective function/constraints respectively. How they are initialized will be clarified later with an example.
Notice that the LPInstance class has two more boolean fields, namely 'max' and 'nonNegative'. The former indicates whether this is a maximization (default option) or a minimization (set max=false for this) problem.
Moreover, it is possible that an LP instance has empty (technically any constant, usually 0) objective function and in this case it is often called constraint satisfaction problem. This happens when you only need to know if there is a solution that simultaneously satisfies all the constraints, without having any quantity to maximize/minimize. This case is handled by the line:
objFunc ? objFunc(x) : new LPFunction([:], 0)
which means "if the closure objFunc is set (i.e. not null), call it, else use the empty objective function (i.e. f(x)=0)".
Step (2) is implemented by 'solve' method, which takes as input an LP instance y and returns its solution. If no solver is set (default option), the Apache Math3 Simplex is used to solve y, else the user provided solver is used. To allow this LPSolver interface provides the proper abstraction:
The solution is modeled by Solution class:
The whole reduction is implemented by 'reduce' method. Obviously, it runs the first two steps as explained above and then it calls the 'result' closure to map the solution of the LP instance y to the solution of P for x.
Notice that the LPInstance class has two more boolean fields, namely 'max' and 'nonNegative'. The former indicates whether this is a maximization (default option) or a minimization (set max=false for this) problem.
Moreover, it is possible that an LP instance has empty (technically any constant, usually 0) objective function and in this case it is often called constraint satisfaction problem. This happens when you only need to know if there is a solution that simultaneously satisfies all the constraints, without having any quantity to maximize/minimize. This case is handled by the line:
objFunc ? objFunc(x) : new LPFunction([:], 0)
which means "if the closure objFunc is set (i.e. not null), call it, else use the empty objective function (i.e. f(x)=0)".
Step (2) is implemented by 'solve' method, which takes as input an LP instance y and returns its solution. If no solver is set (default option), the Apache Math3 Simplex is used to solve y, else the user provided solver is used. To allow this LPSolver interface provides the proper abstraction:
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
@CompileStatic | |
interface LPSolver { | |
/** | |
* Solves the specified LP instance | |
* | |
* @param lpi the LP instance to solve | |
* @return the solution of the LP instance | |
*/ | |
Solution solve(LPInstance lpi) | |
} |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
@TupleConstructor | |
@CompileStatic | |
class Solution<V> { | |
Double optimum | |
Map<V, Double> at | |
LPIState state | |
} | |
enum LPIState { | |
FEASIBLE, | |
NOT_FEASIBLE, | |
UNBOUNDED | |
} |
The whole reduction is implemented by 'reduce' method. Obviously, it runs the first two steps as explained above and then it calls the 'result' closure to map the solution of the LP instance y to the solution of P for x.
The Max Flow Reduction Example
The code below exercises the lp-reduction framework for reducing the Max Flow problem to LP. You can read here what the Max Flow problem is and how an instance of it can be mapped to an LP instance.
As shown above, you can initialize an LPReduction object fields with the 'with' Groovy keyword. After defining the LPReduction object, you can exercise it as shown in MaxFlowSpec.groovy. Notice that the code inside the closures vars/objFunc/constraints/result can be written in Java instead of Groovy. However, you can not replace the entire closure with a Java 8 lambda as Groovy is not compatible with lamdas.
The full code can be found github under the lp-reductions repository.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
static LPReduction reduction() { | |
def reduction = new LPReduction() | |
reduction.with { | |
vars = { Graph g -> | |
def vars = g.edges.collect { it.key } | |
vars << [g.t, g.s] | |
} | |
objFunc = { Graph g -> | |
def coeffMap = g.outEdges[g.s].collectEntries { [(it): 1] } | |
new LPFunction(coeffMap, 0) | |
} | |
constraints = { Graph g -> | |
def constraints = g.edges.collect { e, l -> | |
LPFunction f = new LPFunction([(e):1], 0) | |
new LPConstraint(f, ConstraintRelationship.LEQ, l) | |
} | |
constraints += g.nodes.findAll { | |
it != g.s && it != g.t | |
}.collect { n -> | |
def coeffMap = [:] | |
g.inEdges[n].each { e -> coeffMap[e] = 1 } | |
g.outEdges[n].each { e -> coeffMap[e] = -1 } | |
LPFunction f = new LPFunction(coeffMap, 0) | |
new LPConstraint(f, ConstraintRelationship.EQ, 0) | |
} | |
} | |
result = { Solution s -> s.state == LPIState.FEASIBLE ? new Result(s.optimum, s.at) : null } | |
nonNegative = true | |
} | |
reduction | |
} |
The full code can be found github under the lp-reductions repository.
The Graph Realization Reduction Example
The following code exercises the lp-reduction framework for reducing the Graph Realization problem to LP.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
static LPReduction reduction(boolean decisionVersion) { | |
def reduction = new LPReduction() | |
reduction.with { | |
vars = { List<Integer> seq -> | |
[0..seq.size()-1, 0..seq.size()-1].combinations() | |
.findAll( { u, v -> u < v } ) | |
.collect( { u, v -> [u, v]} ) | |
} | |
constraints = { List<Integer> seq -> | |
def constraints = (0..seq.size()-1).collect { u -> | |
def coeffMap = variables.findAll{ e -> e[0]==u || e[1]==u }.collectEntries { e -> [(e) : 1] } | |
LPFunction f = new LPFunction(coeffMap, 0) | |
new LPConstraint(f, ConstraintRelationship.EQ, seq[u]) | |
} | |
constraints += variables.collect { e -> | |
LPFunction f = new LPFunction([(e):1], 0) | |
new LPConstraint(f, ConstraintRelationship.LEQ, 1) | |
} | |
} | |
if (decisionVersion) | |
result = { Solution s-> s.state==LPIState.FEASIBLE } | |
else | |
result = { Solution s -> s.state==LPIState.FEASIBLE ? s.at : null } | |
} | |
reduction | |
} |
Notice that in its original version the Graph Realization is a decision problem, i.e. is there any graph with the given sequence, YES or NO? In practice, it is often useful to ask for a graph with the given sequence if such a graph exists. To define the decision version you need to define the result closure as:
result = { Solution s-> s.state==LPIState.FEASIBLE }
while for the modified version you need to define it as
result = { Solution s -> s.state==LPIState.FEASIBLE ? s.at : null }
After defining the LPReduction object, you can exercise it as shown in GraphRealizationSpec.groovy.
result = { Solution s-> s.state==LPIState.FEASIBLE }
while for the modified version you need to define it as
result = { Solution s -> s.state==LPIState.FEASIBLE ? s.at : null }
After defining the LPReduction object, you can exercise it as shown in GraphRealizationSpec.groovy.