diff options
Diffstat (limited to 'Solvers/VIATRA-Solver/hu.bme.mit.inf.dslreasoner.viatrasolver.logic2viatra/src/hu/bme/mit/inf/dslreasoner/viatrasolver/logic2viatra/cardinality/CbcPolyhedronSolver.xtend')
-rw-r--r-- | Solvers/VIATRA-Solver/hu.bme.mit.inf.dslreasoner.viatrasolver.logic2viatra/src/hu/bme/mit/inf/dslreasoner/viatrasolver/logic2viatra/cardinality/CbcPolyhedronSolver.xtend | 241 |
1 files changed, 241 insertions, 0 deletions
diff --git a/Solvers/VIATRA-Solver/hu.bme.mit.inf.dslreasoner.viatrasolver.logic2viatra/src/hu/bme/mit/inf/dslreasoner/viatrasolver/logic2viatra/cardinality/CbcPolyhedronSolver.xtend b/Solvers/VIATRA-Solver/hu.bme.mit.inf.dslreasoner.viatrasolver.logic2viatra/src/hu/bme/mit/inf/dslreasoner/viatrasolver/logic2viatra/cardinality/CbcPolyhedronSolver.xtend new file mode 100644 index 00000000..708f93dc --- /dev/null +++ b/Solvers/VIATRA-Solver/hu.bme.mit.inf.dslreasoner.viatrasolver.logic2viatra/src/hu/bme/mit/inf/dslreasoner/viatrasolver/logic2viatra/cardinality/CbcPolyhedronSolver.xtend | |||
@@ -0,0 +1,241 @@ | |||
1 | package hu.bme.mit.inf.dslreasoner.viatrasolver.logic2viatra.cardinality | ||
2 | |||
3 | import com.google.common.collect.ImmutableList | ||
4 | import com.google.common.collect.ImmutableMap | ||
5 | import hu.bme.mit.inf.dslreasoner.ilp.cbc.CbcResult | ||
6 | import hu.bme.mit.inf.dslreasoner.ilp.cbc.CbcSolver | ||
7 | import java.util.HashSet | ||
8 | import java.util.List | ||
9 | import java.util.Map | ||
10 | import java.util.Set | ||
11 | import org.eclipse.xtend.lib.annotations.FinalFieldsConstructor | ||
12 | |||
13 | @FinalFieldsConstructor | ||
14 | class CbcPolyhedronSolver implements PolyhedronSolver { | ||
15 | val boolean lpRelaxation | ||
16 | val double timeoutSeconds | ||
17 | val boolean silent | ||
18 | |||
19 | new() { | ||
20 | this(false, -1, true) | ||
21 | } | ||
22 | |||
23 | override createSaturationOperator(Polyhedron polyhedron) { | ||
24 | new CbcSaturationOperator(polyhedron, lpRelaxation, timeoutSeconds, silent) | ||
25 | } | ||
26 | } | ||
27 | |||
28 | class CbcSaturationOperator extends AbstractPolyhedronSaturationOperator { | ||
29 | static val EPSILON = 1e-6 | ||
30 | |||
31 | val boolean lpRelaxation | ||
32 | val double timeoutSeconds | ||
33 | val boolean silent | ||
34 | val double[] columnLowerBounds | ||
35 | val double[] columnUpperBounds | ||
36 | val double[] objective | ||
37 | val Map<Dimension, Integer> dimensionsToIndicesMap | ||
38 | |||
39 | new(Polyhedron polyhedron, boolean lpRelaxation, double timeoutSeconds, boolean silent) { | ||
40 | super(polyhedron) | ||
41 | this.lpRelaxation = lpRelaxation | ||
42 | this.timeoutSeconds = timeoutSeconds | ||
43 | this.silent = silent | ||
44 | val numDimensions = polyhedron.dimensions.size | ||
45 | columnLowerBounds = newDoubleArrayOfSize(numDimensions) | ||
46 | columnUpperBounds = newDoubleArrayOfSize(numDimensions) | ||
47 | objective = newDoubleArrayOfSize(numDimensions) | ||
48 | dimensionsToIndicesMap = ImmutableMap.copyOf(polyhedron.dimensions.indexed.toMap([value], [key])) | ||
49 | } | ||
50 | |||
51 | override doSaturate() { | ||
52 | val numDimensions = polyhedron.dimensions.size | ||
53 | for (var int j = 0; j < numDimensions; j++) { | ||
54 | val dimension = polyhedron.dimensions.get(j) | ||
55 | columnLowerBounds.set(j, dimension.lowerBound.toDouble(Double.NEGATIVE_INFINITY)) | ||
56 | columnUpperBounds.set(j, dimension.upperBound.toDouble(Double.POSITIVE_INFINITY)) | ||
57 | } | ||
58 | val constraints = nonTrivialConstraints | ||
59 | val numConstraints = constraints.size | ||
60 | val rowStarts = newIntArrayOfSize(numConstraints + 1) | ||
61 | val rowLowerBounds = newDoubleArrayOfSize(numConstraints) | ||
62 | val rowUpperBounds = newDoubleArrayOfSize(numConstraints) | ||
63 | val numEntries = constraints.map[coefficients.size].reduce[a, b|a + b] ?: 0 | ||
64 | rowStarts.set(numConstraints, numEntries) | ||
65 | val columnIndices = newIntArrayOfSize(numEntries) | ||
66 | val entries = newDoubleArrayOfSize(numEntries) | ||
67 | val unconstrainedDimensions = new HashSet | ||
68 | for (dimension : polyhedron.dimensions) { | ||
69 | if (dimension.lowerBound === null && dimension.upperBound === null) { | ||
70 | unconstrainedDimensions += dimension | ||
71 | } | ||
72 | } | ||
73 | var int index = 0 | ||
74 | for (var int i = 0; i < numConstraints; i++) { | ||
75 | rowStarts.set(i, index) | ||
76 | val constraint = constraints.get(i) | ||
77 | rowLowerBounds.set(i, constraint.lowerBound.toDouble(Double.NEGATIVE_INFINITY)) | ||
78 | rowUpperBounds.set(i, constraint.upperBound.toDouble(Double.POSITIVE_INFINITY)) | ||
79 | if (!dimensionsToIndicesMap.keySet.containsAll(constraint.coefficients.keySet)) { | ||
80 | throw new IllegalArgumentException("Constrains has unknown dimensions") | ||
81 | } | ||
82 | for (var int j = 0; j < numDimensions; j++) { | ||
83 | val dimension = polyhedron.dimensions.get(j) | ||
84 | val coefficient = constraint.coefficients.get(dimension) | ||
85 | if (coefficient !== null && coefficient != 0) { | ||
86 | unconstrainedDimensions -= dimension | ||
87 | columnIndices.set(index, j) | ||
88 | entries.set(index, coefficient) | ||
89 | index++ | ||
90 | } | ||
91 | } | ||
92 | } | ||
93 | if (index != numEntries) { | ||
94 | throw new AssertionError("Last entry does not equal the number of entries in the constraint matrix") | ||
95 | } | ||
96 | for (expressionToSaturate : polyhedron.expressionsToSaturate) { | ||
97 | val result = saturate(expressionToSaturate, rowStarts, columnIndices, entries, rowLowerBounds, | ||
98 | rowUpperBounds, unconstrainedDimensions, constraints) | ||
99 | if (result != PolyhedronSaturationResult.SATURATED) { | ||
100 | return result | ||
101 | } | ||
102 | } | ||
103 | PolyhedronSaturationResult.SATURATED | ||
104 | } | ||
105 | |||
106 | protected def saturate(LinearBoundedExpression expressionToSaturate, int[] rowStarts, int[] columnIndices, | ||
107 | double[] entries, double[] rowLowerBounds, double[] rowUpperBounds, Set<Dimension> unconstrainedDimensions, | ||
108 | ImmutableList<LinearConstraint> constraints) { | ||
109 | val numDimensions = objective.size | ||
110 | for (var int j = 0; j < numDimensions; j++) { | ||
111 | objective.set(j, 0) | ||
112 | } | ||
113 | switch (expressionToSaturate) { | ||
114 | Dimension: { | ||
115 | // CBC will return nonsensical results or call free() with invalid arguments if | ||
116 | // it is passed a fully unconstrained (-Inf lower and +Int upper bound, no inequalities) variable | ||
117 | // in the objective function. | ||
118 | if (unconstrainedDimensions.contains(expressionToSaturate)) { | ||
119 | return PolyhedronSaturationResult.SATURATED | ||
120 | } | ||
121 | val j = getIndex(expressionToSaturate) | ||
122 | objective.set(j, 1) | ||
123 | } | ||
124 | LinearConstraint: { | ||
125 | for (pair : expressionToSaturate.coefficients.entrySet) { | ||
126 | val dimension = pair.key | ||
127 | // We also have to check for unconstrained dimensions here to avoid crashing. | ||
128 | if (unconstrainedDimensions.contains(dimension)) { | ||
129 | expressionToSaturate.lowerBound = null | ||
130 | expressionToSaturate.upperBound = null | ||
131 | return PolyhedronSaturationResult.SATURATED | ||
132 | } | ||
133 | val j = getIndex(dimension) | ||
134 | objective.set(j, pair.value) | ||
135 | } | ||
136 | } | ||
137 | default: | ||
138 | throw new IllegalArgumentException("Unknown expression: " + expressionToSaturate) | ||
139 | } | ||
140 | val minimizationResult = CbcSolver.solve(columnLowerBounds, columnUpperBounds, rowStarts, columnIndices, | ||
141 | entries, rowLowerBounds, rowUpperBounds, objective, lpRelaxation, timeoutSeconds, silent) | ||
142 | switch (minimizationResult) { | ||
143 | CbcResult.SolutionBounded: { | ||
144 | val doubleValue = minimizationResult.value | ||
145 | val roundedValue = Math.ceil(doubleValue - EPSILON) | ||
146 | val intValue = roundedValue as int | ||
147 | val oldBound = expressionToSaturate.lowerBound | ||
148 | if (oldBound === null || intValue >= oldBound) { | ||
149 | expressionToSaturate.lowerBound = intValue | ||
150 | setBound(expressionToSaturate, constraints, roundedValue, columnLowerBounds, rowLowerBounds) | ||
151 | } else { | ||
152 | throw new IllegalStateException("Unexpected decrease of lower bound by " + (oldBound - doubleValue)) | ||
153 | } | ||
154 | } | ||
155 | case CbcResult.SOLUTION_UNBOUNDED: { | ||
156 | if (expressionToSaturate.lowerBound !== null) { | ||
157 | throw new IllegalStateException("Finite lower bound became infinite") | ||
158 | } | ||
159 | setBound(expressionToSaturate, constraints, Double.NEGATIVE_INFINITY, columnLowerBounds, rowLowerBounds) | ||
160 | } | ||
161 | case CbcResult.UNSAT: | ||
162 | return PolyhedronSaturationResult.EMPTY | ||
163 | case CbcResult.ABANDONED, | ||
164 | case CbcResult.TIMEOUT: | ||
165 | return PolyhedronSaturationResult.UNKNOWN | ||
166 | default: | ||
167 | throw new RuntimeException("Unknown CbcResult: " + minimizationResult) | ||
168 | } | ||
169 | for (var int j = 0; j < numDimensions; j++) { | ||
170 | val objectiveCoefficient = objective.get(j) | ||
171 | objective.set(j, -objectiveCoefficient) | ||
172 | } | ||
173 | val maximizationResult = CbcSolver.solve(columnLowerBounds, columnUpperBounds, rowStarts, columnIndices, | ||
174 | entries, rowLowerBounds, rowUpperBounds, objective, lpRelaxation, timeoutSeconds, silent) | ||
175 | switch (maximizationResult) { | ||
176 | CbcResult.SolutionBounded: { | ||
177 | val doubleValue = -maximizationResult.value | ||
178 | val roundedValue = Math.floor(doubleValue + EPSILON) | ||
179 | val intValue = roundedValue as int | ||
180 | val oldBound = expressionToSaturate.upperBound | ||
181 | if (oldBound === null || intValue <= oldBound) { | ||
182 | expressionToSaturate.upperBound = intValue | ||
183 | setBound(expressionToSaturate, constraints, roundedValue, columnUpperBounds, rowUpperBounds) | ||
184 | } else { | ||
185 | throw new IllegalStateException("Unexpected increase of upper bound by " + (doubleValue - oldBound)) | ||
186 | } | ||
187 | } | ||
188 | case CbcResult.SOLUTION_UNBOUNDED: { | ||
189 | if (expressionToSaturate.lowerBound !== null) { | ||
190 | throw new IllegalStateException("Finite upper bound became infinite") | ||
191 | } | ||
192 | expressionToSaturate.upperBound = null | ||
193 | setBound(expressionToSaturate, constraints, Double.POSITIVE_INFINITY, columnUpperBounds, rowUpperBounds) | ||
194 | } | ||
195 | case CbcResult.UNSAT: | ||
196 | if (lpRelaxation) { | ||
197 | return PolyhedronSaturationResult.EMPTY | ||
198 | } else { | ||
199 | throw new RuntimeException("Minimization was SAT, but maximization is UNSAT") | ||
200 | } | ||
201 | case CbcResult.ABANDONED, | ||
202 | case CbcResult.TIMEOUT: | ||
203 | return PolyhedronSaturationResult.UNKNOWN | ||
204 | default: | ||
205 | throw new RuntimeException("Unknown CbcResult: " + maximizationResult) | ||
206 | } | ||
207 | return PolyhedronSaturationResult.SATURATED | ||
208 | } | ||
209 | |||
210 | private def toDouble(Integer nullableInt, double defaultValue) { | ||
211 | if (nullableInt === null) { | ||
212 | defaultValue | ||
213 | } else { | ||
214 | nullableInt.doubleValue | ||
215 | } | ||
216 | } | ||
217 | |||
218 | private def int getIndex(Dimension dimension) { | ||
219 | val index = dimensionsToIndicesMap.get(dimension) | ||
220 | if (index === null) { | ||
221 | throw new IllegalArgumentException("Unknown dimension: " + dimension) | ||
222 | } | ||
223 | index | ||
224 | } | ||
225 | |||
226 | private def void setBound(LinearBoundedExpression expression, List<LinearConstraint> constraints, double bound, | ||
227 | double[] columnBounds, double[] rowBounds) { | ||
228 | switch (expression) { | ||
229 | Dimension: { | ||
230 | val j = getIndex(expression) | ||
231 | columnBounds.set(j, bound) | ||
232 | } | ||
233 | LinearConstraint: { | ||
234 | val i = constraints.indexOf(expression) | ||
235 | if (i >= 0) { | ||
236 | rowBounds.set(i, bound) | ||
237 | } | ||
238 | } | ||
239 | } | ||
240 | } | ||
241 | } | ||