svg-layout-designer-react/src/utils/simplex.ts

/**
 * @module {Simplex} Apply the simplex algorithm
 * https://www.imse.iastate.edu/files/2015/08/Explanation-of-Simplex-Method.docx
 */

/**
 * Apply the simplex algorithms to the minimum widths
 *
 * Note: Some optimizations were made to improve performance in order to solve
 *  with max(minWidths). In point of fact most coefficient are equal to 1 or -1.
 *
 * Let the following format be the linear problem :
 * x >= b are the minimum widths constraint
 * sum(x) <= b is the maximum width constraint
 * s are slack variables
 * @param minWidths
 * @param requiredMaxWidth
 * @returns
 */
export function Simplex(minWidths: number[], maxWidths: number[], requiredMaxWidth: number): number[] {
  /// 1) standardized the equations
  // add the min widths constraints
  const constraints = minWidths.map(minWidth => minWidth * -1);

  /// 2) Create the initial matrix
  // get row length (nVariables + nConstraints + 1 (z) + 1 (b))
  const nVariables = minWidths.length;
  const nConstraints = constraints.length;
  const rowlength =
    minWidths.length + // min constraints
    maxWidths.length + // max constraints
    nConstraints + 1 + // slack variables
    1 + // z
    1; // b
  const matrix = GetInitialMatrix(constraints, maxWidths, requiredMaxWidth, rowlength);

  /// Apply the algorithm
  const finalMatrix = ApplyMainLoop(matrix, rowlength);

  // 5) read the solutions
  const solutions: number[] = GetSolutions(nVariables + nConstraints + 1, finalMatrix);
  return solutions;
}

/**
 * Specific to min widths algorithm
 * Get the initial matrix from the maximum constraints
 * and the number of variables
 * @param minConstraints
 * @param rowlength
 * @param nVariables
 * @returns
 */
function GetInitialMatrix(
  minConstraints: number[],
  maxConstraints: number[],
  objectiveConstraint: number,
  rowlength: number
): number[][] {
  const nVariables = maxConstraints.length;
  const constraints = minConstraints.concat(maxConstraints);
  constraints.push(objectiveConstraint);
  const matrix = constraints.map((constraint, index) => {
    const row: number[] = Array(rowlength).fill(0);

    // insert the variable coefficient a of a*x
    if (index < nVariables) {
      // insert the the variable coefficient of the minimum/maximum widths constraints (negative identity matrix)
      row[index] = -1;
    } else if (index < (2 * nVariables)) {
      row[index - (nVariables)] = 1;
    } else {
      // insert the the variable coefficient of the maximum desired width constraint
      row.fill(1, 0, nVariables);
    }

    // insert the slack variable coefficient b of b*s (identity matrix)
    row[index + nVariables] = 1;

    // insert the constraint coefficient (b)
    row[rowlength - 1] = constraint;
    return row;
  });

  // add objective function in the last row
  const row: number[] = Array(rowlength).fill(0);

  // insert z coefficient
  row[rowlength - 2] = 1;

  // insert variable coefficients
  row.fill(-1, 0, nVariables);
  matrix.push(row);
  return matrix;
}

function GetAllIndexes(arr: number[], val: number): number[] {
  const indexes = []; let i = -1;
  while ((i = arr.indexOf(val, i + 1)) !== -1) {
    indexes.push(i);
  }
  return indexes;
}

/**
 * Apply the main loop of the simplex algorithm and return the final matrix:
 * - While the last row of the matrix has negative values :
 *  - 1) find the column with the smallest negative coefficient in the last row
 *  - 2) in that column, find the pivot by selecting the row with the smallest ratio
 *      such as ratio = constraint of last column / coefficient of the selected row of the selected column
 *  - 3) create the new matrix such as:
 *  - 4) the selected column must have 1 in the pivot and zeroes in the other rows
 *  - 5) in the selected rows other columns (other than the selected column)
 *      must be divided by that pivot: coef / pivot
 *  - 6) for the others cells, apply the pivot: new value = (-coefficient in the old col) * (coefficient in the new row) + old value
 *  - 7) if in the new matrix there are still negative values in the last row,
 *      redo the algorithm with the new matrix as the base matrix
 *  - 8) otherwise returns the basic variable such as
 *      a basic variable is defined by a single 1 and only zeroes in its column
 *      other variables are equal to zeroes
 * @param oldMatrix
 * @param rowlength
 * @returns
 */
