diff --git a/.gitignore b/.gitignore
index 299e39ac618542a86dda8dd4615047b4d795328b..b0a9905575783f18d92dd5f505a52d187802980d 100644
--- a/.gitignore
+++ b/.gitignore
@@ -4,4 +4,5 @@ build/
target/
*.iml
.classpath
-.project
\ No newline at end of file
+.project
+**/.DS_Store
\ No newline at end of file
diff --git a/pom.xml b/pom.xml
index 83c98708069ca9c1b11e44218e9219c2810a53c2..18fcc240c7888a84cc7b677b3c3cbdb759e63f20 100644
--- a/pom.xml
+++ b/pom.xml
@@ -1,19 +1,17 @@
<?xml version="1.0" encoding="UTF-8"?>
<project xmlns="http://maven.apache.org/POM/4.0.0"
- xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
- xsi:schemaLocation="http://maven.apache.org/POM/4.0.0 http://maven.apache.org/xsd/maven-4.0.0.xsd">
+ xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
+ xsi:schemaLocation="http://maven.apache.org/POM/4.0.0 http://maven.apache.org/xsd/maven-4.0.0.xsd">
<modelVersion>4.0.0</modelVersion>
<parent>
<groupId>org.bioimageanalysis.icy</groupId>
- <artifactId>parent-pom-plugin</artifactId>
- <version>1.0.3</version>
+ <artifactId>pom-icy</artifactId>
+ <version>2.2.0</version>
</parent>
<artifactId>linear-programming</artifactId>
- <version>1.0.2</version>
-
- <packaging>jar</packaging>
+ <version>2.0.0</version>
<name>Linear Programming</name>
<description>
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalProgramParameters.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalProgramParameters.java
index 2161372aebe27fcbc652c4229a7e644be1b82f90..2c85ae237efeae731d28ffcea9fe379133b9b20c 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalProgramParameters.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalProgramParameters.java
@@ -1,16 +1,34 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
/**
* Parameters used for the simplex algorithm
- *
+ * <p>
* By definition, the optimization problem is of the form:
- *
+ *
* <pre>
* min_x c'*x or max_x c'*x
* A*x <= b
* with x >=0
* </pre>
- *
+ * <p>
* where x, b, c are column vectors and c' is the transpose of c. * stands for
* the matrix multiplication.
* <p>
@@ -21,39 +39,38 @@ package plugins.nchenouard.linearprogrammingfullsimplex;
* <p>
* A[i]*x <= b if equalityConstraints[i]
* <p>
- *
- *
+ * <p>
+ * <p>
* Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
-public class CanonicalProgramParameters
-{
- /**
- * The objective function coefficient: c'*x or sum_i c[i]*x[i] for the variable x.
- * */
- public double[] c;
- /**
- * Defines whether the problem is a maximization or a minimization with respect to x.
- * */
- public boolean maximization;
+public class CanonicalProgramParameters {
+ /**
+ * The objective function coefficient: c'*x or sum_i c[i]*x[i] for the variable x.
+ */
+ public double[] c;
+ /**
+ * Defines whether the problem is a maximization or a minimization with respect to x.
+ */
+ public boolean maximization;
+
+ /**
+ * The constraint matrix structure. Each line defines a linear combination of the input variable.
+ */
+ public double[][] A;
+ /**
+ * The value for each constraint.
+ */
+ public double[] b;
- /**
- * The constraint matrix structure. Each line defines a linear combination of the input variable.
- */
- public double[][] A;
- /**
- * The value for each constraint.
- */
- public double[] b;
-
- /**
- * Defines the type of constraint
- * <p>
- * A[i]*x == b if equalityConstraints[i]
- * <p>
- * A[i]*x <= b if equalityConstraints[i]
- */
- public boolean[] equalityConstraints;
+ /**
+ * Defines the type of constraint
+ * <p>
+ * A[i]*x == b if equalityConstraints[i]
+ * <p>
+ * A[i]*x <= b if equalityConstraints[i]
+ */
+ public boolean[] equalityConstraints;
}
\ No newline at end of file
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalSimplexProgram.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalSimplexProgram.java
index 451e9c12f39a09e2cff73125a30a3e67f472b34b..d54cfb9fa1386b403badc4e3b720f29f5a4c546c 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalSimplexProgram.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/CanonicalSimplexProgram.java
@@ -1,24 +1,37 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
-import java.io.BufferedReader;
-import java.io.BufferedWriter;
-import java.io.File;
-import java.io.FileInputStream;
-import java.io.FileNotFoundException;
-import java.io.FileWriter;
-import java.io.IOException;
-import java.io.InputStreamReader;
+import icy.system.IcyExceptionHandler;
+
+import java.io.*;
/**
* Linear program for the Simplex algorithm. Inheriting classes should implement
* the pivoting rule.
- *
+ * <p>
* By definition, the optimization problem is of the form: min_x c'*x or max_x
* c'*x A*x <= b with x >=0 where x, b, c are column vectors and c' is the
* transpose of c.
- *
+ * <p>
* * stands for the matrix multiplication.
- *
+ * <p>
* Here each constraint <code>A[i]*x</code> can be either "equal" or "lower than
* or equal" constraints based on the boolean value in the equality constraints
* vector:
@@ -27,231 +40,207 @@ import java.io.InputStreamReader;
* <p>
* A[i]*x <= b if equalityConstraints[i]
* <p>
- *
+ * <p>
* Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
-public abstract class CanonicalSimplexProgram
-{
- /**
- * Parameters defining the linear programming problem
- * */
- CanonicalProgramParameters parameters;
- /**
- * Number of constraints
- * */
- int numConstraints;
- /**
- * Number of variables
- * */
- int numVariables;
- /**
- * Solution obtained with the Simplex algorithm
- * */
- double[] solution;
- /**
- * Intermediary tableau used for some types of problems
- * */
- TableauWithSlackVariables tableau0 = null;
- /**
- * Numerical accuracy limit under which a number is considered to be 0
- * */
- double tol = 1e-12;
- /**
- * Display detail of steps
- * */
- public boolean verbose = false;
+public abstract class CanonicalSimplexProgram {
+ /**
+ * Parameters defining the linear programming problem
+ */
+ CanonicalProgramParameters parameters;
+ /**
+ * Number of constraints
+ */
+ int numConstraints;
+ /**
+ * Number of variables
+ */
+ int numVariables;
+ /**
+ * Solution obtained with the Simplex algorithm
+ */
+ double[] solution;
+ /**
+ * Intermediary tableau used for some types of problems
+ */
+ TableauWithSlackVariables tableau0 = null;
+ /**
+ * Numerical accuracy limit under which a number is considered to be 0
+ */
+ double tol = 1e-12;
+ /**
+ * Display detail of steps
+ */
+ public boolean verbose = false;
+
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ public CanonicalSimplexProgram(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ this.parameters = new CanonicalProgramParameters();
+ this.parameters.A = A;
+ this.parameters.b = b;
+ this.parameters.c = c;
+ this.parameters.maximization = maximization;
+ this.parameters.equalityConstraints = equality;
+ this.numVariables = c.length;
+ this.numConstraints = b.length;
+ }
- /**
- * Constructor specifying parameters of the linear programming problem.
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- public CanonicalSimplexProgram(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality)
- {
- this.parameters = new CanonicalProgramParameters();
- this.parameters.A = A;
- this.parameters.b = b;
- this.parameters.c = c;
- this.parameters.maximization = maximization;
- this.parameters.equalityConstraints = equality;
- this.numVariables = c.length;
- this.numConstraints = b.length;
- }
- /**
- * Constructor specifying parameters of the linear programming problem.
- * @param p Parameters of the linear programming problem.
- */
- public CanonicalSimplexProgram(final CanonicalProgramParameters p)
- {
- this.parameters = p;
- this.numVariables = p.c.length;
- this.numConstraints = p.b.length;
- }
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param p Parameters of the linear programming problem.
+ */
+ public CanonicalSimplexProgram(final CanonicalProgramParameters p) {
+ this.parameters = p;
+ this.numVariables = p.c.length;
+ this.numConstraints = p.b.length;
+ }
- /**
- * Solve the linear programming problem using the Simplex algorithm for the primal problem
- * @return true if optimization succeeded, false otherwise.
- * */
- public boolean solvePrimalSimplex()
- {
- final double[][] A = parameters.A;
- final double[] b = parameters.b;
- final boolean[] eq = parameters.equalityConstraints;
- final double[] c = parameters.c;
+ /**
+ * Solve the linear programming problem using the Simplex algorithm for the primal problem
+ *
+ * @return true if optimization succeeded, false otherwise.
+ */
+ public boolean solvePrimalSimplex() {
+ final double[][] A = parameters.A;
+ final double[] b = parameters.b;
+ final boolean[] eq = parameters.equalityConstraints;
+ final double[] c = parameters.c;
- boolean onlyInequality = true;
- for (final boolean e:eq)
- {
- if (e)
- {
- onlyInequality = false;
- break;
- }
- }
- boolean hasNegativeRHS = false;
- for (int i = 0; i < numConstraints; i++)
- {
- if (b[i] < 0)
- {
- if (eq[i]) // just change the signs for the constraint
- {
- b[i] = -b[i];
- for (int j = 0; j < A[i].length; j++)
- A[i][j] = -A[i][j];
- }
- else
- hasNegativeRHS = true;
- }
- }
- if (onlyInequality)
- {
- if(hasNegativeRHS)
- {
- // add slack variables and multiply by -1 when corresponding to a negative right-hand-side
- final double[][] A2 = new double[numConstraints][numConstraints + numVariables];
- final double[] c2 = new double[numConstraints + numVariables];
- System.arraycopy(c, 0, c2, 0, numVariables); // no cost for slack variables
- final boolean[] equality = new boolean[numConstraints];
- for (int i = 0; i < numConstraints; i++)
- {
- equality[i] = true;
- if (b[i] >= 0)
- {
- System.arraycopy(A[i], 0, A2[i], 0, numVariables);
- A2[i][numVariables + i] = 1;
- }
- else
- {
- for (int j = 0; j < numVariables; j++)
- A2[i][j] = -A[i][j];
- A2[i][numVariables + i] = -1;
- b[i] = -b[i];
- }
- }
- parameters.A = A2;
- parameters.c = c2;
- parameters.equalityConstraints = equality;
- // we have now a problem with only equality constraints and non-negative right-hand-side
- if(solveWithEqualities(false))
- {
- final double[] tmpSolution = this.solution;
- this.solution = new double[numVariables];
- System.arraycopy(tmpSolution, 0, solution, 0, numVariables);
- }
- else
- return false;
- }
- else
- solveWithInequalities(false);
- }
- else
- {
- if(hasNegativeRHS)
- {
- // we need to solve first a problem to find the initial tableau
- // we will add slack variables and minimize the cost of the sum of the slack variables for the equality constraints
- final boolean maximization = false;
- final boolean[] eq2 = new boolean[numConstraints];
- final boolean[] isNegativeConstraint = new boolean[numConstraints];
- for (int i = 0; i < numConstraints; i++)
- {
- if (b[i] < 0)
- {
- isNegativeConstraint[i] = true;
- b[i] = -b[i];
- for (int j = 0; j < numVariables; j++)
- A[i][j] = -A[i][j];
- }
- }
- final CanonicalSimplexProgram modifiedProgram = createNewProgram(A, b, c, maximization, eq2);
- if (modifiedProgram.solvePhaseICanonicalSimplex())
- {
- modifiedProgram.tableau0.setUnitScoreToSlackVars(eq2);
- if (modifiedProgram.solvePhaseIICanonicalSimplex(modifiedProgram.tableau0, maximization))
- {
- if (Math.abs(modifiedProgram.tableau0.scoreValue) < tol) // optimal score should be 0
- {
- // we have found an initial tableau for the main problem
- // we should modify the scores now for the tableau
- this.tableau0 = modifiedProgram.tableau0;
- tableau0.modifyScores(c);
- return true;
- }
- else
- return false;
- }
- else
- return false;
- }
- }
- else
- {
- return solveWithEqualities(false);
- }
- }
- return true;
- }
+ boolean onlyInequality = true;
+ for (final boolean e : eq) {
+ if (e) {
+ onlyInequality = false;
+ break;
+ }
+ }
+ boolean hasNegativeRHS = false;
+ for (int i = 0; i < numConstraints; i++) {
+ if (b[i] < 0) {
+ if (eq[i]) // just change the signs for the constraint
+ {
+ b[i] = -b[i];
+ for (int j = 0; j < A[i].length; j++)
+ A[i][j] = -A[i][j];
+ }
+ else
+ hasNegativeRHS = true;
+ }
+ }
+ if (onlyInequality) {
+ if (hasNegativeRHS) {
+ // add slack variables and multiply by -1 when corresponding to a negative right-hand-side
+ final double[][] A2 = new double[numConstraints][numConstraints + numVariables];
+ final double[] c2 = new double[numConstraints + numVariables];
+ System.arraycopy(c, 0, c2, 0, numVariables); // no cost for slack variables
+ final boolean[] equality = new boolean[numConstraints];
+ for (int i = 0; i < numConstraints; i++) {
+ equality[i] = true;
+ if (b[i] >= 0) {
+ System.arraycopy(A[i], 0, A2[i], 0, numVariables);
+ A2[i][numVariables + i] = 1;
+ }
+ else {
+ for (int j = 0; j < numVariables; j++)
+ A2[i][j] = -A[i][j];
+ A2[i][numVariables + i] = -1;
+ b[i] = -b[i];
+ }
+ }
+ parameters.A = A2;
+ parameters.c = c2;
+ parameters.equalityConstraints = equality;
+ // we have now a problem with only equality constraints and non-negative right-hand-side
+ if (solveWithEqualities(false)) {
+ final double[] tmpSolution = this.solution;
+ this.solution = new double[numVariables];
+ System.arraycopy(tmpSolution, 0, solution, 0, numVariables);
+ }
+ else
+ return false;
+ }
+ else
+ solveWithInequalities(false);
+ }
+ else {
+ if (hasNegativeRHS) {
+ // we need to solve first a problem to find the initial tableau
+ // we will add slack variables and minimize the cost of the sum of the slack variables for the equality constraints
+ final boolean maximization = false;
+ final boolean[] eq2 = new boolean[numConstraints];
+ final boolean[] isNegativeConstraint = new boolean[numConstraints];
+ for (int i = 0; i < numConstraints; i++) {
+ if (b[i] < 0) {
+ isNegativeConstraint[i] = true;
+ b[i] = -b[i];
+ for (int j = 0; j < numVariables; j++)
+ A[i][j] = -A[i][j];
+ }
+ }
+ final CanonicalSimplexProgram modifiedProgram = createNewProgram(A, b, c, maximization, eq2);
+ if (modifiedProgram.solvePhaseICanonicalSimplex()) {
+ modifiedProgram.tableau0.setUnitScoreToSlackVars(eq2);
+ if (modifiedProgram.solvePhaseIICanonicalSimplex(modifiedProgram.tableau0, maximization)) {
+ if (Math.abs(modifiedProgram.tableau0.scoreValue) < tol) // optimal score should be 0
+ {
+ // we have found an initial tableau for the main problem
+ // we should modify the scores now for the tableau
+ this.tableau0 = modifiedProgram.tableau0;
+ tableau0.modifyScores(c);
+ return true;
+ }
+ else
+ return false;
+ }
+ else
+ return false;
+ }
+ }
+ else {
+ return solveWithEqualities(false);
+ }
+ }
+ return true;
+ }
- /**
- * Solve a problem that contains only inequality constraints
- * @return true if optimization succeeded, false otherwise.
- * */
- private boolean solveWithInequalities(final boolean compareWithLPSolve)
- {
- // Canonical form: initial tableau is built by introducing unit vectors
- solvePhaseICanonicalSimplex();
- if (verbose)
- {
- System.out.println("Tableau after phase 1 :");
- tableau0.printTableau();
- }
- // Phase II: optimisation by Gauss-Jordan replacement of non-basic columns
- return solvePhaseIICanonicalSimplex(tableau0, parameters.maximization);
- }
+ /**
+ * Solve a problem that contains only inequality constraints
+ *
+ * @return true if optimization succeeded, false otherwise.
+ */
+ private boolean solveWithInequalities(final boolean compareWithLPSolve) {
+ // Canonical form: initial tableau is built by introducing unit vectors
+ solvePhaseICanonicalSimplex();
+ if (verbose) {
+ System.out.println("Tableau after phase 1 :");
+ tableau0.printTableau();
+ }
+ // Phase II: optimisation by Gauss-Jordan replacement of non-basic columns
+ return solvePhaseIICanonicalSimplex(tableau0, parameters.maximization);
+ }
- /**
- * Solve a problem that contains some equality constraints
- * @return true if optimization succeeded, false otherwise.
- * */
- private boolean solveWithEqualities(final boolean compareWithLPSolve)
- {
- tableau0 = new TableauMinSlackObjective(this);
- if (verbose)
- tableau0.printTableau();
- // check for some duplicates in tableau0 columns
+ /**
+ * Solve a problem that contains some equality constraints
+ *
+ * @return true if optimization succeeded, false otherwise.
+ */
+ private boolean solveWithEqualities(final boolean compareWithLPSolve) {
+ tableau0 = new TableauMinSlackObjective(this);
+ if (verbose)
+ tableau0.printTableau();
+ // check for some duplicates in tableau0 columns
// for (int i = 0; i < tableau0.columnVectors.length; i++)
// {
// for (int ii = 0; ii < tableau0.columnVectors.length; ii++)
@@ -269,777 +258,699 @@ public abstract class CanonicalSimplexProgram
// }
// }
// }
- final boolean success = solvePhaseIICanonicalSimplex(tableau0, false); // minimization of the init tableau
- // {
- // double c[] = new double[tableau0.numCol];
- // for (int i = 0; i < numConstraints; i ++)
- // if (parameters.equalityConstraints[i])
- // c[i] = 1;//put weights to slack variables which correspond to equalities
- //// tableau0.solverWithLPSolve(c, b0, false);
- // }
- if (!success)
- return false;
- if (verbose)
- {
- System.out.println("Phase 1 ");
- tableau0.printTableau();
- }
- if (Math.abs(tableau0.scoreValue) < tol) // optimal score should be 0
- {
- // check that the initial solution is valid
- // for (int i = 0; i < A0.length; i++)
- // {
- // double v = 0;
- // double[] s = tableau0.getSolution();
- // // compute the constraint
- // for (int j = 0; j < A0[i].length; j++)
- // v += A0[i][j]*s[j];
- // if (equality[i])
- // System.out.println(v+" == "+b0[i]);
- // else
- // System.out.println(v+" <= "+b0[i]);
- // }
- // modifies the score for the objective function
- tableau0.modifyScoresBis(parameters.c);
- }
- else
- return false;
- // phase II: optimization by Gauss-Jordan replacement of non-basic columns
- return solvePhaseIICanonicalSimplex(tableau0, parameters.maximization);
- }
+ final boolean success = solvePhaseIICanonicalSimplex(tableau0, false); // minimization of the init tableau
+ // {
+ // double c[] = new double[tableau0.numCol];
+ // for (int i = 0; i < numConstraints; i ++)
+ // if (parameters.equalityConstraints[i])
+ // c[i] = 1;//put weights to slack variables which correspond to equalities
+ //// tableau0.solverWithLPSolve(c, b0, false);
+ // }
+ if (!success)
+ return false;
+ if (verbose) {
+ System.out.println("Phase 1 ");
+ tableau0.printTableau();
+ }
+ if (Math.abs(tableau0.scoreValue) < tol) // optimal score should be 0
+ {
+ // check that the initial solution is valid
+ // for (int i = 0; i < A0.length; i++)
+ // {
+ // double v = 0;
+ // double[] s = tableau0.getSolution();
+ // // compute the constraint
+ // for (int j = 0; j < A0[i].length; j++)
+ // v += A0[i][j]*s[j];
+ // if (equality[i])
+ // System.out.println(v+" == "+b0[i]);
+ // else
+ // System.out.println(v+" <= "+b0[i]);
+ // }
+ // modifies the score for the objective function
+ tableau0.modifyScoresBis(parameters.c);
+ }
+ else
+ return false;
+ // phase II: optimization by Gauss-Jordan replacement of non-basic columns
+ return solvePhaseIICanonicalSimplex(tableau0, parameters.maximization);
+ }
+
+ /**
+ * Solve the phase I of the Simplex algorithm for canonical problem: no equality constraints nor negative constraints
+ *
+ * @return true if optimization succeeded, false otherwise.
+ */
+ public boolean solvePhaseICanonicalSimplex() {
+ // phase I: initial solution with slack variables
+ tableau0 = new TableauWithSlackVariables(this.parameters);
+ return true;
+ }
+
+ /**
+ * Solve the phase II of the Simplex algorithm for canonical problem: no equality constraints nor negative constraints
+ *
+ * @return true if optimization succeeded, false otherwise.
+ */
+ public boolean solvePhaseIICanonicalSimplex(final TableauWithSlackVariables tableau, final boolean max) {
+ if (tableau == null)
+ return false;
+
+ // Gauss-Jordan simplex method
+ boolean feasible = true;
+ boolean reachedOptimum = false;
+
+ while (feasible && !reachedOptimum) {
+ // start optimization
+ final int colIdx = getPivotColumn(tableau, max);
+ if (colIdx < 0)// optimal solution achieved
+ {
+ reachedOptimum = true;
+ break;
+ }
+ else {
+ final int rowIdx = getRowidx(tableau, colIdx);
+ if (rowIdx < 0) // not feasible
+ {
+ feasible = false;
+ break;
+ }
+ else {
+ if (verbose)
+ System.out.println("replacement: row = " + rowIdx + ", column = " + colIdx);
+ tableau.pivot(colIdx, rowIdx);
+ }
+ }
+ if (verbose) {
+ System.out.println("Tableau at iteration:");
+ tableau.printTableau();
+ }
+ }
+ if (verbose) {
+ System.out.println("Final tableau:");
+ tableau.printTableau();
+ }
+ if (!feasible) {
+ if (verbose)
+ System.out.println("UNFEASIBLE problem");
+ return false;
+ }
+ else {
+ solution = tableau.getSolution();
+ }
+ return true;
+ }
+
+
+ /**
+ * Duplicate the program
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ protected abstract CanonicalSimplexProgram createNewProgram(double[][] A, double[] b, double[] c, boolean maximization, boolean[] equality);
+
+ /**
+ * Get the row index for pivoting
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param colIdx index of the pivoting column in the tableau
+ * @return the row index in the tableau
+ */
+ protected abstract int getRowidx(TableauWithSlackVariables tableau, int colIdx);
+
+ /**
+ * Get the column index for pivoting
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param max true if maximization problem, false if minization problem
+ * @return the column index in the tableau
+ */
+ protected abstract int getPivotColumn(TableauWithSlackVariables tableau, boolean max);
- /**
- * Solve the phase I of the Simplex algorithm for canonical problem: no equality constraints nor negative constraints
- * @return true if optimization succeeded, false otherwise.
- * */
- public boolean solvePhaseICanonicalSimplex()
- {
- // phase I: initial solution with slack variables
- tableau0 = new TableauWithSlackVariables(this.parameters);
- return true;
- }
+ /**
+ * Display whether the computed solution satisfies the constraints
+ */
+ public void checkSolution() {
+ final double[][] A0 = parameters.A;
+ final double[] b0 = parameters.b;
+ if (solution != null) {
+ for (int i = 0; i < numConstraints; i++) {
+ double v = 0;
+ for (int j = 0; j < solution.length; j++) {
+ v += A0[i][j] * solution[j];
+ }
+ System.out.println(v + "<?>" + b0[i]);
+ }
+ }
+ }
- /**
- * Solve the phase II of the Simplex algorithm for canonical problem: no equality constraints nor negative constraints
- * @return true if optimization succeeded, false otherwise.
- * */
- public boolean solvePhaseIICanonicalSimplex(final TableauWithSlackVariables tableau, final boolean max)
- {
- if (tableau == null)
- return false;
+ /**
+ * Returns the solution of the program, once already computed.
+ *
+ * @return computed solution for the linear program
+ */
+ public double[] getSolution() {
+ return solution;
+ }
- // Gauss-Jordan simplex method
- boolean feasible = true;
- boolean reachedOptimum = false;
+ /**
+ * Returns the score of the solution of the program, once already computed.
+ *
+ * @return score for the computed solution for the linear program
+ */
+ public double getScore() {
+ if (solution != null)
+ return tableau0.scoreValue;
+ else
+ return 0;
+ }
- while(feasible && !reachedOptimum)
- {
- // start optimization
- final int colIdx = getPivotColumn(tableau, max);
- if (colIdx < 0)// optimal solution achieved
- {
- reachedOptimum = true;
- break;
- }
- else
- {
- final int rowIdx = getRowidx(tableau, colIdx);
- if (rowIdx < 0) // not feasible
- {
- feasible = false;
- break;
- }
- else
- {
- if (verbose)
- System.out.println("replacement: row = "+ rowIdx + ", column = " + colIdx);
- tableau.pivot(colIdx, rowIdx);
- }
- }
- if (verbose)
- {
- System.out.println("Tableau at iteration:");
- tableau.printTableau();
- }
- }
- if (verbose)
- {
- System.out.println("Final tableau:");
- tableau.printTableau();
- } if (!feasible)
- {
- if (verbose)
- System.out.println("UNFEASIBLE problem");
- return false;
- }
- else
- {
- solution = tableau.getSolution();
- }
- return true;
- }
+ /**
+ * Returns the number of variable for the optimization problem
+ *
+ * @return number of variables
+ */
+ public int getNumVariables() {
+ return numVariables;
+ }
+ /**
+ * Returns the number of constraints for the optimization problem
+ *
+ * @return number of constraints
+ */
+ public int getNumConstraints() {
+ return numConstraints;
+ }
- /**
- * Duplicate the program
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- protected abstract CanonicalSimplexProgram createNewProgram(double[][] A, double[] b, double[] c, boolean maximization, boolean[] equality);
- /**
- * Get the row index for pivoting
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param colIdx index of the pivoting column in the tableau
- * @return the row index in the tableau
- * */
- protected abstract int getRowidx(TableauWithSlackVariables tableau, int colIdx);
+ /**
+ * @return Parameters of the linear programming problem
+ */
+ public CanonicalProgramParameters getParameters() {
+ return parameters;
+ }
- /**
- * Get the column index for pivoting
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param max true if maximization problem, false if minization problem
- * @return the column index in the tableau
- * */
- protected abstract int getPivotColumn(TableauWithSlackVariables tableau, boolean max);
+ /**
+ * Set the parameters of the linear programming problem
+ */
+ public void setParameters(final CanonicalProgramParameters parameters) {
+ this.parameters = parameters;
+ }
- /**
- * Display whether the computed solution satisfies the constraints
- * */
- public void checkSolution()
- {
- final double[][] A0 = parameters.A;
- final double[] b0 = parameters.b;
- if (solution != null)
- {
- for (int i = 0; i < numConstraints; i++)
- {
- double v = 0;
- for (int j = 0; j < solution.length; j++)
- {
- v += A0[i][j]*solution[j];
- }
- System.out.println(v + "<?>" + b0[i]);
- }
- }
- }
- /**
- * Returns the solution of the program, once already computed.
- * @return computed solution for the linear program
- * */
- public double[] getSolution()
- {
- return solution;
- }
- /**
- * Returns the score of the solution of the program, once already computed.
- * @return score for the computed solution for the linear program
- * */
- public double getScore()
- {
- if (solution != null)
- return tableau0.scoreValue;
- else
- return 0;
- }
- /**
- * Returns the number of variable for the optimization problem
- * @return number of variables
- * */
- public int getNumVariables() {
- return numVariables;
- }
- /**
- * Returns the number of constraints for the optimization problem
- * @return number of constraints
- * */
- public int getNumConstraints()
- {
- return numConstraints;
- }
+ public void displayProblem() {
+ System.out.println("Linear Programming problem max c'x st. Ax <= b, x >= 0 and some equality contraints defined by eq");
+ if (parameters.maximization)
+ System.out.println("Maximization problem");
+ else
+ System.out.println("Minimization problem");
+ System.out.print("c = [");
+ for (int i = 0; i < parameters.c.length; i++)
+ System.out.print(parameters.c[i] + ", ");
+ System.out.println("]");
+ System.out.print("b = [");
+ for (int i = 0; i < parameters.b.length; i++)
+ System.out.print(parameters.b[i] + ", ");
+ System.out.println("]");
+ System.out.print("eq = [");
+ for (int i = 0; i < parameters.equalityConstraints.length; i++)
+ System.out.print(parameters.equalityConstraints[i] + ", ");
+ System.out.println("]");
+ System.out.print("A = [");
+ for (int i = 0; i < parameters.A.length; i++) {
+ System.out.print("[");
+ for (int j = 0; j < parameters.A[i].length; j++)
+ System.out.print(parameters.A[i][j] + ", ");
+ System.out.println("],");
+ }
+ System.out.println("]");
+ }
+ public void displayProblemTutorial() {
+ System.out.println("Linear Programming problem:");
+ System.out.println();
+ StringBuilder line = new StringBuilder();
+ if (parameters.maximization)
+ line.append("max_x");
+ else
+ line.append("min_x");
+ line.append(" ");
+ for (int i = 0; i < parameters.c.length; i++) {
+ if (parameters.c[i] >= 0 && i > 0)
+ line.append(" + ").append(parameters.c[i]).append("*x[").append(i).append("]");
+ else
+ line.append(" - ").append(Math.abs(parameters.c[i])).append("*x[").append(i).append("]");
+ }
+ System.out.println(line);
+ System.out.println();
+ System.out.println("Such that x >= 0 and:");
+ System.out.println();
+ for (int i = 0; i < parameters.b.length; i++) {
+ line = new StringBuilder();
+ for (int j = 0; j < parameters.A[i].length; j++) {
+ if (parameters.A[i][j] >= 0 && j > 0)
+ line.append(" + ").append(parameters.A[i][j]).append("*x[").append(i).append("]");
+ else
+ line.append(" - ").append(Math.abs(parameters.A[i][j])).append("*x[").append(i).append("]");
+ }
+ if (parameters.equalityConstraints[i])
+ line.append(" = ");
+ else
+ line.append(" <= ");
+ line.append(parameters.b[i]);
+ System.out.println(line);
+ }
+ }
- /**
- * @return Parameters of the linear programming problem
- * */
- public CanonicalProgramParameters getParameters() {
- return parameters;
- }
+ /**
+ * @return Manual of the library for stand-alone use through command line
+ */
+ public static String getHelp() {
+ return "Run the solver for example problems or for a user-defined system\n"
+ + "Empty arguments to use the default example.\n"
+ + "Enter an integer from 0 to 3 for different example scenarios.\n"
+ + "Enter -1 for a user-defined scenario through text files.\n "
+ + "For custom scenarios, enter filenames (text files) for different parameters of the problem preceded by the appropriate prefix:\n"
+ + "-c for the objective function file.\n"
+ + "-A for constraint matrix file.\n"
+ + "-b for the constraint value file.\n"
+ + "-e for the equality constraint file.\n"
+ + "-max or -min to indicate a maximization or minimization problem, respectively. Default in minimization.\n"
+ + "Example arguments: java -jar linearProgrammingICY.jar -c c.txt -A A.txt -b b.txt -e eq.txt -max -o solution.txt"
+ + "Each text file must contain a series of double values separated by ',' in a single line, except for the constraint file which contains one line per constraint.\n"
+ + "For the equality file '0' stands for 'false' and '1' for true.\n\n"
+ + "Version 1.0. April 2014. Author: Nicolas Chenouard. nicolas.chenouard.dev@gmail.com. Licence GPL V3.0";
+ }
- /**
- * Set the parameters of the linear programming problem
- * */
- public void setParameters(final CanonicalProgramParameters parameters) {
- this.parameters = parameters;
- }
+ /**
+ * Display the manual of the library for stand-alone use through command line
+ */
+ public static void displayHelp() {
+ System.out.println(getHelp());
+ }
- public void displayProblem()
- {
- System.out.println("Linear Programming problem max c'x st. Ax <= b, x >= 0 and some equality contraints defined by eq");
- if (parameters.maximization)
- System.out.println("Maximization problem");
- else
- System.out.println("Minimization problem");
- System.out.print("c = [");
- for (int i = 0; i < parameters.c.length; i++)
- System.out.print(parameters.c[i]+", ");
- System.out.println("]");
- System.out.print("b = [");
- for (int i = 0; i < parameters.b.length; i++)
- System.out.print(parameters.b[i]+", ");
- System.out.println("]");
- System.out.print("eq = [");
- for (int i = 0; i < parameters.equalityConstraints.length; i++)
- System.out.print(parameters.equalityConstraints[i]+", ");
- System.out.println("]");
- System.out.print("A = [");
- for (int i = 0; i < parameters.A.length; i++)
- {
- System.out.print("[");
- for (int j = 0; j < parameters.A[i].length; j++)
- System.out.print(parameters.A[i][j]+", ");
- System.out.println("],");
- }
- System.out.println("]");
- }
+ /**
+ * Run the solver for example problems or for a user-defined system
+ *
+ * @param args empty to used the default example.
+ * Enter an integer from 0 to 3 for different example scenarios.
+ * Enter -1 for a user-defined scenario through text files
+ * <p>
+ * For custom scenarios, enter filenames (text files) for different parameters of the problem preceded by the appropriate prefix:
+ * -c for the objective function file
+ * -A for constraint matrix file
+ * -b for the constraint value file
+ * -e for the equality constraint file
+ * -max or -min to indicate a maximization or minimization problem, respectively. Default in minimization.
+ * -v for verbose mode
+ * -o for an output file containing the solution
+ * Example arguments: -c c.txt -A A.txt -b b.txt -e eq.txt -max -o solution.txt
+ * Each text file must contain a series of double values separated by ',' in a single line, except for the constraint file which contains one line per constraint.
+ * For the equality file '0' stands for 'false' and '1' for true.
+ */
+ public static void main(final String[] args) {
+ int scenario = 0;
+ if (args.length > 0) {
+ if (args[0].compareTo("-h") == 0 || args[0].compareTo("-help") == 0 || args[0].compareTo("--h") == 0 || args[0].compareTo("--help") == 0) {
+ displayHelp();
+ return;
+ }
+ else {
+ try {
+ scenario = Integer.parseInt(args[0]);
+ }
+ catch (final NumberFormatException e) {
+ System.out.println("Invalid input argument");
+ System.out.println("Use input option -help to display the manual of the software");
+ IcyExceptionHandler.showErrorMessage(e, true);
+ return;
+ }
+ }
+ }
+ else {
+ System.out.println("Solving a default simple linear problem with the Simplex algorithm.");
+ System.out.println("Enter -help for the manual explicating the customization options.");
+ }
+ if (scenario < -1 || scenario > 3) {
+ System.out.println("Invalid input argument");
+ System.out.println("Use input option -help to display the manual of the software");
+ return;
+ }
+ boolean maximization = false;
+ double[] b = null;
+ double[] c = null;
+ boolean[] equality = null;
+ double[][] A = null;
+ String outputFile = null;
+ boolean verbose = true;
+ switch (scenario) {
+ case -1: {
+ verbose = false;
+ for (final String s : args) {
+ if (s.compareTo("-v") == 0) {
+ verbose = true;
+ break;
+ }
+ }
+ for (final String arg : args) {
+ if (arg.compareTo("-max") == 0) {
+ maximization = true;
+ break;
+ }
+ }
+ String fileA = null;
+ for (int i = 0; i < args.length - 1; i++) {
+ if (args[i].compareTo("-A") == 0 || args[i].compareTo("-a") == 0) {
+ fileA = args[i + 1];
+ break;
+ }
+ }
+ if (fileA == null) {
+ System.out.println("Missing -A input argument");
+ System.out.println("Use input option -help to display the manual of the software");
+ return;
+ }
+ String fileB = null;
+ for (int i = 0; i < args.length - 1; i++) {
+ if (args[i].compareTo("-b") == 0) {
+ fileB = args[i + 1];
+ break;
+ }
+ }
+ if (fileB == null) {
+ System.out.println("Missing -b input argument");
+ System.out.println("Use input option -help to display the manual of the software");
+ return;
+ }
- public void displayProblemTutorial()
- {
- System.out.println("Linear Programming problem:");
- System.out.println();
- String line = "";
- if (parameters.maximization)
- line = line + "max_x";
- else
- line = line + "min_x";
- line = line + " ";
- for (int i = 0; i < parameters.c.length; i++)
- {
- if (parameters.c[i] >= 0 && i > 0)
- line = line + " + " + parameters.c[i] + "*x[" + i+"]";
- else
- line = line + " - " + Math.abs(parameters.c[i]) + "*x[" + i+"]";
- }
- System.out.println(line);
- System.out.println();
- System.out.println("Such that x >= 0 and:");
- System.out.println();
- for (int i = 0; i < parameters.b.length; i++)
- {
- line = "";
- for (int j = 0; j < parameters.A[i].length; j++)
- {
- if (parameters.A[i][j] >= 0 && j > 0)
- line = line + " + " + parameters.A[i][j] + "*x[" + i+"]";
- else
- line = line + " - "+ Math.abs(parameters.A[i][j]) + "*x[" + i+"]";
- }
- if (parameters.equalityConstraints[i])
- line = line + " = ";
- else
- line = line + " <= ";
- line = line + parameters.b[i];
- System.out.println(line);
- }
- }
-
- /**
- * @return Manual of the library for stand-alone use through command line
- * */
- public static String getHelp()
- {
- final String helpString = "Run the solver for example problems or for a user-defined system\n"
- + "Empty arguments to use the default example.\n"
- + "Enter an integer from 0 to 3 for different example scenarios.\n"
- + "Enter -1 for a user-defined scenario through text files.\n "
- + "For custom scenarios, enter filenames (text files) for different parameters of the problem preceded by the appropriate prefix:\n"
- + "-c for the objective function file.\n"
- + "-A for constraint matrix file.\n"
- + "-b for the constraint value file.\n"
- + "-e for the equality constraint file.\n"
- + "-max or -min to indicate a maximization or minimization problem, respectively. Default in minimization.\n"
- + "Example arguments: java -jar linearProgrammingICY.jar -c c.txt -A A.txt -b b.txt -e eq.txt -max -o solution.txt"
- + "Each text file must contain a series of double values separated by ',' in a single line, except for the constraint file which contains one line per constraint.\n"
- + "For the equality file '0' stands for 'false' and '1' for true.\n\n"
- + "Version 1.0. April 2014. Author: Nicolas Chenouard. nicolas.chenouard.dev@gmail.com. Licence GPL V3.0";
- return helpString;
- }
+ String fileC = null;
+ for (int i = 0; i < args.length - 1; i++) {
+ if (args[i].compareTo("-c") == 0) {
+ fileC = args[i + 1];
+ break;
+ }
+ }
+ if (fileC == null) {
+ System.out.println("Missing -c input argument");
+ System.out.println("Use input option -help to display the manual of the software");
+ return;
+ }
+ for (int i = 0; i < args.length - 1; i++) {
+ if (args[i].compareTo("-o") == 0) {
+ outputFile = args[i + 1];
+ break;
+ }
+ }
+ String fileE = null;
+ for (int i = 0; i < args.length - 1; i++) {
+ if (args[i].compareTo("-e") == 0) {
+ fileE = args[i + 1];
+ break;
+ }
+ }
+ if (fileE == null) {
+ System.out.println("Missing -e input argument");
+ System.out.println("Use input option -help to display the manual of the software");
+ return;
+ }
+ try {
+ File file = new File(fileB);
+ FileInputStream fis;
+ fis = new FileInputStream(file);
+ BufferedReader reader = new BufferedReader(new InputStreamReader(fis));
+ String line = reader.readLine();
+ String[] tab = line.split(",");
+ b = new double[tab.length - 1];
+ for (int i = 0; i < tab.length - 1; i++)
+ b[i] = Double.parseDouble(tab[i]);
+ reader.close();
+ fis.close();
- /**
- * Display the manual of the library for stand-alone use through command line
- * */
- public static void displayHelp()
- {
- System.out.println(getHelp());
- }
-
-
- /**
- * Run the solver for example problems or for a user-defined system
- *
- * @param args empty to used the default example.
- * Enter an integer from 0 to 3 for different example scenarios.
- * Enter -1 for a user-defined scenario through text files
- *
- * For custom scenarios, enter filenames (text files) for different parameters of the problem preceded by the appropriate prefix:
- * -c for the objective function file
- * -A for constraint matrix file
- * -b for the constraint value file
- * -e for the equality constraint file
- * -max or -min to indicate a maximization or minimization problem, respectively. Default in minimization.
- * -v for verbose mode
- * -o for an output file containing the solution
- * Example arguments: -c c.txt -A A.txt -b b.txt -e eq.txt -max -o solution.txt
- * Each text file must contain a series of double values separated by ',' in a single line, except for the constraint file which contains one line per constraint.
- * For the equality file '0' stands for 'false' and '1' for true.
- *
- * */
- public static void main(final String [] args)
- {
- int scenario = 0;
- if (args.length > 0)
- {
- if(args[0].compareTo("-h") == 0 || args[0].compareTo("-help") == 0 || args[0].compareTo("--h") == 0 || args[0].compareTo("--help") == 0)
- {
- displayHelp();
- return;
- }
- else
- {
- try {
- scenario = Integer.parseInt(args[0]);
- }
- catch (final NumberFormatException e)
- {
- System.out.println("Invalid input argument");
- System.out.println("Use input option -help to display the manual of the software");
- e.printStackTrace();
- return;
- }
- }
- }
- else
- {
- System.out.println("Solving a default simple linear problem with the Simplex algorithm.");
- System.out.println("Enter -help for the manual explicating the customization options.");
- }
- if (scenario < -1 || scenario > 3)
- {
- System.out.println("Invalid input argument");
- System.out.println("Use input option -help to display the manual of the software");
- return;
- }
- boolean maximization = false;
- double[] b = null;
- double[] c = null;
- boolean[] equality = null;
- double[][] A = null;
- String outputFile = null;
- boolean verbose = true;
- switch (scenario) {
- case -1:
- {
- verbose = false;
- for (int i = 0; i < args.length; i++)
- {
- if (args[i].compareTo("-v") == 0)
- {
- verbose = true;
- break;
- }
- }
- for (int i = 0; i < args.length; i++)
- {
- if (args[i].compareTo("-max") == 0)
- {
- maximization = true;
- break;
- }
- }
- String fileA = null;
- for (int i = 0; i < args.length -1; i++)
- {
- if (args[i].compareTo("-A") == 0 || args[i].compareTo("-a") == 0)
- {
- fileA = args[i + 1];
- break;
- }
- }
- if (fileA == null)
- {
- System.out.println("Missing -A input argument");
- System.out.println("Use input option -help to display the manual of the software");
- return;
- }
- String fileB = null;
- for (int i = 0; i < args.length -1; i++)
- {
- if (args[i].compareTo("-b") == 0)
- {
- fileB = args[i + 1];
- break;
- }
- }
- if (fileB == null)
- {
- System.out.println("Missing -b input argument");
- System.out.println("Use input option -help to display the manual of the software");
- return;
- }
-
- String fileC = null;
- for (int i = 0; i < args.length -1; i++)
- {
- if (args[i].compareTo("-c") == 0)
- {
- fileC = args[i + 1];
- break;
- }
- }
- if (fileC == null)
- {
- System.out.println("Missing -c input argument");
- System.out.println("Use input option -help to display the manual of the software");
- return;
- }
- for (int i = 0; i < args.length -1; i++)
- {
- if (args[i].compareTo("-o") == 0)
- {
- outputFile = args[i + 1];
- break;
- }
- }
- String fileE = null;
- for (int i = 0; i < args.length -1; i++)
- {
- if (args[i].compareTo("-e") == 0)
- {
- fileE = args[i + 1];
- break;
- }
- }
- if (fileE == null)
- {
- System.out.println("Missing -e input argument");
- System.out.println("Use input option -help to display the manual of the software");
- return;
- }
- try {
- File file = new File(fileB);
- FileInputStream fis;
- fis = new FileInputStream(file);
- BufferedReader reader = new BufferedReader(new InputStreamReader(fis));
- String line = reader.readLine();
- String[] tab = line.split(",");
- b = new double[tab.length-1];
- for (int i = 0; i < tab.length-1; i++)
- b[i] = Double.parseDouble(tab[i]);
- reader.close();
- fis.close();
+ file = new File(fileC);
+ fis = new FileInputStream(file);
+ reader = new BufferedReader(new InputStreamReader(fis));
+ line = reader.readLine();
+ tab = line.split(",");
+ c = new double[tab.length - 1];
+ for (int i = 0; i < tab.length - 1; i++)
+ c[i] = Double.parseDouble(tab[i]);
+ reader.close();
+ fis.close();
- file = new File(fileC);
- fis = new FileInputStream(file);
- reader = new BufferedReader(new InputStreamReader(fis));
- line = reader.readLine();
- tab = line.split(",");
- c = new double[tab.length-1];
- for (int i = 0; i < tab.length-1; i++)
- c[i] = Double.parseDouble(tab[i]);
- reader.close();
- fis.close();
+ file = new File(fileE);
+ fis = new FileInputStream(file);
+ reader = new BufferedReader(new InputStreamReader(fis));
+ line = reader.readLine();
+ tab = line.split(",");
+ equality = new boolean[tab.length - 1];
+ for (int i = 0; i < tab.length - 1; i++) {
+ final double d = Double.parseDouble(tab[i]);
+ equality[i] = d > 0;
+ }
+ reader.close();
+ fis.close();
- file = new File(fileE);
- fis = new FileInputStream(file);
- reader = new BufferedReader(new InputStreamReader(fis));
- line = reader.readLine();
- tab = line.split(",");
- equality = new boolean[tab.length-1];
- for (int i = 0; i < tab.length-1; i++)
- {
- final double d = Double.parseDouble(tab[i]);
- equality[i] = d > 0;
- }
- reader.close();
- fis.close();
+ final int numRow = b.length;
+ final int numCol = c.length;
- final int numRow = b.length;
- final int numCol = c.length;
+ file = new File(fileA);
+ fis = new FileInputStream(file);
+ reader = new BufferedReader(new InputStreamReader(fis));
+ A = new double[numRow][numCol];
+ line = reader.readLine();
+ int j = 0;
+ while (line != null) {
+ tab = line.split(",");
+ for (int i = 0; i < tab.length - 1; i++) {
+ A[j][i] = Double.parseDouble(tab[i]);
+ }
+ line = reader.readLine();
+ j++;
+ }
+ reader.close();
+ fis.close();
+ }
+ catch (final IOException e) {
+ IcyExceptionHandler.showErrorMessage(e, true);
+ return;
+ }
+ break;
+ }
+ case 0: {
+ System.out.println("Basic example with 3 pivots");
+ System.out.println("Optimal is 6.5");
+ System.out.println("Solution is [1.0, 1.0, 0.5, 0.0]");
+ A = new double[][]{{2, 1, 0, 0},
+ {0, 1, 4, 1},
+ {1, 3, 0, 1}
+ };
+ b = new double[]{3, 3, 4};
+ c = new double[]{2, 4, 1, 1};
+ equality = new boolean[b.length];
+ maximization = true;
+ break;
+ }
+ case 1: {
+ System.out.println("Example that would cycle with Bland rule.");
+ System.out.println("Optimal score is 1/20.");
+ System.out.println("Solution is [1/25, 0, 1, 0].");
+ A = new double[][]{{1d / 4, -60, -1d / 25, 9},
+ {1d / 2, -90, -1d / 50, 3},
+ {0, 0, 1, 0}
+ };
+ b = new double[]{0, 0, 1};
+ c = new double[]{3d / 4, -150, 1d / 50, -6};
+ equality = new boolean[b.length];
+ maximization = true;
+ break;
+ }
+ case 2: {
+ System.out.println("An example with equality constraints");
+ System.out.println("Solution is [0, 4/7, 1 + 5/7, 0, 0].");
+ A = new double[][]
+ {{5, -4, 13, -2, 1},
+ {1, -1, 5, -1, 1},
+ };
+ b = new double[]{20, 8};
+ c = new double[]{1, 6, -7, 1, 5};
+ equality = new boolean[]{true, true};
+ maximization = false;
+ break;
+ }
+ case 3: {
+ System.out.println("An example with equality constraints and negative right-hand-side");
+ System.out.println("Optimal score is -3 - 1/16.");
+ System.out.println("Solution is [3/16, 1 + 1/4, 0, 5/16].");
+ A = new double[][]
+ {{1, 2, 1, 1},
+ {1, -2, 2, 1},
+ {3, -1, 0, -1}
+ };
+ b = new double[]{3, -2, -1};
+ c = new double[]{2, -3, 1, 1};
+ equality = new boolean[]{true, true, true};
+ maximization = true;
+ }
+ }
+ if (A == null || b == null || c == null || equality == null) {
+ System.out.println("Missing arguments");
+ System.out.println("Use input option -help to display the manual of the software");
+ return;
+ }
+ final CanonicalSimplexProgram program = new SimplexLEXICO(A, b, c, maximization, equality);
+ if (verbose) {
+ program.displayProblem();
+ }
+ if (program.solvePrimalSimplex()) {
+ final double[] solution = program.solution;
+ if (verbose) {
+ System.out.print("Computed solution = [");
+ for (final double value : solution)
+ System.out.print(value + ", ");
+ System.out.println("]");
+ System.out.println("Computed score = " + program.tableau0.scoreValue);
+ System.out.println("Computed constraint values:");
+ for (int i = 0; i < A.length; i++) {
+ double v = 0;
+ // compute the constraint
+ for (int j = 0; j < A[i].length; j++)
+ v += A[i][j] * solution[j];
+ if (equality[i])
+ System.out.println(v + " == " + b[i]);
+ else
+ System.out.println(v + " <= " + b[i]);
+ }
+ }
+ if (outputFile != null) {
+ final BufferedWriter out;
+ try {
+ out = new BufferedWriter(new FileWriter(outputFile));
+ out.write("solution\n");
+ for (final double v : solution) {
+ out.write(v + ", ");
+ }
+ out.write("\nscore\n");
+ out.write(Double.toString(program.tableau0.scoreValue));
+ out.close();
+ }
+ catch (final IOException e) {
+ IcyExceptionHandler.showErrorMessage(e, true);
+ }
+ }
+ }
+ else {
+ System.out.println("Solve simplex failed");
+ }
+ }
- file = new File(fileA);
- fis = new FileInputStream(file);
- reader = new BufferedReader(new InputStreamReader(fis));
- A = new double[numRow][numCol];
- line = reader.readLine();
- int j = 0;
- while (line!=null)
- {
- tab = line.split(",");
- for (int i = 0; i < tab.length-1; i++)
- {
- A[j][i] = Double.parseDouble(tab[i]);
- }
- line = reader.readLine();
- j++;
- }
- reader.close();
- fis.close();
- } catch (final FileNotFoundException e) {
- e.printStackTrace();
- return;
- } catch (final IOException e) {
- e.printStackTrace();
- return;
- }
- break;
- }
- case 0:
- {
- System.out.println("Basic example with 3 pivots");
- System.out.println("Optimal is 6.5");
- System.out.println("Solution is [1.0, 1.0, 0.5, 0.0]");
- A = new double[][]{{2, 1, 0, 0},
- {0, 1, 4, 1},
- {1, 3, 0, 1}
- };
- b = new double[]{3, 3, 4};
- c = new double[]{2, 4, 1, 1};
- equality = new boolean[b.length];
- maximization = true;
- break;
- }
- case 1:
- {
- System.out.println("Example that would cycle with Bland rule.");
- System.out.println("Optimal score is 1/20.");
- System.out.println("Solution is [1/25, 0, 1, 0].");
- A = new double[][]{{1d/4, -60, -1d/25, 9},
- {1d/2, -90, -1d/50, 3},
- {0, 0, 1, 0}
- };
- b = new double[]{0, 0, 1};
- c = new double[]{3d/4, -150, 1d/50, -6};
- equality = new boolean[b.length];
- maximization = true;
- break;
- }
- case 2:
- {
- System.out.println("An example with equality constraints");
- System.out.println("Solution is [0, 4/7, 1 + 5/7, 0, 0].");
- A = new double[][]
- {{5, -4, 13, -2, 1},
- {1, -1, 5, -1, 1},
- };
- b = new double[]{20, 8};
- c = new double[]{1, 6, -7, 1, 5};
- equality = new boolean[]{true, true};
- maximization = false;
- break;
- }
- case 3:
- {
- System.out.println("An example with equality constraints and negative right-hand-side");
- System.out.println("Optimal score is -3 - 1/16.");
- System.out.println("Solution is [3/16, 1 + 1/4, 0, 5/16].");
- A = new double[][]
- {{1, 2, 1, 1},
- {1, -2, 2, 1},
- {3, -1, 0, -1}
- };
- b = new double[]{3, -2, -1};
- c = new double[]{2, -3, 1, 1};
- equality = new boolean[]{true, true, true};
- maximization = true;
- }
- }
- if (A == null || b == null || c == null || equality == null)
- {
- System.out.println("Missing arguments");
- System.out.println("Use input option -help to display the manual of the software");
- return;
- }
- final CanonicalSimplexProgram program = new SimplexLEXICO(A, b, c, maximization, equality);
- if (verbose)
- {
- program.displayProblem();
- }
- if(program.solvePrimalSimplex())
- {
- final double[] solution = program.solution;
- if (verbose)
- {
- System.out.print("Computed solution = [");
- for (int i = 0; i < solution.length; i++)
- System.out.print(solution[i] + ", ");
- System.out.println("]");
- System.out.println("Computed score = " + program.tableau0.scoreValue);
- System.out.println("Computed constraint values:");
- for (int i = 0; i < A.length; i++)
- {
- double v = 0;
- // compute the constraint
- for (int j = 0; j < A[i].length; j++)
- v += A[i][j]*solution[j];
- if (equality[i])
- System.out.println(v+" == "+b[i]);
- else
- System.out.println(v+" <= "+b[i]);
- }
- }
- if (outputFile != null)
- {
- BufferedWriter out;
- try {
- out = new BufferedWriter(new FileWriter(outputFile));
- out.write("solution\n");
- for (int i = 0; i < solution.length; i++)
- {
- out.write(solution[i]+", ");
- }
- out.write("\nscore\n");
- out.write(Double.toString(program.tableau0.scoreValue));
- out.close();
- } catch (final IOException e) {
- e.printStackTrace();
- return;
- }
- }
- return;
- }
- else
- {
- System.out.println("Solve simplex failed");
- return;
- }
- }
-
- public static void runExampleScenario(final int scenario)
- {
- if (scenario < 0 || scenario > 3)
- {
- System.out.println("Invalid input argument");
- System.out.println("Scenario indices are: 0, 1, 2 and 3");
- return;
- }
- boolean maximization = false;
- double[] b = null;
- double[] c = null;
- boolean[] equality = null;
- double[][] A = null;
- final boolean verbose = true;
- switch (scenario) {
- case 0:
- {
- System.out.println("Basic example with 3 pivots");
- System.out.println("Optimal score is 6.5");
- System.out.println("Solution is [1.0, 1.0, 0.5, 0.0]");
- A = new double[][]{{2, 1, 0, 0},
- {0, 1, 4, 1},
- {1, 3, 0, 1}
- };
- b = new double[]{3, 3, 4};
- c = new double[]{2, 4, 1, 1};
- equality = new boolean[b.length];
- maximization = true;
- break;
- }
- case 1:
- {
- System.out.println("Example that would cycle with Bland rule.");
- System.out.println("Optimal score is 1/20.");
- System.out.println("Solution is [1/25, 0, 1, 0].");
- A = new double[][]{{1d/4, -60, -1d/25, 9},
- {1d/2, -90, -1d/50, 3},
- {0, 0, 1, 0}
- };
- b = new double[]{0, 0, 1};
- c = new double[]{3d/4, -150, 1d/50, -6};
- equality = new boolean[b.length];
- maximization = true;
- break;
- }
- case 2:
- {
- System.out.println("An example with equality constraints");
- System.out.println("Solution is [0, 4/7, 1 + 5/7, 0, 0].");
- A = new double[][]
- {{5, -4, 13, -2, 1},
- {1, -1, 5, -1, 1},
- };
- b = new double[]{20, 8};
- c = new double[]{1, 6, -7, 1, 5};
- equality = new boolean[]{true, true};
- maximization = false;
- break;
- }
- case 3:
- {
- System.out.println("An example with equality constraints and negative right-hand-side");
- System.out.println("Optimal score is -3 - 1/16.");
- System.out.println("Solution is [3/16, 1 + 1/4, 0, 5/16].");
- A = new double[][]
- {{1, 2, 1, 1},
- {1, -2, 2, 1},
- {3, -1, 0, -1}
- };
- b = new double[]{3, -2, -1};
- c = new double[]{2, -3, 1, 1};
- equality = new boolean[]{true, true, true};
- maximization = true;
- break;
- }
- }
- final CanonicalSimplexProgram program = new SimplexLEXICO(A, b, c, maximization, equality);
- if (verbose)
- {
- System.out.println();
- program.displayProblemTutorial();
- }
- if(program.solvePrimalSimplex())
- {
- final double[] solution = program.solution;
- if (verbose)
- {
- System.out.println();
- System.out.print("Computed solution = [");
- for (int i = 0; i < solution.length; i++)
- System.out.print(solution[i] + ", ");
- System.out.println("]");
- System.out.println();
- System.out.println("Computed score = " + program.tableau0.scoreValue);
- System.out.println();
- System.out.println("Computed constraint values:");
- for (int i = 0; i < A.length; i++)
- {
- double v = 0;
- // compute the constraint
- for (int j = 0; j < A[i].length; j++)
- v += A[i][j]*solution[j];
- if (equality[i])
- System.out.println(v+" == "+b[i]);
- else
- System.out.println(v+" <= "+b[i]);
- }
- }
- return;
- }
- else
- {
- System.out.println("Solve simplex failed");
- return;
- }
- }
+ public static void runExampleScenario(final int scenario) {
+ if (scenario < 0 || scenario > 3) {
+ System.out.println("Invalid input argument");
+ System.out.println("Scenario indices are: 0, 1, 2 and 3");
+ return;
+ }
+ boolean maximization = false;
+ double[] b = null;
+ double[] c = null;
+ boolean[] equality = null;
+ double[][] A = null;
+ final boolean verbose = true;
+ switch (scenario) {
+ case 0: {
+ System.out.println("Basic example with 3 pivots");
+ System.out.println("Optimal score is 6.5");
+ System.out.println("Solution is [1.0, 1.0, 0.5, 0.0]");
+ A = new double[][]{{2, 1, 0, 0},
+ {0, 1, 4, 1},
+ {1, 3, 0, 1}
+ };
+ b = new double[]{3, 3, 4};
+ c = new double[]{2, 4, 1, 1};
+ equality = new boolean[b.length];
+ maximization = true;
+ break;
+ }
+ case 1: {
+ System.out.println("Example that would cycle with Bland rule.");
+ System.out.println("Optimal score is 1/20.");
+ System.out.println("Solution is [1/25, 0, 1, 0].");
+ A = new double[][]{{1d / 4, -60, -1d / 25, 9},
+ {1d / 2, -90, -1d / 50, 3},
+ {0, 0, 1, 0}
+ };
+ b = new double[]{0, 0, 1};
+ c = new double[]{3d / 4, -150, 1d / 50, -6};
+ equality = new boolean[b.length];
+ maximization = true;
+ break;
+ }
+ case 2: {
+ System.out.println("An example with equality constraints");
+ System.out.println("Solution is [0, 4/7, 1 + 5/7, 0, 0].");
+ A = new double[][]
+ {{5, -4, 13, -2, 1},
+ {1, -1, 5, -1, 1},
+ };
+ b = new double[]{20, 8};
+ c = new double[]{1, 6, -7, 1, 5};
+ equality = new boolean[]{true, true};
+ maximization = false;
+ break;
+ }
+ case 3: {
+ System.out.println("An example with equality constraints and negative right-hand-side");
+ System.out.println("Optimal score is -3 - 1/16.");
+ System.out.println("Solution is [3/16, 1 + 1/4, 0, 5/16].");
+ A = new double[][]
+ {{1, 2, 1, 1},
+ {1, -2, 2, 1},
+ {3, -1, 0, -1}
+ };
+ b = new double[]{3, -2, -1};
+ c = new double[]{2, -3, 1, 1};
+ equality = new boolean[]{true, true, true};
+ maximization = true;
+ break;
+ }
+ }
+ final CanonicalSimplexProgram program = new SimplexLEXICO(A, b, c, maximization, equality);
+ if (verbose) {
+ System.out.println();
+ program.displayProblemTutorial();
+ }
+ if (program.solvePrimalSimplex()) {
+ final double[] solution = program.solution;
+ if (verbose) {
+ System.out.println();
+ System.out.print("Computed solution = [");
+ for (final double value : solution)
+ System.out.print(value + ", ");
+ System.out.println("]");
+ System.out.println();
+ System.out.println("Computed score = " + program.tableau0.scoreValue);
+ System.out.println();
+ System.out.println("Computed constraint values:");
+ for (int i = 0; i < A.length; i++) {
+ double v = 0;
+ // compute the constraint
+ for (int j = 0; j < A[i].length; j++)
+ v += A[i][j] * solution[j];
+ if (equality[i])
+ System.out.println(v + " == " + b[i]);
+ else
+ System.out.println(v + " <= " + b[i]);
+ }
+ }
+ }
+ else {
+ System.out.println("Solve simplex failed");
+ }
+ }
}
\ No newline at end of file
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/LinearProgrammingICYPlugin.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/LinearProgrammingICYPlugin.java
index 48969e6fb948b0f5eee1dca9b95e315233e90bdb..728b638cc0ec9a6309b0c23985b02aca5b14b130 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/LinearProgrammingICYPlugin.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/LinearProgrammingICYPlugin.java
@@ -1,23 +1,37 @@
-package plugins.nchenouard.linearprogrammingfullsimplex;
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+package plugins.nchenouard.linearprogrammingfullsimplex;
import icy.gui.frame.progress.AnnounceFrame;
import icy.plugin.abstract_.PluginActionable;
-public class LinearProgrammingICYPlugin extends PluginActionable
-{
- @Override
- public void run()
- {
- new AnnounceFrame("It is now running a series of example optimization problems. See the console for the output.");
- new AnnounceFrame("This is a utility plugin, intended to be used by other ICY plugins.");
- for (int i = 0; i < 3; i++)
- {
- System.out.println();
- System.out.println("=== Linear Programming test scenario "+i+" ===");
- System.out.println();
- CanonicalSimplexProgram.runExampleScenario(i);
- }
- }
-
+public class LinearProgrammingICYPlugin extends PluginActionable {
+ @Override
+ public void run() {
+ new AnnounceFrame("It is now running a series of example optimization problems. See the console for the output.");
+ new AnnounceFrame("This is a utility plugin, intended to be used by other ICY plugins.");
+ for (int i = 0; i < 3; i++) {
+ System.out.println();
+ System.out.println("=== Linear Programming test scenario " + i + " ===");
+ System.out.println();
+ CanonicalSimplexProgram.runExampleScenario(i);
+ }
+ }
+
}
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexBLAND.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexBLAND.java
index 706a0d531e43b98d76f9ee70ab20c9df91cc316c..a92558acbc8a79bb5ef7742d619125b645ef3813 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexBLAND.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexBLAND.java
@@ -1,126 +1,128 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
/**
* Implements the Bland rule for pivoting for the Simplex algorithm.
- *
+ * <p>
* Warning: this rule shows high probability of cycling.
- *
+ * <p>
* Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
-public class SimplexBLAND extends CanonicalSimplexProgram
-{
- /**
- * Constructor specifying parameters of the linear programming problem.
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- public SimplexBLAND(final double[][] A, final double[] b, final double[] c,
- final boolean maximization, final boolean[] equality) {
- super(A, b, c, maximization, equality);
- }
- /**
- * Constructor specifying parameters of the linear programming problem.
- * @param p Parameters of the linear programming problem.
- */
- public SimplexBLAND(final CanonicalProgramParameters p) {
- super(p);
- }
- /**
- * Get the row index for pivoting
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param colIdx index of the pivoting column in the tableau
- * @return the row index in the tableau
- * */
- @Override
- protected int getRowidx(final TableauWithSlackVariables tableau, final int colIdx) {
- boolean found = false;
- double maxVal = -1;
- int rowIdx = -1;
- for (int i = 0; i < tableau.xCol.length; i++)
- {
- final double val = tableau.tableau[i][colIdx];
- if (val > tol )
- {
- if (!found)
- {
- found = true;
- maxVal = tableau.xCol[i]/val;
- rowIdx = i;
- }
- else
- {
- if (tableau.xCol[i]/val < maxVal)
- {
- maxVal = tableau.xCol[i]/val;
- rowIdx = i;
- }
- }
- }
- }
- return rowIdx;
- }
- /**
- * Get the column index for pivoting using the Bland rule
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param max true if maximization problem, false if minization problem
- * @return the column index in the tableau
- * */
- @Override
- protected int getPivotColumn(final TableauWithSlackVariables tableau, final boolean max) {
- if (max)
- {
- int idx = -1;
- for (int i = 0; i <tableau.scoreRow.length; i++)
- if (tableau.scoreRow[i] < 0 && tableau.originalColOrder[i] >= 0)
- {
- idx = i;
- break;
- }
- return idx;
- }
- else
- {
- int idx = -1;
- for (int i = 0; i <tableau.scoreRow.length; i++)
- if (tableau.scoreRow[i] > 0 && tableau.originalColOrder[i] >= 0)
- {
- idx = i;
- break;
- }
- return idx;
- }
- }
-
- /**
- * Duplicate the program
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- @Override
- protected CanonicalSimplexProgram createNewProgram(final double[][] A,
- final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
- return new SimplexBLAND(A, b, c, maximization, equality);
- }
+public class SimplexBLAND extends CanonicalSimplexProgram {
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ public SimplexBLAND(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ super(A, b, c, maximization, equality);
+ }
+
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param p Parameters of the linear programming problem.
+ */
+ public SimplexBLAND(final CanonicalProgramParameters p) {
+ super(p);
+ }
+
+ /**
+ * Get the row index for pivoting
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param colIdx index of the pivoting column in the tableau
+ * @return the row index in the tableau
+ */
+ @Override
+ protected int getRowidx(final TableauWithSlackVariables tableau, final int colIdx) {
+ boolean found = false;
+ double maxVal = -1;
+ int rowIdx = -1;
+ for (int i = 0; i < tableau.xCol.length; i++) {
+ final double val = tableau.tableau[i][colIdx];
+ if (val > tol) {
+ if (!found) {
+ found = true;
+ maxVal = tableau.xCol[i] / val;
+ rowIdx = i;
+ }
+ else {
+ if (tableau.xCol[i] / val < maxVal) {
+ maxVal = tableau.xCol[i] / val;
+ rowIdx = i;
+ }
+ }
+ }
+ }
+ return rowIdx;
+ }
+
+ /**
+ * Get the column index for pivoting using the Bland rule
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param max true if maximization problem, false if minization problem
+ * @return the column index in the tableau
+ */
+ @Override
+ protected int getPivotColumn(final TableauWithSlackVariables tableau, final boolean max) {
+ if (max) {
+ int idx = -1;
+ for (int i = 0; i < tableau.scoreRow.length; i++)
+ if (tableau.scoreRow[i] < 0 && tableau.originalColOrder[i] >= 0) {
+ idx = i;
+ break;
+ }
+ return idx;
+ }
+ else {
+ int idx = -1;
+ for (int i = 0; i < tableau.scoreRow.length; i++)
+ if (tableau.scoreRow[i] > 0 && tableau.originalColOrder[i] >= 0) {
+ idx = i;
+ break;
+ }
+ return idx;
+ }
+ }
+
+ /**
+ * Duplicate the program
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ @Override
+ protected CanonicalSimplexProgram createNewProgram(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ return new SimplexBLAND(A, b, c, maximization, equality);
+ }
}
\ No newline at end of file
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexLEXICO.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexLEXICO.java
index 690331f06554ab781991b24688b7322c5c72627e..172848cf4ed0a44182c754f0e4ddb80f91270506 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexLEXICO.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexLEXICO.java
@@ -1,180 +1,172 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
import java.util.ArrayList;
/**
* Implements the lexicographic rule for pivoting for the Simplex algorithm.
- *
+ * <p>
* Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
-public class SimplexLEXICO extends CanonicalSimplexProgram
-{
- /**
- * Constructor specifying parameters of the linear programming problem.
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- public SimplexLEXICO(final double[][] A, final double[] b, final double[] c,
- final boolean maximization, final boolean[] equality) {
- super(A, b, c, maximization, equality);
- }
- /**
- * Constructor specifying parameters of the linear programming problem.
- * @param p Parameters of the linear programming problem.
- */
- public SimplexLEXICO(final CanonicalProgramParameters p) {
- super(p);
- }
+public class SimplexLEXICO extends CanonicalSimplexProgram {
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ public SimplexLEXICO(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ super(A, b, c, maximization, equality);
+ }
+
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param p Parameters of the linear programming problem.
+ */
+ public SimplexLEXICO(final CanonicalProgramParameters p) {
+ super(p);
+ }
+
+ /**
+ * Get the column index for pivoting using the lexicographic order rule
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param max true if maximization problem, false if minization problem
+ * @return the column index in the tableau
+ */
+ @Override
+ protected int getPivotColumn(final TableauWithSlackVariables tableau, final boolean max) {
+ // max score choice strategy
+ if (max) {
+ int idx = -1;
+ double score = 0;
+ for (int i = 0; i < tableau.scoreRow.length; i++)
+ if (tableau.scoreRow[i] < score) {
+ if (tableau.originalColOrder[i] >= 0 || !parameters.equalityConstraints[i]) {
+ score = tableau.scoreRow[i];
+ idx = i;
+ }
+ }
+ if (idx >= 0)
+ return idx;
+ else
+ return -1;
+ }
+ else {
+ int idx = -1;
+ double score = 0;
+ for (int i = 0; i < tableau.scoreRow.length; i++)
+ if (tableau.scoreRow[i] > score)
+ if (tableau.originalColOrder[i] >= 0 || !parameters.equalityConstraints[i]) {
+ score = tableau.scoreRow[i];
+ idx = i;
+ }
+ if (idx >= 0)
+ return idx;
+ else
+ return -1;
+ }
+ }
+
+ /**
+ * Get the row index for pivoting
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param colIdx index of the pivoting column in the tableau
+ * @return the row index in the tableau
+ */
+ @Override
+ protected int getRowidx(final TableauWithSlackVariables tableau, final int colIdx) {
+ int rowIdx = -1;
+ int numBasicRow = 0;
+ final ArrayList<Integer> basicCols = new ArrayList<>();
+ for (int i = 0; i < tableau.originalColOrder.length; i++)
+ if (tableau.originalColOrder[i] < 0) {
+ basicCols.add(i);
+ numBasicRow++;
+ }
+ final double[] minLexico = new double[numBasicRow + 1];
+ boolean found = false;
+ for (int i = 0; i < tableau.xCol.length; i++) {
+ final double val = tableau.tableau[i][colIdx];
+ if (val > tol) {
+ int k = 0;
+ boolean isMinLexico = false;
+ if (!found) {
+ found = true;
+ isMinLexico = true;
+ }
+ else {
+ if (tableau.xCol[i] / val == minLexico[k]) {
+ k++;
+ for (final Integer ii : basicCols) {
+ final double t = tableau.tableau[i][ii] / val;
+ if (t > minLexico[k]) {
+ isMinLexico = false;
+ break;
+ }
+ else if (t < minLexico[k]) {
+ isMinLexico = true;
+ break;
+ }
+ k++;
+ }
+ }
+ else if (tableau.xCol[i] / val < minLexico[k])
+ isMinLexico = true;
+ }
+ if (isMinLexico) {
+ if (k == 0) {
+ minLexico[k] = tableau.xCol[i] / val;
+ k++;
+ }
+ while (k <= basicCols.size()) {
+ minLexico[k] = tableau.tableau[i][basicCols.get(k - 1)] / val;
+ k++;
+ }
+ rowIdx = i;
+ }
+ }
+ }
+ return rowIdx;
+ }
- /**
- * Get the column index for pivoting using the lexicographic order rule
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param max true if maximization problem, false if minization problem
- * @return the column index in the tableau
- * */
- @Override
- protected int getPivotColumn(final TableauWithSlackVariables tableau, final boolean max) {
- // max score choice strategy
- if (max)
- {
- int idx = -1;
- double score = 0;
- for (int i = 0; i <tableau.scoreRow.length; i++)
- if (tableau.scoreRow[i] < score)
- {
- if (tableau.originalColOrder[i] >= 0 || !parameters.equalityConstraints[i])
- {
- score = tableau.scoreRow[i];
- idx = i;
- }
- }
- if (idx >= 0)
- return idx;
- else
- return -1;
- }
- else
- {
- int idx = -1;
- double score = 0;
- for (int i = 0; i <tableau.scoreRow.length; i++)
- if (tableau.scoreRow[i] > score)
- if (tableau.originalColOrder[i] >= 0 || !parameters.equalityConstraints[i])
- {
- score = tableau.scoreRow[i];
- idx = i;
- }
- if (idx >= 0)
- return idx;
- else
- return -1;
- }
- }
- /**
- * Get the row index for pivoting
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param colIdx index of the pivoting column in the tableau
- * @return the row index in the tableau
- * */
- @Override
- protected int getRowidx(final TableauWithSlackVariables tableau, final int colIdx)
- {
- int rowIdx = -1;
- int numBasicRow = 0;
- final ArrayList<Integer> basicCols = new ArrayList<Integer>();
- for (int i = 0; i < tableau.originalColOrder.length; i++)
- if (tableau.originalColOrder[i] < 0)
- {
- basicCols.add(i);
- numBasicRow ++;
- }
- final double[] minLexico = new double[numBasicRow + 1];
- boolean found = false;
- for (int i = 0; i < tableau.xCol.length; i++)
- {
- final double val = tableau.tableau[i][colIdx];
- if (val > tol)
- {
- int k = 0;
- boolean isMinLexico = false;
- if (!found)
- {
- found = true;
- isMinLexico = true;
- }
- else
- {
- if (tableau.xCol[i]/val == minLexico[k])
- {
- k++;
- for (final Integer ii:basicCols)
- {
- final double t = tableau.tableau[i][ii]/val;
- if (t > minLexico[k])
- {
- isMinLexico = false;
- break;
- }
- else if (t < minLexico[k])
- {
- isMinLexico = true;
- break;
- }
- k++;
- }
- }
- else if (tableau.xCol[i]/val < minLexico[k])
- isMinLexico = true;
- }
- if (isMinLexico)
- {
- if (k == 0)
- {
- minLexico[k] = tableau.xCol[i]/val;
- k++;
- }
- while (k <= basicCols.size() )
- {
- minLexico[k] = tableau.tableau[i][basicCols.get(k - 1)]/val;
- k++;
- }
- rowIdx = i;
- }
- }
- }
- return rowIdx;
- }
-
- /**
- * Duplicate the program
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- @Override
- protected CanonicalSimplexProgram createNewProgram(final double[][] A,
- final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
- return new SimplexLEXICO(A, b, c, maximization, equality);
- }
+ /**
+ * Duplicate the program
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ @Override
+ protected CanonicalSimplexProgram createNewProgram(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ return new SimplexLEXICO(A, b, c, maximization, equality);
+ }
}
\ No newline at end of file
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexMAX.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexMAX.java
index 705ddab138910d7f4802c08a59432769f9030aae..a65f066db8104b58b7a399262f754a9ed97d6d6e 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexMAX.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/SimplexMAX.java
@@ -1,138 +1,136 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
/**
* Implements the maximum score rule for pivoting for the Simplex algorithm.
- *
+ * <p>
* Warning: this rule shows high probability of cycling.
- *
+ * <p>
* Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
-public class SimplexMAX extends CanonicalSimplexProgram
-{
- /**
- * Constructor specifying parameters of the linear programming problem.
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- public SimplexMAX(final double[][] A, final double[] b, final double[] c,
- final boolean maximization, final boolean[] equality) {
- super(A, b, c, maximization, equality);
- }
- /**
- * Constructor specifying parameters of the linear programming problem.
- * @param p Parameters of the linear programming problem.
- */
- public SimplexMAX(final CanonicalProgramParameters p)
- {
- super(p);
- }
+public class SimplexMAX extends CanonicalSimplexProgram {
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ public SimplexMAX(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ super(A, b, c, maximization, equality);
+ }
+
+ /**
+ * Constructor specifying parameters of the linear programming problem.
+ *
+ * @param p Parameters of the linear programming problem.
+ */
+ public SimplexMAX(final CanonicalProgramParameters p) {
+ super(p);
+ }
+
+ /**
+ * Get the column index for pivoting using the maximum score rule
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param max true if maximization problem, false if minization problem
+ * @return the column index in the tableau
+ */
+ @Override
+ protected int getPivotColumn(final TableauWithSlackVariables tableau, final boolean max) {
+ if (max) {
+ int idx = -1;
+ double score = 0;
+ for (int i = 0; i < tableau.scoreRow.length; i++)
+ if (tableau.scoreRow[i] < score && tableau.originalColOrder[i] >= 0) {
+ score = tableau.scoreRow[i];
+ idx = i;
+ }
+ if (idx >= 0)
+ return idx;
+ else
+ return -1;
+ }
+ else {
+ int idx = -1;
+ double score = 0;
+ for (int i = 0; i < tableau.scoreRow.length; i++)
+ if (tableau.scoreRow[i] > score && tableau.originalColOrder[i] >= 0) {
+ score = tableau.scoreRow[i];
+ idx = i;
+ }
+ if (idx >= 0)
+ return idx;
+ else
+ return -1;
+ }
+ }
+
+ /**
+ * Get the row index for pivoting
+ *
+ * @param tableau the tableau for Gauss-Jordan elimination
+ * @param colIdx index of the pivoting column in the tableau
+ * @return the row index in the tableau
+ */
+ @Override
+ protected int getRowidx(final TableauWithSlackVariables tableau, final int colIdx) {
+ boolean found = false;
+ double maxVal = -1;
+ int rowIdx = -1;
+ for (int i = 0; i < tableau.xCol.length; i++) {
+ final double val = tableau.tableau[i][colIdx];
+ if (val > tol) {
+ if (!found) {
+ found = true;
+ maxVal = tableau.xCol[i] / val;
+ rowIdx = i;
+ }
+ else {
+ if (tableau.xCol[i] / val < maxVal) {
+ maxVal = tableau.xCol[i] / val;
+ rowIdx = i;
+ }
+ }
+ }
+ }
+ return rowIdx;
+ }
- /**
- * Get the column index for pivoting using the maximum score rule
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param max true if maximization problem, false if minization problem
- * @return the column index in the tableau
- * */
- @Override
- protected int getPivotColumn(final TableauWithSlackVariables tableau, final boolean max)
- {
- if (max)
- {
- int idx = -1;
- double score = 0;
- for (int i = 0; i <tableau.scoreRow.length; i++)
- if (tableau.scoreRow[i] < score && tableau.originalColOrder[i] >= 0)
- {
- score = tableau.scoreRow[i];
- idx = i;
- }
- if (idx >= 0)
- return idx;
- else
- return -1;
- }
- else
- {
- int idx = -1;
- double score = 0;
- for (int i = 0; i <tableau.scoreRow.length; i++)
- if (tableau.scoreRow[i] > score && tableau.originalColOrder[i] >= 0)
- {
- score = tableau.scoreRow[i];
- idx = i;
- }
- if (idx >= 0)
- return idx;
- else
- return -1;
- }
- }
- /**
- * Get the row index for pivoting
- * @param tableau the tableau for Gauss-Jordan elimination
- * @param colIdx index of the pivoting column in the tableau
- * @return the row index in the tableau
- * */
- @Override
- protected int getRowidx(final TableauWithSlackVariables tableau, final int colIdx)
- {
- boolean found = false;
- double maxVal = -1;
- int rowIdx = -1;
- for (int i = 0; i < tableau.xCol.length; i++)
- {
- final double val = tableau.tableau[i][colIdx];
- if (val > tol)
- {
- if (!found)
- {
- found = true;
- maxVal = tableau.xCol[i]/val;
- rowIdx = i;
- }
- else
- {
- if (tableau.xCol[i]/val < maxVal)
- {
- maxVal = tableau.xCol[i]/val;
- rowIdx = i;
- }
- }
- }
- }
- return rowIdx;
- }
-
- /**
- * Duplicate the program
- *
- * @param A
- * constraint matrix
- * @param b
- * constraint values
- * @param c
- * linear objective function weights
- * @param maximization
- * maximization problem if true, minimization else.
- * @param equality
- * indicates whether each constraint is an equality. Else it is
- * <=.
- */
- @Override
- protected CanonicalSimplexProgram createNewProgram(final double[][] A,
- final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
- return new SimplexMAX(A, b, c, maximization, equality);
- }
+ /**
+ * Duplicate the program
+ *
+ * @param A constraint matrix
+ * @param b constraint values
+ * @param c linear objective function weights
+ * @param maximization maximization problem if true, minimization else.
+ * @param equality indicates whether each constraint is an equality. Else it is
+ * <=.
+ */
+ @Override
+ protected CanonicalSimplexProgram createNewProgram(final double[][] A, final double[] b, final double[] c, final boolean maximization, final boolean[] equality) {
+ return new SimplexMAX(A, b, c, maximization, equality);
+ }
}
\ No newline at end of file
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauMinSlackObjective.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauMinSlackObjective.java
index 229fc310425747f79541fbb6dfec99f370a43f03..55a2ac2a29a356aac9ebad6bc04e02803273cd77 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauMinSlackObjective.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauMinSlackObjective.java
@@ -1,88 +1,94 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
/**
* Tableau used for Gauss-Jordan elimination in which slack variables are given a 1 score depending on whether constraints are equalities.
* Slack variables are added to the set of column vectors in a leading position.
- *
- * Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ * <p>
+ * Part of the Linear Programming plugin for ICY <a href="https://icy.bioimageanalysis.org">https://icy.bioimageanalysis.org</a>
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
-public class TableauMinSlackObjective extends TableauWithSlackVariables
-{
+public class TableauMinSlackObjective extends TableauWithSlackVariables {
+
+ /**
+ * Create the tableau based on a linear program. The scores in the tableau
+ * will not be initialized to their standard value but set to handle
+ * equality constraints.
+ *
+ * @param simplexProgram parameters defining the linear optimization problem.
+ */
+ public TableauMinSlackObjective(final CanonicalSimplexProgram simplexProgram) {
+ super(simplexProgram.parameters.A);
+ // now make the slacks feasible
+ xCol = simplexProgram.parameters.b.clone(); // equal to b since we start with unit vectors as initial left vectors
+ // build the score row
+ scoreRow = new double[numCol];
+ // add a score for slacks corresponding to equality constraints
+ final double[] slackScores = new double[simplexProgram.parameters.equalityConstraints.length];
+ this.scoreValue = 0;
+ for (int i = 0; i < slackScores.length; i++)
+ if (simplexProgram.parameters.equalityConstraints[i]) {
+ scoreRow[i] = 1;
+ scoreValue += xCol[i];
+ }
+ for (int i = 0; i < numCol; i++) {
+ if (originalColOrder[i] >= 0) {
+ double pi = 0;
+ for (int j = 0; j < numRows; j++) {
+ if (simplexProgram.parameters.equalityConstraints[j])
+ pi += tableau[j][i];
+ }
+ scoreRow[i] = pi;
+ }
+ }
+ }
- /**
- * Create the tableau based on a linear program. The scores in the tableau
- * will not be initialized to their standard value but set to handle
- * equality constraints.
- *
- * @param simplexProgram
- * parameters defining the linear optimization problem.
- */
- public TableauMinSlackObjective(final CanonicalSimplexProgram simplexProgram) {
- super(simplexProgram.parameters.A);
- // now make the slacks feasible
- xCol = simplexProgram.parameters.b.clone(); // equal to b since we start with unit vectors as initial left vectors
- // build the score row
- scoreRow = new double[numCol];
- // add a score for slacks corresponding to equality constraints
- final double[] slackScores = new double[simplexProgram.parameters.equalityConstraints.length];
- this.scoreValue = 0;
- for (int i = 0; i < slackScores.length; i++)
- if (simplexProgram.parameters.equalityConstraints[i])
- {
- scoreRow[i] = 1;
- scoreValue += xCol[i];
- }
- for (int i = 0; i < numCol; i++)
- {
- if (originalColOrder[i] >= 0)
- {
- double pi = 0;
- for (int j = 0; j < numRows; j++)
- {
- if (simplexProgram.parameters.equalityConstraints[j])
- pi += tableau[j][i];
- }
- scoreRow[i] = pi;
- }
- }
- }
-
- /**
- * Initialize scores for the tableau.
- * The scores are not set to their standard value.
- * Instead score 1 is given to slack variables corresponding to equality constraints and 0 to others.
- *
- * @param c standard score for variables.
- * @param slackScore score values for slack variables.
- */
- public void initScores(final double[] c, final double[] slackScore)
- {
- final double[] scoreValues = new double[numCol];
- for (int i = 0; i < slackScore.length; i++)
- scoreValues[i] = slackScore[i];
- for (int i = 0; i < numCol; i++)
- {
- double pi = 0;
- for (int j = 0; j < numRows; j++)
- {
- if (leftVectorsIdx[j] >= 0)
- pi += tableau[j][i]*scoreValues[leftVectorsIdx[j]];
- // else do nothing as unit vectors have no cost
- }
- if (originalColOrder[i] >= 0)
- scoreRow[i] = pi - c[originalColOrder[i]];
- else
- scoreRow[i] = pi;
- }
- scoreValue = 0;
- for (int i = 0; i < numRows; i++)
- {
- if (leftVectorsIdx[i] >= 0 && originalColOrder[leftVectorsIdx[i]] >= 0)
- {
- scoreValue += c[originalColOrder[leftVectorsIdx[i]]]*xCol[i];
- }
- }
- }
+ /**
+ * Initialize scores for the tableau.
+ * The scores are not set to their standard value.
+ * Instead score 1 is given to slack variables corresponding to equality constraints and 0 to others.
+ *
+ * @param c standard score for variables.
+ * @param slackScore score values for slack variables.
+ */
+ public void initScores(final double[] c, final double[] slackScore) {
+ final double[] scoreValues = new double[numCol];
+ System.arraycopy(slackScore, 0, scoreValues, 0, slackScore.length);
+ for (int i = 0; i < numCol; i++) {
+ double pi = 0;
+ for (int j = 0; j < numRows; j++) {
+ if (leftVectorsIdx[j] >= 0)
+ pi += tableau[j][i] * scoreValues[leftVectorsIdx[j]];
+ // else do nothing as unit vectors have no cost
+ }
+ if (originalColOrder[i] >= 0)
+ scoreRow[i] = pi - c[originalColOrder[i]];
+ else
+ scoreRow[i] = pi;
+ }
+ scoreValue = 0;
+ for (int i = 0; i < numRows; i++) {
+ if (leftVectorsIdx[i] >= 0 && originalColOrder[leftVectorsIdx[i]] >= 0) {
+ scoreValue += c[originalColOrder[leftVectorsIdx[i]]] * xCol[i];
+ }
+ }
+ }
}
diff --git a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauWithSlackVariables.java b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauWithSlackVariables.java
index bb5799272aef235bc45b337a000f0325bf79d295..110f17d2a770324cb1ce761b9a7555e8f24817b6 100644
--- a/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauWithSlackVariables.java
+++ b/src/main/java/plugins/nchenouard/linearprogrammingfullsimplex/TableauWithSlackVariables.java
@@ -1,342 +1,336 @@
+/*
+ * Copyright (c) 2010-2023. Institut Pasteur.
+ *
+ * This file is part of Icy.
+ * Icy is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * Icy is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with Icy. If not, see <https://www.gnu.org/licenses/>.
+ */
+
package plugins.nchenouard.linearprogrammingfullsimplex;
/**
* Tableau used for Gauss-Jordan elimination.
* Initial basis is built by using slack variables.
* Slack variables are added to the set of column vectors in a leading position.
- *
- * Part of the Linear Programming plugin for ICY http://icy.bioimageanalysis.org
- *
+ * <p>
+ * Part of the Linear Programming plugin for ICY <a href="https://icy.bioimageanalysis.org">https://icy.bioimageanalysis.org</a>
+ *
* @author Nicolas Chenouard (nicolas.chenouard.dev@gmail.com)
*/
public class TableauWithSlackVariables {
- /**
- * Set of column vectors for the tableau. The first numVariables ones are unit vectors corresponding to slack variables.
- * */
- double[][] columnVectors;
- /**
- * Basis vectors for the tableau. In the initial state they are unit vectors corresponding to slack variables.
- * */
- double[][] basisVectors;
- /**
- * Indices of basis vectors in the set of column vectors columnVectors.
- * */
- int[] leftVectorsIdx;
- /**
- * Index of column vectors in the original constraint matrix. Unit vectors added for slack variables have index -1.
- * */
- int[] originalColOrder;
- /**
- * Row of the tableau corresponding to scores
- * */
- double[] scoreRow;
- /**
- * Number of columns of the tableau
- * */
- int numCol;
- /**
- * Number of rows of the tableau
- * */
- int numRows;
- /**
- * Tableau for Gauss-Jordan elimination.
- * */
- double[][] tableau;
- /**
- * Solution column. First rows correspond to slack variables.
- * */
- double[] xCol;
- /**
- * Number of variables
- * */
- int numVariables;
- /**
- * Score value for the current solution
- * */
- double scoreValue;
+ /**
+ * Set of column vectors for the tableau. The first numVariables ones are unit vectors corresponding to slack variables.
+ */
+ double[][] columnVectors;
+ /**
+ * Basis vectors for the tableau. In the initial state they are unit vectors corresponding to slack variables.
+ */
+ double[][] basisVectors;
+ /**
+ * Indices of basis vectors in the set of column vectors columnVectors.
+ */
+ int[] leftVectorsIdx;
+ /**
+ * Index of column vectors in the original constraint matrix. Unit vectors added for slack variables have index -1.
+ */
+ int[] originalColOrder;
+ /**
+ * Row of the tableau corresponding to scores
+ */
+ double[] scoreRow;
+ /**
+ * Number of columns of the tableau
+ */
+ int numCol;
+ /**
+ * Number of rows of the tableau
+ */
+ int numRows;
+ /**
+ * Tableau for Gauss-Jordan elimination.
+ */
+ double[][] tableau;
+ /**
+ * Solution column. First rows correspond to slack variables.
+ */
+ double[] xCol;
+ /**
+ * Number of variables
+ */
+ int numVariables;
+ /**
+ * Score value for the current solution
+ */
+ double scoreValue;
+
+ /**
+ * Create the tableau based on a constraint matrix.
+ * Each element A[i] defines a linear combination of the variable corresponding to a constraint.
+ *
+ * @param A table of linear combinations of the variable, each defining a constraint.
+ */
+ public TableauWithSlackVariables(final double[][] A) {
+ createVectors(A);
+ }
+
+ /**
+ * Create the tableau based on a linear program.
+ *
+ * @param parameters parameters defining the linear optimization problem.
+ */
+ public TableauWithSlackVariables(final CanonicalProgramParameters parameters) {
+ createVectors(parameters.A); // constraint value column
+ xCol = parameters.b.clone(); // equal to b since we start with unit vectors as initial left vectors
+ // build the score row
+ scoreRow = new double[numCol];
+ modifyScores(parameters.c);
+ }
+
+ /**
+ * Change the score of slack variables to 1
+ *
+ * @param changeSlack indicate whether the score of each slack variable should be changed
+ */
+ protected void setUnitScoreToSlackVars(final boolean[] changeSlack) {
+ final double[] scoreValues = new double[numCol];
+ // recompute scoreValues based on c
+ int k = 0;
+ for (int i = 0; i < originalColOrder.length; i++) {
+ if (originalColOrder[i] < 0) {
+ if (changeSlack[k])
+ scoreValues[i] = 1d;
+ k++;
+ }
+ }
+ scoreRow = new double[numCol];
+ for (int i = 0; i < numCol; i++) {
+ double pi = 0;
+ for (int j = 0; j < numRows; j++) {
+ if (leftVectorsIdx[j] >= 0)
+ pi += tableau[j][i] * scoreValues[leftVectorsIdx[j]];
+ }
+ scoreRow[i] = pi - scoreValues[i];
+ }
+ scoreValue = 0;
+ for (int i = 0; i < leftVectorsIdx.length; i++) {
+ scoreValue += scoreValues[leftVectorsIdx[i]] * xCol[i];
+ }
+
+ k = 0;
+ for (int i = 0; i < originalColOrder.length; i++) {
+ if (originalColOrder[i] < 0) {
+ if (changeSlack[k])
+ scoreRow[i] = 1;
+ k++;
+ }
+ }
+ }
+
+ /**
+ * Create column and basis vectors based on the constraint matrix
+ *
+ * @param A Constraint matrix. Each element A[i] defines a linear combination of the variable corresponding to a constraint.
+ */
+ protected void createVectors(final double[][] A) {
+ final int numConstraints = A.length;
+ numVariables = A[0].length;
+ numCol = numConstraints + numVariables;
+ numRows = numConstraints;
+ tableau = new double[numRows][numCol];
+
+ columnVectors = new double[numConstraints + numVariables][numConstraints];
+ originalColOrder = new int[numConstraints + numVariables];
+ basisVectors = new double[numConstraints][];
+ leftVectorsIdx = new int[numConstraints];
+ // set unit vectors as left vectors of the tableau, and add them to the column vector list
+ for (int i = 0; i < numConstraints; i++) {
+ columnVectors[i][i] = 1d;
+ tableau[i][i] = 1d;
+ basisVectors[i] = columnVectors[i];
+ originalColOrder[i] = -1;
+ leftVectorsIdx[i] = i;
+ }
+ // build the set of columns corresponding to initial variables
+ for (int i = 0; i < numVariables; i++) {
+ for (int j = 0; j < numConstraints; j++)
+ columnVectors[i + numConstraints][j] = A[j][i];
+ originalColOrder[i + numConstraints] = i;
+ }
+ for (int j = 0; j < numConstraints; j++)
+ System.arraycopy(A[j], 0, tableau[j], numConstraints, numVariables);
+ }
- /**
- * Create the tableau based on a constraint matrix.
- * Each element A[i] defines a linear combination of the variable corresponding to a constraint.
- *
- * @param A table of linear combinations of the variable, each defining a constraint.
- */
- public TableauWithSlackVariables(final double[][] A)
- {
- createVectors(A);
- }
-
- /**
- * Create the tableau based on a linear program.
- *
- * @param parameters parameters defining the linear optimization problem.
- */
- public TableauWithSlackVariables(final CanonicalProgramParameters parameters)
- {
- createVectors(parameters.A); // constraint value column
- xCol = parameters.b.clone(); // equal to b since we start with unit vectors as initial left vectors
- // build the score row
- scoreRow = new double[numCol];
- modifyScores(parameters.c);
- }
- /**
- * Change the score of slack variables to 1
- *
- * @param changeSlack indicate whether the score of each slack variable should be changed
- */
- protected void setUnitScoreToSlackVars(final boolean[] changeSlack)
- {
- final double[] scoreValues = new double[numCol];
- // recompute scoreValues based on c
- int k = 0;
- for (int i = 0; i < originalColOrder.length; i++)
- {
- if (originalColOrder[i] < 0)
- {
- if (changeSlack[k])
- scoreValues[i] = 1d;
- k++;
- }
- }
- scoreRow = new double[numCol];
- for (int i = 0; i < numCol; i++)
- {
- double pi = 0;
- for (int j = 0; j < numRows; j++)
- {
- if (leftVectorsIdx[j] >= 0)
- pi += tableau[j][i]*scoreValues[leftVectorsIdx[j]];
- }
- scoreRow[i] = pi - scoreValues[i];
- }
- scoreValue = 0;
- for (int i = 0; i < leftVectorsIdx.length; i++)
- {
- scoreValue += scoreValues[leftVectorsIdx[i]]*xCol[i];
- }
+ /**
+ * Print the tableau
+ */
+ public void printTableau() {
+ System.out.println("Column vectors");
+ for (int i = 0; i < columnVectors.length; i++) {
+ final double[] vector = columnVectors[i];
+ System.out.print("a" + i + " = [");
+ for (final double v : vector)
+ System.out.print(v + ", ");
+ System.out.println("]'");
+ }
+ System.out.println("Basis vectors");
+ for (int i = 0; i < leftVectorsIdx.length; i++) {
+ final double[] vector = columnVectors[leftVectorsIdx[i]];
+ System.out.print("v" + i + " = [");
+ for (final double v : vector)
+ System.out.print(v + ", ");
+ System.out.println("]'");
+ }
+ System.out.println("Tableau");
+ for (final double[] doubles : tableau) {
+ for (final double aDouble : doubles)
+ System.out.print(aDouble + " ");
+ System.out.println();
+ }
- k = 0;
- for (int i = 0; i < originalColOrder.length; i++)
- {
- if (originalColOrder[i] < 0)
- {
- if (changeSlack[k])
- scoreRow[i] = 1;
- k++;
- }
- }
- }
- /**
- * Create column and basis vectors based on the constraint matrix
- * @param A Constraint matrix. Each element A[i] defines a linear combination of the variable corresponding to a constraint.
- * */
- protected void createVectors(final double[][] A)
- {
- final int numConstraints = A.length;
- numVariables = A[0].length;
- numCol = numConstraints + numVariables;
- numRows = numConstraints;
- tableau = new double[numRows][numCol];
+ System.out.println();
+ System.out.print("x = [");
+ for (final double v : xCol)
+ System.out.print(v + ", ");
+ System.out.println("]'");
- columnVectors = new double[numConstraints + numVariables][numConstraints];
- originalColOrder = new int[numConstraints + numVariables];
- basisVectors = new double[numConstraints][];
- leftVectorsIdx = new int[numConstraints];
- // set unit vectors as left vectors of the tableau, and add them to the column vector list
- for (int i = 0; i < numConstraints; i++)
- {
- columnVectors[i][i] = 1d;
- tableau[i][i] = 1d;
- basisVectors[i] = columnVectors[i];
- originalColOrder[i] = -1;
- leftVectorsIdx[i] = i;
- }
- // build the set of columns corresponding to initial variables
- for (int i = 0; i < numVariables; i++)
- {
- for (int j = 0; j < numConstraints; j++)
- columnVectors[i + numConstraints][j] = A[j][i];
- originalColOrder[i + numConstraints] = i;
- }
- for (int j = 0; j < numConstraints; j++)
- System.arraycopy(A[j], 0, tableau[j], numConstraints, numVariables);
- }
- /**
- * Print the tableau
- * */
- public void printTableau()
- {
- System.out.println("Column vectors");
- for (int i = 0; i < columnVectors.length; i++)
- {
- final double[] vector = columnVectors[i];
- System.out.print("a"+i+" = [");
- for (int j = 0; j < vector.length; j++)
- System.out.print(vector[j]+", ");
- System.out.println("]'");
- }
- System.out.println("Basis vectors");
- for (int i = 0; i < leftVectorsIdx.length; i++)
- {
- final double[] vector = columnVectors[leftVectorsIdx[i]];
- System.out.print("v"+i+" = [");
- for (int j = 0; j < vector.length; j++)
- System.out.print(vector[j]+", ");
- System.out.println("]'");
- }
- System.out.println("Tableau");
- for (int i = 0; i < tableau.length; i++)
- {
- for (int j = 0; j <tableau[i].length; j++)
- System.out.print(tableau[i][j]+" ");
- System.out.println();
- }
+ System.out.println();
+ System.out.print("score row = [");
+ for (final double v : scoreRow)
+ System.out.print(v + ", ");
+ System.out.println("]'");
- System.out.println();
- System.out.print("x = [");
- for (int i = 0; i < xCol.length; i++)
- System.out.print(xCol[i]+", ");
- System.out.println("]'");
+ System.out.println("Score = " + scoreValue);
+ }
- System.out.println();
- System.out.print("score row = [");
- for (int i = 0; i < scoreRow.length; i++)
- System.out.print(scoreRow[i]+", ");
- System.out.println("]'");
+ /**
+ * Pivot at a given position in the tableau
+ *
+ * @param colIdx index of the pivoting column in the tableau
+ * @param rowIdx index of the pivoting row in the tableau
+ */
+ public void pivot(final int colIdx, final int rowIdx) {
+ leftVectorsIdx[rowIdx] = colIdx;
- System.out.println("Score = " + scoreValue);
- }
-
- /**
- * Pivot at a given position in the tableau
- * @param colIdx index of the pivoting column in the tableau
- * @param rowIdx index of the pivoting row in the tableau
- */
- public void pivot(final int colIdx, final int rowIdx)
- {
- leftVectorsIdx[rowIdx] = colIdx;
+ final double pivotVal = tableau[rowIdx][colIdx];
+ final double[] pivotRow = tableau[rowIdx];
+ final double[] pivotCol = new double[numRows];
+ for (int i = 0; i < numRows; i++)
+ pivotCol[i] = tableau[i][colIdx];
+ final double scorePivot = scoreRow[colIdx];
- final double pivotVal = tableau[rowIdx][colIdx];
- final double[] pivotRow = tableau[rowIdx];
- final double[] pivotCol = new double[numRows];
- for (int i = 0; i < numRows; i++)
- pivotCol[i] = tableau[i][colIdx];
- final double scorePivot = scoreRow[colIdx];
+ // update the score row
+ for (int j = 0; j < colIdx; j++)
+ scoreRow[j] = scoreRow[j] - pivotRow[j] * scorePivot / pivotVal;
+ for (int j = colIdx + 1; j < numCol; j++)
+ scoreRow[j] = scoreRow[j] - pivotRow[j] * scorePivot / pivotVal;
+ scoreRow[colIdx] = 0;
- // update the score row
- for (int j = 0; j < colIdx; j++)
- scoreRow[j] = scoreRow[j] - pivotRow[j]*scorePivot/pivotVal;
- for (int j = colIdx + 1; j < numCol; j++)
- scoreRow[j] = scoreRow[j] - pivotRow[j]*scorePivot/pivotVal;
- scoreRow[colIdx] = 0;
+ scoreValue = scoreValue - xCol[rowIdx] * scorePivot / pivotVal;
- scoreValue = scoreValue - xCol[rowIdx]*scorePivot/pivotVal;
+ // update the variable column
+ for (int i = 0; i < rowIdx; i++)
+ xCol[i] = xCol[i] - xCol[rowIdx] * pivotCol[i] / pivotVal;
+ for (int i = rowIdx + 1; i < numRows; i++)
+ xCol[i] = xCol[i] - xCol[rowIdx] * pivotCol[i] / pivotVal;
+ xCol[rowIdx] = xCol[rowIdx] / pivotVal;
- // update the variable column
- for (int i = 0; i < rowIdx; i++)
- xCol[i] = xCol[i] - xCol[rowIdx]*pivotCol[i]/pivotVal;
- for (int i = rowIdx + 1; i < numRows; i++)
- xCol[i] = xCol[i] - xCol[rowIdx]*pivotCol[i]/pivotVal;
- xCol[rowIdx] = xCol[rowIdx]/pivotVal;
+ // update the tableau
+ for (int i = 0; i < rowIdx; i++) {
+ for (int j = 0; j < colIdx; j++)
+ tableau[i][j] = tableau[i][j] - pivotCol[i] * pivotRow[j] / pivotVal; // tij - til*tkj/tkl for col l and row k
+ for (int j = colIdx + 1; j < numCol; j++)
+ tableau[i][j] = tableau[i][j] - pivotCol[i] * pivotRow[j] / pivotVal;
+ }
+ for (int i = rowIdx + 1; i < numRows; i++) {
+ for (int j = 0; j < colIdx; j++)
+ tableau[i][j] = tableau[i][j] - pivotCol[i] * pivotRow[j] / pivotVal;
+ for (int j = colIdx + 1; j < numCol; j++)
+ tableau[i][j] = tableau[i][j] - pivotCol[i] * pivotRow[j] / pivotVal;
+ }
+ for (int j = 0; j < numCol; j++)
+ tableau[rowIdx][j] = pivotRow[j] / pivotVal;
+ for (int i = 0; i < numRows; i++)
+ tableau[i][colIdx] = 0;
+ tableau[rowIdx][colIdx] = 1;
+ }
- // update the tableau
- for (int i = 0; i < rowIdx; i++)
- {
- for (int j = 0; j < colIdx; j++)
- tableau[i][j] = tableau[i][j] - pivotCol[i]*pivotRow[j]/pivotVal; // tij - til*tkj/tkl for col l and row k
- for (int j = colIdx + 1; j < numCol; j++)
- tableau[i][j] = tableau[i][j] - pivotCol[i]*pivotRow[j]/pivotVal;
- }
- for (int i = rowIdx + 1; i < numRows; i++)
- {
- for (int j = 0; j < colIdx; j++)
- tableau[i][j] = tableau[i][j] - pivotCol[i]*pivotRow[j]/pivotVal;
- for (int j = colIdx + 1; j < numCol; j++)
- tableau[i][j] = tableau[i][j] - pivotCol[i]*pivotRow[j]/pivotVal;
- }
- for (int j = 0; j < numCol; j++)
- tableau[rowIdx][j] = pivotRow[j]/pivotVal;
- for (int i = 0; i < numRows; i++)
- tableau[i][colIdx] = 0;
- tableau[rowIdx][colIdx] = 1;
- }
+ /**
+ * Get the solution. Slack variables are removed.
+ *
+ * @return a table of variable values defining the current solution
+ */
+ public double[] getSolution() {
+ final double[] solution = new double[numVariables];
+ for (int i = 0; i < leftVectorsIdx.length; i++) {
+ final int originalIdx = originalColOrder[leftVectorsIdx[i]];
+ if (originalIdx >= 0)
+ solution[originalIdx] = xCol[i];
- /**
- * Get the solution. Slack variables are removed.
- * @return a table of variable values defining the current solution
- * */
- public double[] getSolution()
- {
- final double[] solution = new double[numVariables];
- for (int i = 0; i < leftVectorsIdx.length; i++)
- {
- final int originalIdx = originalColOrder[leftVectorsIdx[i]];
- if (originalIdx >= 0)
- solution[originalIdx] = xCol[i];
+ }
+ return solution;
+ }
- }
- return solution;
- }
- /**
- * Modify scores in the tableau according to a set of scores for the constraints
- * */
- public void modifyScores(final double[] c)
- {
- final double[] scoreValues = new double[numCol];
- for (int i = 0; i < numCol; i++)
- {
- if (originalColOrder[i] >= 0)
- scoreValues[i] = c[originalColOrder[i]];
- else
- scoreValues[i] = 0;
- }
- for (int i = 0; i < numCol; i++)
- {
- double pi = 0;
- for (int j = 0; j < numRows; j++)
- {
- if (leftVectorsIdx[j] >= 0)
- pi += tableau[j][i]*scoreValues[leftVectorsIdx[j]];
- // else do nothing as unit vectors have no cost
- }
- if (originalColOrder[i] >= 0)
- scoreRow[i] = pi - c[originalColOrder[i]];
- else
- scoreRow[i] = pi;
- }
- scoreValue = 0;
- for (int i = 0; i < numRows; i++)
- {
- if (leftVectorsIdx[i] >= 0 && originalColOrder[leftVectorsIdx[i]] >= 0)
- {
- scoreValue += c[originalColOrder[leftVectorsIdx[i]]]*xCol[i];
- }
- }
- }
+ /**
+ * Modify scores in the tableau according to a set of scores for the constraints
+ */
+ public void modifyScores(final double[] c) {
+ final double[] scoreValues = new double[numCol];
+ for (int i = 0; i < numCol; i++) {
+ if (originalColOrder[i] >= 0)
+ scoreValues[i] = c[originalColOrder[i]];
+ else
+ scoreValues[i] = 0;
+ }
+ for (int i = 0; i < numCol; i++) {
+ double pi = 0;
+ for (int j = 0; j < numRows; j++) {
+ if (leftVectorsIdx[j] >= 0)
+ pi += tableau[j][i] * scoreValues[leftVectorsIdx[j]];
+ // else do nothing as unit vectors have no cost
+ }
+ if (originalColOrder[i] >= 0)
+ scoreRow[i] = pi - c[originalColOrder[i]];
+ else
+ scoreRow[i] = pi;
+ }
+ scoreValue = 0;
+ for (int i = 0; i < numRows; i++) {
+ if (leftVectorsIdx[i] >= 0 && originalColOrder[leftVectorsIdx[i]] >= 0) {
+ scoreValue += c[originalColOrder[leftVectorsIdx[i]]] * xCol[i];
+ }
+ }
+ }
- /**
- * Modify scores in the tableau according to a set of scores for the constraints
- * Same as modifyScoresBis. Just a backup function. Still working on it.
- * */
- public void modifyScoresBis(final double[] c)
- {
- final double[] scoreValues = new double[numCol];
- for (int i = 0; i < numCol; i++)
- {
- if (originalColOrder[i] >= 0)
- scoreValues[i] = c[originalColOrder[i]];
- else
- scoreValues[i] = 0;
- }
- for (int i = 0; i < numCol; i++)
- {
- double pi = 0;
- for (int j = 0; j < numRows; j++)
- pi += tableau[j][i]*scoreValues[leftVectorsIdx[j]];
- scoreRow[i] = pi - scoreValues[i];
- }
- scoreValue = 0;
- for (int i = 0; i < numRows; i++)
- scoreValue += scoreValues[leftVectorsIdx[i]]*xCol[i];
- }
+ /**
+ * Modify scores in the tableau according to a set of scores for the constraints
+ * Same as modifyScoresBis. Just a backup function. Still working on it.
+ */
+ public void modifyScoresBis(final double[] c) {
+ final double[] scoreValues = new double[numCol];
+ for (int i = 0; i < numCol; i++) {
+ if (originalColOrder[i] >= 0)
+ scoreValues[i] = c[originalColOrder[i]];
+ else
+ scoreValues[i] = 0;
+ }
+ for (int i = 0; i < numCol; i++) {
+ double pi = 0;
+ for (int j = 0; j < numRows; j++)
+ pi += tableau[j][i] * scoreValues[leftVectorsIdx[j]];
+ scoreRow[i] = pi - scoreValues[i];
+ }
+ scoreValue = 0;
+ for (int i = 0; i < numRows; i++)
+ scoreValue += scoreValues[leftVectorsIdx[i]] * xCol[i];
+ }
}
\ No newline at end of file