1 /*
2 * Licensed to the Apache Software Foundation (ASF) under one or more
3 * contributor license agreements. See the NOTICE file distributed with
4 * this work for additional information regarding copyright ownership.
5 * The ASF licenses this file to You under the Apache License, Version 2.0
6 * (the "License"); you may not use this file except in compliance with
7 * the License. You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 */
17
18 package org.apache.commons.math.optimization.linear;
19
20 import org.apache.commons.math.optimization.OptimizationException;
21 import org.apache.commons.math.optimization.RealPointValuePair;
22 import org.apache.commons.math.util.MathUtils;
23
24
25 /**
26 * Solves a linear problem using the Two-Phase Simplex Method.
27 * @version $Revision: 797801 $ $Date: 2009-07-25 13:22:06 -0400 (Sat, 25 Jul 2009) $
28 * @since 2.0
29 */
30 public class SimplexSolver extends AbstractLinearOptimizer {
31
32 /** Default amount of error to accept in floating point comparisons. */
33 private static final double DEFAULT_EPSILON = 1.0e-6;
34
35 /** Amount of error to accept in floating point comparisons. */
36 protected final double epsilon;
37
38 /**
39 * Build a simplex solver with default settings.
40 */
41 public SimplexSolver() {
42 this(DEFAULT_EPSILON);
43 }
44
45 /**
46 * Build a simplex solver with a specified accepted amount of error
47 * @param epsilon the amount of error to accept in floating point comparisons
48 */
49 public SimplexSolver(final double epsilon) {
50 this.epsilon = epsilon;
51 }
52
53 /**
54 * Returns the column with the most negative coefficient in the objective function row.
55 * @param tableau simple tableau for the problem
56 * @return column with the most negative coefficient
57 */
58 private Integer getPivotColumn(SimplexTableau tableau) {
59 double minValue = 0;
60 Integer minPos = null;
61 for (int i = tableau.getNumObjectiveFunctions(); i < tableau.getWidth() - 1; i++) {
62 if (MathUtils.compareTo(tableau.getEntry(0, i), minValue, epsilon) < 0) {
63 minValue = tableau.getEntry(0, i);
64 minPos = i;
65 }
66 }
67 return minPos;
68 }
69
70 /**
71 * Returns the row with the minimum ratio as given by the minimum ratio test (MRT).
72 * @param tableau simple tableau for the problem
73 * @param col the column to test the ratio of. See {@link #getPivotColumn(SimplexTableau)}
74 * @return row with the minimum ratio
75 */
76 private Integer getPivotRow(final int col, final SimplexTableau tableau) {
77 double minRatio = Double.MAX_VALUE;
78 Integer minRatioPos = null;
79 for (int i = tableau.getNumObjectiveFunctions(); i < tableau.getHeight(); i++) {
80 double rhs = tableau.getEntry(i, tableau.getWidth() - 1);
81 if (MathUtils.compareTo(tableau.getEntry(i, col), 0, epsilon) >= 0) {
82 double ratio = rhs / tableau.getEntry(i, col);
83 if (ratio < minRatio) {
84 minRatio = ratio;
85 minRatioPos = i;
86 }
87 }
88 }
89 return minRatioPos;
90 }
91
92
93 /**
94 * Runs one iteration of the Simplex method on the given model.
95 * @param tableau simple tableau for the problem
96 * @throws OptimizationException if the maximal iteration count has been
97 * exceeded or if the model is found not to have a bounded solution
98 */
99 protected void doIteration(final SimplexTableau tableau)
100 throws OptimizationException {
101
102 incrementIterationsCounter();
103
104 Integer pivotCol = getPivotColumn(tableau);
105 Integer pivotRow = getPivotRow(pivotCol, tableau);
106 if (pivotRow == null) {
107 throw new UnboundedSolutionException();
108 }
109
110 // set the pivot element to 1
111 double pivotVal = tableau.getEntry(pivotRow, pivotCol);
112 tableau.divideRow(pivotRow, pivotVal);
113
114 // set the rest of the pivot column to 0
115 for (int i = 0; i < tableau.getHeight(); i++) {
116 if (i != pivotRow) {
117 double multiplier = tableau.getEntry(i, pivotCol);
118 tableau.subtractRow(i, pivotRow, multiplier);
119 }
120 }
121 }
122
123 /**
124 * Checks whether Phase 1 is solved.
125 * @param tableau simple tableau for the problem
126 * @return whether Phase 1 is solved
127 */
128 private boolean isPhase1Solved(final SimplexTableau tableau) {
129 if (tableau.getNumArtificialVariables() == 0) {
130 return true;
131 }
132 for (int i = tableau.getNumObjectiveFunctions(); i < tableau.getWidth() - 1; i++) {
133 if (MathUtils.compareTo(tableau.getEntry(0, i), 0, epsilon) < 0) {
134 return false;
135 }
136 }
137 return true;
138 }
139
140 /**
141 * Returns whether the problem is at an optimal state.
142 * @param tableau simple tableau for the problem
143 * @return whether the model has been solved
144 */
145 public boolean isOptimal(final SimplexTableau tableau) {
146 if (tableau.getNumArtificialVariables() > 0) {
147 return false;
148 }
149 for (int i = tableau.getNumObjectiveFunctions(); i < tableau.getWidth() - 1; i++) {
150 if (MathUtils.compareTo(tableau.getEntry(0, i), 0, epsilon) < 0) {
151 return false;
152 }
153 }
154 return true;
155 }
156
157 /**
158 * Solves Phase 1 of the Simplex method.
159 * @param tableau simple tableau for the problem
160 * @exception OptimizationException if the maximal number of iterations is
161 * exceeded, or if the problem is found not to have a bounded solution, or
162 * if there is no feasible solution
163 */
164 protected void solvePhase1(final SimplexTableau tableau)
165 throws OptimizationException {
166 // make sure we're in Phase 1
167 if (tableau.getNumArtificialVariables() == 0) {
168 return;
169 }
170
171 while (!isPhase1Solved(tableau)) {
172 doIteration(tableau);
173 }
174
175 // if W is not zero then we have no feasible solution
176 if (!MathUtils.equals(tableau.getEntry(0, tableau.getRhsOffset()), 0, epsilon)) {
177 throw new NoFeasibleSolutionException();
178 }
179 }
180
181 /** {@inheritDoc} */
182 @Override
183 public RealPointValuePair doOptimize()
184 throws OptimizationException {
185 final SimplexTableau tableau =
186 new SimplexTableau(f, constraints, goalType, restrictToNonNegative, epsilon);
187 solvePhase1(tableau);
188 tableau.discardArtificialVariables();
189 while (!isOptimal(tableau)) {
190 doIteration(tableau);
191 }
192 return tableau.getSolution();
193 }
194
195 }