function ApplyMainLoop(oldMatrix: number[][], rowlength: number): number[][] {
  let matrix = oldMatrix;
  const maxTries = oldMatrix.length * 2;
  let tries = maxTries;
  const indexesTried: Record<number, number> = {};
  while (matrix[matrix.length - 1].some((v: number) => v < 0) && tries > 0) {
    // 1) find the index with smallest coefficient (O(n)+)
    const lastRow = matrix[matrix.length - 1];
    const min = Math.min(...lastRow);
    const indexes = GetAllIndexes(lastRow, min);
    // to avoid infinite loop try to select the least used selected index
    const pivotColIndex = GetLeastUsedIndex(indexes, indexesTried);
    // record the usage of index by incrementing
    indexesTried[pivotColIndex] = indexesTried[pivotColIndex] !== undefined ? indexesTried[pivotColIndex] + 1 : 1;

    // 2) find the smallest non negative non null ratio bi/xij (O(m))
    const ratios = [];
    for (let i = 0; i <= matrix.length - 2; i++) {
      const coefficient = matrix[i][pivotColIndex];
      const constraint = matrix[i][rowlength - 1];
      if (coefficient === 0) {
        ratios.push(Infinity);
        continue;
      }
      const ratio = constraint / coefficient;
      if (ratio < 0) {
        ratios.push(Infinity);
        continue;
      }
      ratios.push(ratio);
    }
    const minRatio = Math.min(...ratios);
    const pivotRowIndex = ratios.indexOf(minRatio); // i

    /// Init the new matrix
    const newMatrix: number[][] = structuredClone(matrix);
    const pivot = matrix[pivotRowIndex][pivotColIndex];

    // 3) apply on the pivot row the inverse of the pivot
    const newPivotRow = newMatrix[pivotRowIndex];
    newPivotRow.forEach((coef, colIndex) => {
      newPivotRow[colIndex] = coef / pivot;
    });

    // 4) update all values
    newMatrix.forEach((row, rowIndex) => {
      if (rowIndex === pivotRowIndex) {
        return;
      }

      row.forEach((coef, colIndex) => {
        if (colIndex === pivotColIndex) {
          // set zeroes on pivot col
          row[colIndex] = 0;
          return;
        }

        // update value = old value + ((-old coef of pivot column) * (new coef of pivot row))
        row[colIndex] = coef + (-matrix[rowIndex][pivotColIndex] * newMatrix[pivotRowIndex][colIndex]);
      });
    });

    matrix = newMatrix;
    tries--;
  }

  if (tries === 0) {
    console.table(matrix);
    throw new Error('[Flex] Simplexe: Could not find a solution');
  }
  console.debug(`Simplex was solved in ${maxTries - tries} tries`);
  return matrix;
}

/**
 * Get the solutions from the final matrix
 *
 * @param {number} nCols Number of solutions that you want to obtain
 * @param {number[][]} finalMatrix Final matrix after the algorithm is applied
 * @return {*}  {number[]} A list of solutions of the final matrix
 */
function GetSolutions(nCols: number, finalMatrix: number[][]): number[] {
  const solutions: number[] = Array(nCols).fill(0);
  for (let i = 0; i < nCols; i++) {
    const counts: Record<number, number> = {};
    const col: number[] = [];
    for (let j = 0; j < finalMatrix.length; j++) {
      const row = finalMatrix[j];
      counts[row[i]] = counts[row[i]] !== undefined ? counts[row[i]] + 1 : 1;
      col.push(row[i]);
    }

    // a basic variable has a single 1 and only zeroes in the column
    const nRows = finalMatrix.length;
    const isBasic = counts[1] === 1 && counts[0] === (nRows - 1);
    if (isBasic) {
      const oneIndex = col.indexOf(1);
      const row = finalMatrix[oneIndex];
      solutions[i] = row[row.length - 1];
    } else {
      solutions[i] = 0;
    }
  }
  return solutions;
}

/**
 * Returns the least used index from the indexesTried
 * @param indexes Indexes of all occurences
 * @param indexesTried Record of indexes. Count the number of times the index was used.
 * @returns The least used index
 */
function GetLeastUsedIndex(indexes: number[], indexesTried: Record<number, number>): number {
  let minUsed = Infinity;
  let minIndex = -1;
  for (const index of indexes) {
    const occ = indexesTried[index];
    if (occ === undefined) {
      minIndex = index;
      break;
    }

    if (occ < minUsed) {
      minIndex = index;
      minUsed = occ;
    }
  }
  return minIndex;
}