From 2b0e0a90bb0af43957e565a388c15c282f06bac8 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu Date: Tue, 13 Jun 2023 18:13:09 +0100 Subject: [PATCH 01/13] feat: added simplex.py --- other/simplex.py | 205 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 205 insertions(+) create mode 100644 other/simplex.py diff --git a/other/simplex.py b/other/simplex.py new file mode 100644 index 000000000000..78b0356c5268 --- /dev/null +++ b/other/simplex.py @@ -0,0 +1,205 @@ +# Simplex + +import numpy as np + + +class Tableau: + """Generate and operate on simplex tableaus""" + def __init__(self, tableau, n_vars, n_art_vars): + """Initialise Tableau class + + Args: + lin_prog : list[str] + Line separated string input in list + n_art_vars : int + Number of artificial/ surplus variables needed + """ + + # Initial tableau with no entries + self.tableau = tableau + + self.n_art_vars = n_art_vars + + # 2 if there are >= or == constraints (nonstandard), 1 otherwise (std) + self.n_stages = (n_art_vars > 0) + 1 + + # Number of rows the initial simplex tableau will have + self.n_rows = len(tableau[:, 0]) + + # Number of columns of the initial simplex tableau + self.n_cols = len(tableau[0]) + + # Number of slack variables added to make inequalities into equalities + self.n_slack = self.n_rows - self.n_stages + + # Number of decision variables x1, x2, x3... + self.n_vars = n_vars + + # Values of non-basic variables following iterations + self.output_dict = {} + + # Objectives for each stage + self.objectives = ["max"] + + if n_art_vars: + self.objectives.append("min") + + self.col_titles = [] + + # Index of current pivot row and column + self.row_idx = None + self.col_idx = None + + # Does objective row only contain (non)-negative values? + self.stop_iter = False + + def generate_col_titles(self): + """Simplex tableau contains variable, slack, artificial, and RHS + columns e.g. x_1, x_2, s_1, s_2, a_1, RHS + """ + string_starts = ["x_", "s_", "a_"] + constants = self.n_vars, self.n_slack, self.n_art_vars + titles = [] + for i in range(3): + for j in range(constants[i]): + titles.append(string_starts[i] + str(j + 1)) + titles.append("RHS") + self.col_titles = titles + + def find_pivot(self): + """Finds the pivot row and column.""" + tableau = self.tableau + objective = self.objectives[-1] + + # Find entries of highest magnitude in objective rows + sign = (objective == "min") - (objective == "max") + self.col_idx = np.argmax(sign * tableau[0, :-1]) + + # Check if choice is valid, or if iteration must be stopped + if sign * self.tableau[0, self.col_idx] <= 0: + self.stop_iter = True + return + + # Pivot row is chosen as having the lowest quotient when elements of + # the pivot column divide the right-hand side + + # Slice excluding the objective rows + s = slice(self.n_stages, self.n_rows) + + # RHS + dividend = tableau[s, -1] + + # Elements of pivot column within slice + divisor = tableau[s, self.col_idx] + + # Array filled with nans + nans = np.full(self.n_rows - self.n_stages, np.nan) + + # If element in pivot column is greater than zero, return quotient + # or nan otherwise + quotients = np.divide(dividend, divisor, out=nans, where=divisor > 0) + + # Arg of minimum quotient excluding the nan values. `n_stages` is added + # to compensate for earlier exclusion of objective columns + self.row_idx = np.nanargmin(quotients) + self.n_stages + + def pivot(self): + """Pivots on value on the intersection of pivot row and column.""" + + # Avoid changes to original tableau + piv_row = self.tableau[self.row_idx].copy() + + piv_val = piv_row[self.col_idx] + + # Entry becomes 1 + piv_row *= 1 / piv_val + + # Variable in pivot column becomes basic, ie the only non-zero entry + for idx, coeff in enumerate(self.tableau[:, self.col_idx]): + self.tableau[idx] += -coeff * piv_row + self.tableau[self.row_idx] = piv_row + + def change_stage(self): + """Exits first phase of the two-stage method by deleting artificial + rows and columns, or completes the algorithm if exiting the standard + case.""" + # Objective of original objective row remains + self.objectives.pop() + + if not self.objectives: + return + + # Slice containing ids for artificial columns + s = slice(-self.n_art_vars - 1, -1) + + # Delete the artificial variable columns + self.tableau = np.delete(self.tableau, s, axis=1) + + # Delete the objective row of the first stage + self.tableau = np.delete(self.tableau, 0, axis=0) + + self.n_stages = 1 + self.n_rows -= 1 + self.n_art_vars = 0 + self.stop_iter = False + + def run_simp(self): + """Recursively operate on tableau until objective function cannot be + improved further""" + # If optimal solution reached + if not self.objectives: + self.interpret_tableau() + raise Exception + + self.find_pivot() + if self.stop_iter: + self.change_stage() + self.run_simp() + self.pivot() + self.run_simp() + + def interpret_tableau(self): + """Given the final tableau, add the corresponding values of the basic + variables to the `output_dict`""" + # P = RHS of final tableau + self.output_dict["P"] = self.tableau[0, -1] + + n_current_cols = len(self.tableau[0]) + for i in range(n_current_cols): + + # Gives ids of nonzero entries in the ith column + nonzero = np.nonzero(self.tableau[:, i]) + n_nonzero = len(nonzero[0]) + + # First entry in the nonzero ids + nonzero_rowidx = nonzero[0][0] + nonzero_val = self.tableau[nonzero_rowidx, i] + + # If there is only one nonzero value in column, which is one + if n_nonzero == nonzero_val == 1: + rhs_val = self.tableau[nonzero_rowidx, -1] + self.output_dict[self.col_titles[i]] = rhs_val + +tableau = np.array([ + [2, 0, 0, 0, -1, -1, 0, 0, 20], + [-2, -3, -1, 0, 0, 0, 0, 0, 0], + [1, 1, 1, 1, 0, 0, 0, 0, 40], + [2, 1, -1, 0, -1, 0, 1, 0, 10], + [0, -1, 1, 0, 0, -1, 0, 1, 10] + ], dtype=float) + +n_vars = 3 +n_art_vars = 2 + +def main(): + t = Tableau(tableau, n_vars, n_art_vars) + t.generate_col_titles() + try: + t.run_simp() + except Exception: + # If optimal solution is reached + for key, value in t.output_dict.items(): + print(key, ":", value) + +if __name__ == "__main__": + main() \ No newline at end of file From 2f19c85c779cd0d93bcc2cfbca88f8546e308e16 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu Date: Fri, 16 Jun 2023 21:04:40 +0100 Subject: [PATCH 02/13] added docstrings --- linear_programming/simplex.py | 327 ++++++++++++++++++++++++++++++++++ other/simplex.py | 205 --------------------- 2 files changed, 327 insertions(+), 205 deletions(-) create mode 100644 linear_programming/simplex.py delete mode 100644 other/simplex.py diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py new file mode 100644 index 000000000000..f6c17136da9a --- /dev/null +++ b/linear_programming/simplex.py @@ -0,0 +1,327 @@ +""" +Python implementation of the simplex algorithm for solving linear programs in +tabular form with +- `>=`, `<=`, and `=` constraints and +- each variable `x1, x2, ...>= 0`. + +See https://gist.github.com/imengus/f9619a568f7da5bc74eaf20169a24d98 for how +to convert linear programs to simplex tableaus, and the steps taken in the +simplex algorithm. + +Resources: +https://en.wikipedia.org/wiki/Linear_programming +https://en.wikipedia.org/wiki/Simplex_algorithm + +https://towardsdatascience.com/ \ + a-beginners-guide-to-linear-programming-and-the-simplex-algorithm \ + -87db017e92b4 +""" +from typing import Any + +import numpy as np + + +class Tableau: + """Operate on simplex tableaus""" + + def __init__(self, tableau: np.ndarray, n_vars: int) -> None: + """ + >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) + >>> t = Tableau(tableau, 2) + Traceback (most recent call last): + ... + ValueError: RHS must be > 0 + """ + rhs = tableau[:, -1] + if np.any(rhs, where=rhs < 0): + raise ValueError("RHS must be > 0") + + if tableau is None: + return + + self.tableau = tableau + self.n_rows = len(tableau[:, 0]) + self.n_cols = len(tableau[0]) + + # Number of decision variables x1, x2, x3... + self.n_vars = n_vars + + # Number of artificial variables to be minimised + self.n_art_vars = len(np.where(tableau[self.n_vars : -1] == -1)[0]) + + # 2 if there are >= or == constraints (nonstandard), 1 otherwise (std) + self.n_stages = (self.n_art_vars > 0) + 1 + + # Number of slack variables added to make inequalities into equalities + self.n_slack = self.n_rows - self.n_stages + + # Objectives for each stage + self.objectives = ["max"] + + # In two stage simplex, first minimise then maximise + if self.n_art_vars: + self.objectives.append("min") + + self.col_titles = [""] + + # Index of current pivot row and column + self.row_idx = None + self.col_idx = None + + # Does objective row only contain (non)-negative values? + self.stop_iter = False + + @staticmethod + def generate_col_titles(*args: int) -> list[str]: + """Generate column titles for tableau of specific dimensions + + >>> Tableau.generate_col_titles(2, 3, 1) + ['x1', 'x2', 's1', 's2', 's3', 'a1', 'RHS'] + """ + # decision | slack | artificial + string_starts = ["x", "s", "a"] + titles = [] + for i in range(3): + for j in range(args[i]): + titles.append(string_starts[i] + str(j + 1)) + titles.append("RHS") + return titles + + def find_pivot(self, tableau: np.ndarray) -> tuple[Any, Any]: + """Finds the pivot row and column. + >>> tableau = np.array([ + ... [-2,-1,0,0,0], + ... [3,1,1,0,6], + ... [1,2,0,1,7.] + ... ]) + >>> t = Tableau(tableau, 2) + >>> t.find_pivot(t.tableau) + (1, 0) + """ + objective = self.objectives[-1] + + # Find entries of highest magnitude in objective rows + sign = (objective == "min") - (objective == "max") + col_idx = np.argmax(sign * tableau[0, : self.n_vars]) + + # Choice is only valid if below 0 for maximise, and above for minimise + if sign * self.tableau[0, col_idx] <= 0: + self.stop_iter = True + return 0, 0 + + # Pivot row is chosen as having the lowest quotient when elements of + # the pivot column divide the right-hand side + + # Slice excluding the objective rows + s = slice(self.n_stages, self.n_rows) + + # RHS + dividend = tableau[s, -1] + + # Elements of pivot column within slice + divisor = tableau[s, col_idx] + + # Array filled with nans + nans = np.full(self.n_rows - self.n_stages, np.nan) + + # If element in pivot column is greater than zeron_stages, return + # quotient or nan otherwise + quotients = np.divide(dividend, divisor, out=nans, where=divisor > 0) + + # Arg of minimum quotient excluding the nan values. n_stages is added + # to compensate for earlier exclusion of objective columns + row_idx = np.nanargmin(quotients) + self.n_stages + return row_idx, col_idx + + def pivot(self, tableau: np.ndarray, row_idx: int, col_idx: int) -> np.ndarray: + """Pivots on value on the intersection of pivot row and column. + + >>> tableau = np.array([[-2,-3,0,0,0],[1,3,1,0,4],[3,1,0,1,4.]]) + >>> t = Tableau(tableau, 2) + >>> t.pivot(t.tableau, 1, 0).tolist() + ... # doctest: +NORMALIZE_WHITESPACE + [[0.0, 3.0, 2.0, 0.0, 8.0], + [1.0, 3.0, 1.0, 0.0, 4.0], + [0.0, -8.0, -3.0, 1.0, -8.0]] + """ + # Avoid changes to original tableau + piv_row = tableau[row_idx].copy() + + piv_val = piv_row[col_idx] + + # Entry becomes 1 + piv_row *= 1 / piv_val + + # Variable in pivot column becomes basic, ie the only non-zero entry + for idx, coeff in enumerate(tableau[:, col_idx]): + tableau[idx] += -coeff * piv_row + tableau[row_idx] = piv_row + return tableau + + def change_stage(self, tableau: np.ndarray, objectives: list[Any]) -> np.ndarray: + """Exits first phase of the two-stage method by deleting artificial + rows and columns, or completes the algorithm if exiting the standard + case. + + >>> tableau = np.array([ + ... [3,3,-1,-1,0,0,4], + ... [2,1,0,0,0,0,0.], + ... [1,2,-1,0,1,0,2], + ... [2,1,0,-1,0,1,2] + ... ]) + >>> t = Tableau(tableau, 2) + >>> t.change_stage(t.tableau, t.objectives).tolist() + ... # doctest: +NORMALIZE_WHITESPACE + [[2.0, 1.0, 0.0, 0.0, 0.0, 0.0], + [1.0, 2.0, -1.0, 0.0, 1.0, 2.0], + [2.0, 1.0, 0.0, -1.0, 0.0, 2.0]] + """ + # Objective of original objective row remains + objectives.pop() + + if not objectives: + return tableau + + # Slice containing ids for artificial columns + s = slice(-self.n_art_vars - 1, -1) + + # Delete the artificial variable columns + tableau = np.delete(tableau, s, axis=1) + + # Delete the objective row of the first stage + tableau = np.delete(tableau, 0, axis=0) + + self.n_stages = 1 + self.n_rows -= 1 + self.n_art_vars = 0 + self.stop_iter = False + self.objectives = objectives + return tableau + + def run_simplex(self) -> dict[Any, Any]: + """Operate on tableau until objective function cannot be + improved further. + + # Standard linear program: + Max: x1 + x2 + ST: x1 + 3x2 <= 4 + 3x1 + x2 <= 4 + >>> tableau = np.array([[-1,-1,0,0,0],[1,3,1,0,4],[3,1,0,1,4.]]) + >>> t = Tableau(tableau, 2) + >>> t.run_simplex() + {'P': 2.0, 'x1': 1.0, 'x2': 1.0} + + # Optimal tableau input: + >>> tableau = np.array([ + ... [0,0,0.25,0.25,2], + ... [0,1,0.375,-0.125,1], + ... [1,0,-0.125,0.375,1] + ... ]) + >>> t = Tableau(tableau, 2) + >>> t.run_simplex() + {'P': 2.0, 'x1': 1.0, 'x2': 1.0} + + # Non-standard: >= constraints + Max: 2x1 + 3x2 + x3 + ST: x1 + x2 + x3 <= 40 + 2x1 + x2 - x3 >= 10 + - x2 + x3 >= 10 + >>> tableau = np.array([ + ... [2,0,0,0,-1,-1,0,0,20], + ... [-2,-3,-1,0,0,0,0,0,0], + ... [1,1,1,1,0,0,0,0,40], + ... [2,1,-1,0,-1,0,1,0,10], + ... [0,-1,1,0,0,-1,0,1,10.] + ... ]) + >>> t = Tableau(tableau, 3) + >>> t.run_simplex() + {'P': 70.0, 'x1': 10.0, 'x2': 10.0, 'x3': 20.0} + + # Non standard: minimisation and equalities + Min: x1 + x2 + ST: 2x1 + x2 = 12 + 6x1 + 5x2 = 40 + >>> tableau = np.array([ + ... [8,6,0,-1,0,-1,0,0,52], + ... [1,1,0,0,0,0,0,0,0], + ... [2,1,1,0,0,0,0,0,12], + ... [2,1,0,-1,0,0,1,0,12], + ... [6,5,0,0,1,0,0,0,40], + ... [6,5,0,0,0,-1,0,1,40.] + ... ]) + >>> t = Tableau(tableau, 2) + >>> t.run_simplex() + {'P': 7.0, 'x1': 5.0, 'x2': 2.0} + """ + # Stop simplex algorithm from cycling. + iter_num = 0 + + while iter_num < 100: + # Completion of each stage removes an objective. If both stages + # are complete, then no objectives are left + if not self.objectives: + self.col_titles = self.generate_col_titles( + self.n_vars, self.n_slack, self.n_art_vars + ) + + # Find the values of each variable at optimal solution + return self.interpret_tableau(self.tableau, self.col_titles) + + row_idx, col_idx = self.find_pivot(self.tableau) + + # If there are no more negative values in objective row + if self.stop_iter: + # Delete artificial variable columns and rows. Update attributes + self.tableau = self.change_stage(self.tableau, self.objectives) + + # Pivot again + continue + + self.tableau = self.pivot(self.tableau, row_idx, col_idx) + iter_num += 1 + return {} + + def interpret_tableau( + self, tableau: np.ndarray, col_titles: list[str] + ) -> dict[str, float]: + """Given the final tableau, add the corresponding values of the basic + decision variables to the `output_dict` + >>> tableau = np.array([ + ... [0,0,0.875,0.375,5], + ... [0,1,0.375,-0.125,1], + ... [1,0,-0.125,0.375,1] + ... ]) + >>> t = Tableau(tableau, 2) + >>> t.interpret_tableau(tableau, ["x1", "x2", "s1", "s2", "RHS"]) + {'P': 5.0, 'x1': 1.0, 'x2': 1.0} + """ + output_dict = {} + + # P = RHS of final tableau + output_dict["P"] = abs(tableau[0, -1]) + + for i in range(self.n_vars): + # Gives ids of nonzero entries in the ith column + nonzero = np.nonzero(tableau[:, i]) + n_nonzero = len(nonzero[0]) + + # First entry in the nonzero ids + nonzero_rowidx = nonzero[0][0] + nonzero_val = tableau[nonzero_rowidx, i] + + # If there is only one nonzero value in column, which is one + if n_nonzero == nonzero_val == 1: + rhs_val = tableau[nonzero_rowidx, -1] + output_dict[col_titles[i]] = rhs_val + # Check for basic variables + for title in col_titles: + # Don't add RHS or slack variables to output dict + if title[0] not in "R-s": + output_dict.setdefault(title, 0) + return output_dict + + +if __name__ == "__main__": + import doctest + + doctest.testmod() diff --git a/other/simplex.py b/other/simplex.py deleted file mode 100644 index 78b0356c5268..000000000000 --- a/other/simplex.py +++ /dev/null @@ -1,205 +0,0 @@ -# Simplex - -import numpy as np - - -class Tableau: - """Generate and operate on simplex tableaus""" - def __init__(self, tableau, n_vars, n_art_vars): - """Initialise Tableau class - - Args: - lin_prog : list[str] - Line separated string input in list - n_art_vars : int - Number of artificial/ surplus variables needed - """ - - # Initial tableau with no entries - self.tableau = tableau - - self.n_art_vars = n_art_vars - - # 2 if there are >= or == constraints (nonstandard), 1 otherwise (std) - self.n_stages = (n_art_vars > 0) + 1 - - # Number of rows the initial simplex tableau will have - self.n_rows = len(tableau[:, 0]) - - # Number of columns of the initial simplex tableau - self.n_cols = len(tableau[0]) - - # Number of slack variables added to make inequalities into equalities - self.n_slack = self.n_rows - self.n_stages - - # Number of decision variables x1, x2, x3... - self.n_vars = n_vars - - # Values of non-basic variables following iterations - self.output_dict = {} - - # Objectives for each stage - self.objectives = ["max"] - - if n_art_vars: - self.objectives.append("min") - - self.col_titles = [] - - # Index of current pivot row and column - self.row_idx = None - self.col_idx = None - - # Does objective row only contain (non)-negative values? - self.stop_iter = False - - def generate_col_titles(self): - """Simplex tableau contains variable, slack, artificial, and RHS - columns e.g. x_1, x_2, s_1, s_2, a_1, RHS - """ - string_starts = ["x_", "s_", "a_"] - constants = self.n_vars, self.n_slack, self.n_art_vars - titles = [] - for i in range(3): - for j in range(constants[i]): - titles.append(string_starts[i] + str(j + 1)) - titles.append("RHS") - self.col_titles = titles - - def find_pivot(self): - """Finds the pivot row and column.""" - tableau = self.tableau - objective = self.objectives[-1] - - # Find entries of highest magnitude in objective rows - sign = (objective == "min") - (objective == "max") - self.col_idx = np.argmax(sign * tableau[0, :-1]) - - # Check if choice is valid, or if iteration must be stopped - if sign * self.tableau[0, self.col_idx] <= 0: - self.stop_iter = True - return - - # Pivot row is chosen as having the lowest quotient when elements of - # the pivot column divide the right-hand side - - # Slice excluding the objective rows - s = slice(self.n_stages, self.n_rows) - - # RHS - dividend = tableau[s, -1] - - # Elements of pivot column within slice - divisor = tableau[s, self.col_idx] - - # Array filled with nans - nans = np.full(self.n_rows - self.n_stages, np.nan) - - # If element in pivot column is greater than zero, return quotient - # or nan otherwise - quotients = np.divide(dividend, divisor, out=nans, where=divisor > 0) - - # Arg of minimum quotient excluding the nan values. `n_stages` is added - # to compensate for earlier exclusion of objective columns - self.row_idx = np.nanargmin(quotients) + self.n_stages - - def pivot(self): - """Pivots on value on the intersection of pivot row and column.""" - - # Avoid changes to original tableau - piv_row = self.tableau[self.row_idx].copy() - - piv_val = piv_row[self.col_idx] - - # Entry becomes 1 - piv_row *= 1 / piv_val - - # Variable in pivot column becomes basic, ie the only non-zero entry - for idx, coeff in enumerate(self.tableau[:, self.col_idx]): - self.tableau[idx] += -coeff * piv_row - self.tableau[self.row_idx] = piv_row - - def change_stage(self): - """Exits first phase of the two-stage method by deleting artificial - rows and columns, or completes the algorithm if exiting the standard - case.""" - # Objective of original objective row remains - self.objectives.pop() - - if not self.objectives: - return - - # Slice containing ids for artificial columns - s = slice(-self.n_art_vars - 1, -1) - - # Delete the artificial variable columns - self.tableau = np.delete(self.tableau, s, axis=1) - - # Delete the objective row of the first stage - self.tableau = np.delete(self.tableau, 0, axis=0) - - self.n_stages = 1 - self.n_rows -= 1 - self.n_art_vars = 0 - self.stop_iter = False - - def run_simp(self): - """Recursively operate on tableau until objective function cannot be - improved further""" - # If optimal solution reached - if not self.objectives: - self.interpret_tableau() - raise Exception - - self.find_pivot() - if self.stop_iter: - self.change_stage() - self.run_simp() - self.pivot() - self.run_simp() - - def interpret_tableau(self): - """Given the final tableau, add the corresponding values of the basic - variables to the `output_dict`""" - # P = RHS of final tableau - self.output_dict["P"] = self.tableau[0, -1] - - n_current_cols = len(self.tableau[0]) - for i in range(n_current_cols): - - # Gives ids of nonzero entries in the ith column - nonzero = np.nonzero(self.tableau[:, i]) - n_nonzero = len(nonzero[0]) - - # First entry in the nonzero ids - nonzero_rowidx = nonzero[0][0] - nonzero_val = self.tableau[nonzero_rowidx, i] - - # If there is only one nonzero value in column, which is one - if n_nonzero == nonzero_val == 1: - rhs_val = self.tableau[nonzero_rowidx, -1] - self.output_dict[self.col_titles[i]] = rhs_val - -tableau = np.array([ - [2, 0, 0, 0, -1, -1, 0, 0, 20], - [-2, -3, -1, 0, 0, 0, 0, 0, 0], - [1, 1, 1, 1, 0, 0, 0, 0, 40], - [2, 1, -1, 0, -1, 0, 1, 0, 10], - [0, -1, 1, 0, 0, -1, 0, 1, 10] - ], dtype=float) - -n_vars = 3 -n_art_vars = 2 - -def main(): - t = Tableau(tableau, n_vars, n_art_vars) - t.generate_col_titles() - try: - t.run_simp() - except Exception: - # If optimal solution is reached - for key, value in t.output_dict.items(): - print(key, ":", value) - -if __name__ == "__main__": - main() \ No newline at end of file From face51655e79d5bc28a6404ea352f3d6e5642d80 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu <113149540+imengus@users.noreply.github.com> Date: Sat, 17 Jun 2023 00:29:19 +0100 Subject: [PATCH 03/13] Update linear_programming/simplex.py Co-authored-by: Caeden Perelli-Harris --- linear_programming/simplex.py | 1 - 1 file changed, 1 deletion(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index f6c17136da9a..913803eac22b 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -9,7 +9,6 @@ simplex algorithm. Resources: -https://en.wikipedia.org/wiki/Linear_programming https://en.wikipedia.org/wiki/Simplex_algorithm https://towardsdatascience.com/ \ From fc5c14d36ed6116a9879252e05f7492526d19d41 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu <113149540+imengus@users.noreply.github.com> Date: Sat, 17 Jun 2023 00:33:16 +0100 Subject: [PATCH 04/13] Update linear_programming/simplex.py Co-authored-by: Caeden Perelli-Harris --- linear_programming/simplex.py | 2 -- 1 file changed, 2 deletions(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index 913803eac22b..0037993579bf 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -35,8 +35,6 @@ def __init__(self, tableau: np.ndarray, n_vars: int) -> None: if np.any(rhs, where=rhs < 0): raise ValueError("RHS must be > 0") - if tableau is None: - return self.tableau = tableau self.n_rows = len(tableau[:, 0]) From 119b755b444ee3320c8ec7b75f012ced30ea6a02 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 16 Jun 2023 23:33:42 +0000 Subject: [PATCH 05/13] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- linear_programming/simplex.py | 1 - 1 file changed, 1 deletion(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index 0037993579bf..3664de42fdef 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -35,7 +35,6 @@ def __init__(self, tableau: np.ndarray, n_vars: int) -> None: if np.any(rhs, where=rhs < 0): raise ValueError("RHS must be > 0") - self.tableau = tableau self.n_rows = len(tableau[:, 0]) self.n_cols = len(tableau[0]) From ed34adc8dcfa253cb98d14112b159d09e27ddb2b Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu <113149540+imengus@users.noreply.github.com> Date: Sat, 17 Jun 2023 00:35:47 +0100 Subject: [PATCH 06/13] Update linear_programming/simplex.py Co-authored-by: Caeden Perelli-Harris --- linear_programming/simplex.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index 3664de42fdef..f3ef289830e9 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -21,16 +21,16 @@ class Tableau: - """Operate on simplex tableaus""" + """Operate on simplex tableaus + + >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) + >>> t = Tableau(tableau, 2) + Traceback (most recent call last): + ... + ValueError: RHS must be > 0 + """ def __init__(self, tableau: np.ndarray, n_vars: int) -> None: - """ - >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) - >>> t = Tableau(tableau, 2) - Traceback (most recent call last): - ... - ValueError: RHS must be > 0 - """ rhs = tableau[:, -1] if np.any(rhs, where=rhs < 0): raise ValueError("RHS must be > 0") From c6b4703babebe7885630ca5b1bf91c7391139ca1 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 16 Jun 2023 23:36:14 +0000 Subject: [PATCH 07/13] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- linear_programming/simplex.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index f3ef289830e9..020b49aa9865 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -22,7 +22,7 @@ class Tableau: """Operate on simplex tableaus - + >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) >>> t = Tableau(tableau, 2) Traceback (most recent call last): From d95b4f7b7e21d69c74ec3e315a09165afd4b0d47 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu Date: Sat, 17 Jun 2023 01:15:05 +0100 Subject: [PATCH 08/13] ruff fix Co-authored by: CaedenPH --- linear_programming/README.md | 19 +++++ linear_programming/simplex.py | 137 ++++++++++++++++------------------ 2 files changed, 83 insertions(+), 73 deletions(-) create mode 100644 linear_programming/README.md diff --git a/linear_programming/README.md b/linear_programming/README.md new file mode 100644 index 000000000000..4f520e173493 --- /dev/null +++ b/linear_programming/README.md @@ -0,0 +1,19 @@ +# Linear Programming +Sub-field of mathematical optimisation, where the best solution is found while satisfying linear functions. + +Taking this linear program as an example +``` + Maximise P = 8x1 + 6x2 + subject to: x1 + x2 <= 20 + x1 >= 5 + x2 = 6 + with: + x1 >= 0, x2 >= 0 +``` +x1 could be the number of blue cars and x2 the number of red cars produced. The company would be trying to maximise profit P. + +Resources: + +https://en.wikipedia.org/wiki/Linear_programming +https://brilliant.org/wiki/linear-programming/ +https://byjus.com/maths/linear-programming/ \ No newline at end of file diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index f6c17136da9a..f59433989377 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -4,17 +4,13 @@ - `>=`, `<=`, and `=` constraints and - each variable `x1, x2, ...>= 0`. -See https://gist.github.com/imengus/f9619a568f7da5bc74eaf20169a24d98 for how -to convert linear programs to simplex tableaus, and the steps taken in the -simplex algorithm. +See https://gist.github.com/imengus/f9619a568f7da5bc74eaf20169a24d98 for how to +convert linear programs to simplex tableaus, and the steps taken in the simplex +algorithm. Resources: -https://en.wikipedia.org/wiki/Linear_programming https://en.wikipedia.org/wiki/Simplex_algorithm - -https://towardsdatascience.com/ \ - a-beginners-guide-to-linear-programming-and-the-simplex-algorithm \ - -87db017e92b4 +https://tinyurl.com/simplex4beginners """ from typing import Any @@ -22,23 +18,20 @@ class Tableau: - """Operate on simplex tableaus""" + """Operate on simplex tableaus + + >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) + >>> t = Tableau(tableau, 2) + Traceback (most recent call last): + ... + ValueError: RHS must be > 0 + """ def __init__(self, tableau: np.ndarray, n_vars: int) -> None: - """ - >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) - >>> t = Tableau(tableau, 2) - Traceback (most recent call last): - ... - ValueError: RHS must be > 0 - """ rhs = tableau[:, -1] if np.any(rhs, where=rhs < 0): raise ValueError("RHS must be > 0") - if tableau is None: - return - self.tableau = tableau self.n_rows = len(tableau[:, 0]) self.n_cols = len(tableau[0]) @@ -77,7 +70,22 @@ def generate_col_titles(*args: int) -> list[str]: >>> Tableau.generate_col_titles(2, 3, 1) ['x1', 'x2', 's1', 's2', 's3', 'a1', 'RHS'] + + >>> Tableau.generate_col_titles() + Traceback (most recent call last): + ... + ValueError: Must provide n_vars, n_slack, and n_art_vars + >>> Tableau.generate_col_titles(-2, 3, 1) + Traceback (most recent call last): + ... + ValueError: All arguments must be non-negative integers """ + if len(args) != 3: + raise ValueError("Must provide n_vars, n_slack, and n_art_vars") + + if not all(x >= 0 and isinstance(x, int) for x in args): + raise ValueError("All arguments must be non-negative integers") + # decision | slack | artificial string_starts = ["x", "s", "a"] titles = [] @@ -89,12 +97,7 @@ def generate_col_titles(*args: int) -> list[str]: def find_pivot(self, tableau: np.ndarray) -> tuple[Any, Any]: """Finds the pivot row and column. - >>> tableau = np.array([ - ... [-2,-1,0,0,0], - ... [3,1,1,0,6], - ... [1,2,0,1,7.] - ... ]) - >>> t = Tableau(tableau, 2) + >>> t = Tableau(np.array([[-2,1,0,0,0], [3,1,1,0,6], [1,2,0,1,7.]]), 2) >>> t.find_pivot(t.tableau) (1, 0) """ @@ -136,8 +139,7 @@ def find_pivot(self, tableau: np.ndarray) -> tuple[Any, Any]: def pivot(self, tableau: np.ndarray, row_idx: int, col_idx: int) -> np.ndarray: """Pivots on value on the intersection of pivot row and column. - >>> tableau = np.array([[-2,-3,0,0,0],[1,3,1,0,4],[3,1,0,1,4.]]) - >>> t = Tableau(tableau, 2) + >>> t = Tableau(np.array([[-2,-3,0,0,0],[1,3,1,0,4],[3,1,0,1,4.]]), 2) >>> t.pivot(t.tableau, 1, 0).tolist() ... # doctest: +NORMALIZE_WHITESPACE [[0.0, 3.0, 2.0, 0.0, 8.0], @@ -158,28 +160,27 @@ def pivot(self, tableau: np.ndarray, row_idx: int, col_idx: int) -> np.ndarray: tableau[row_idx] = piv_row return tableau - def change_stage(self, tableau: np.ndarray, objectives: list[Any]) -> np.ndarray: + def change_stage(self, tableau: np.ndarray) -> np.ndarray: """Exits first phase of the two-stage method by deleting artificial rows and columns, or completes the algorithm if exiting the standard case. - >>> tableau = np.array([ - ... [3,3,-1,-1,0,0,4], - ... [2,1,0,0,0,0,0.], - ... [1,2,-1,0,1,0,2], - ... [2,1,0,-1,0,1,2] - ... ]) - >>> t = Tableau(tableau, 2) - >>> t.change_stage(t.tableau, t.objectives).tolist() + >>> t = Tableau(np.array([ + ... [3, 3, -1, -1, 0, 0, 4], + ... [2, 1, 0, 0, 0, 0, 0.], + ... [1, 2, -1, 0, 1, 0, 2], + ... [2, 1, 0, -1, 0, 1, 2] + ... ]), 2) + >>> t.change_stage(t.tableau).tolist() ... # doctest: +NORMALIZE_WHITESPACE [[2.0, 1.0, 0.0, 0.0, 0.0, 0.0], [1.0, 2.0, -1.0, 0.0, 1.0, 2.0], [2.0, 1.0, 0.0, -1.0, 0.0, 2.0]] """ # Objective of original objective row remains - objectives.pop() + self.objectives.pop() - if not objectives: + if not self.objectives: return tableau # Slice containing ids for artificial columns @@ -195,7 +196,6 @@ def change_stage(self, tableau: np.ndarray, objectives: list[Any]) -> np.ndarray self.n_rows -= 1 self.n_art_vars = 0 self.stop_iter = False - self.objectives = objectives return tableau def run_simplex(self) -> dict[Any, Any]: @@ -206,19 +206,16 @@ def run_simplex(self) -> dict[Any, Any]: Max: x1 + x2 ST: x1 + 3x2 <= 4 3x1 + x2 <= 4 - >>> tableau = np.array([[-1,-1,0,0,0],[1,3,1,0,4],[3,1,0,1,4.]]) - >>> t = Tableau(tableau, 2) - >>> t.run_simplex() + >>> Tableau(np.array([[-1,-1,0,0,0],[1,3,1,0,4],[3,1,0,1,4.]]), + ... 2).run_simplex() {'P': 2.0, 'x1': 1.0, 'x2': 1.0} # Optimal tableau input: - >>> tableau = np.array([ - ... [0,0,0.25,0.25,2], - ... [0,1,0.375,-0.125,1], - ... [1,0,-0.125,0.375,1] - ... ]) - >>> t = Tableau(tableau, 2) - >>> t.run_simplex() + >>> Tableau(np.array([ + ... [0, 0, 0.25, 0.25, 2], + ... [0, 1, 0.375, -0.125, 1], + ... [1, 0, -0.125, 0.375, 1] + ... ]), 2).run_simplex() {'P': 2.0, 'x1': 1.0, 'x2': 1.0} # Non-standard: >= constraints @@ -226,31 +223,27 @@ def run_simplex(self) -> dict[Any, Any]: ST: x1 + x2 + x3 <= 40 2x1 + x2 - x3 >= 10 - x2 + x3 >= 10 - >>> tableau = np.array([ - ... [2,0,0,0,-1,-1,0,0,20], - ... [-2,-3,-1,0,0,0,0,0,0], - ... [1,1,1,1,0,0,0,0,40], - ... [2,1,-1,0,-1,0,1,0,10], - ... [0,-1,1,0,0,-1,0,1,10.] - ... ]) - >>> t = Tableau(tableau, 3) - >>> t.run_simplex() + >>> Tableau(np.array([ + ... [2, 0, 0, 0, -1, -1, 0, 0, 20], + ... [-2, -3, -1, 0, 0, 0, 0, 0, 0], + ... [1, 1, 1, 1, 0, 0, 0, 0, 40], + ... [2, 1, -1, 0, -1, 0, 1, 0, 10], + ... [0, -1, 1, 0, 0, -1, 0, 1, 10.] + ... ]), 3).run_simplex() {'P': 70.0, 'x1': 10.0, 'x2': 10.0, 'x3': 20.0} # Non standard: minimisation and equalities Min: x1 + x2 ST: 2x1 + x2 = 12 6x1 + 5x2 = 40 - >>> tableau = np.array([ - ... [8,6,0,-1,0,-1,0,0,52], - ... [1,1,0,0,0,0,0,0,0], - ... [2,1,1,0,0,0,0,0,12], - ... [2,1,0,-1,0,0,1,0,12], - ... [6,5,0,0,1,0,0,0,40], - ... [6,5,0,0,0,-1,0,1,40.] - ... ]) - >>> t = Tableau(tableau, 2) - >>> t.run_simplex() + >>> Tableau(np.array([ + ... [8, 6, 0, -1, 0, -1, 0, 0, 52], + ... [1, 1, 0, 0, 0, 0, 0, 0, 0], + ... [2, 1, 1, 0, 0, 0, 0, 0, 12], + ... [2, 1, 0, -1, 0, 0, 1, 0, 12], + ... [6, 5, 0, 0, 1, 0, 0, 0, 40], + ... [6, 5, 0, 0, 0, -1, 0, 1, 40.] + ... ]), 2).run_simplex() {'P': 7.0, 'x1': 5.0, 'x2': 2.0} """ # Stop simplex algorithm from cycling. @@ -272,7 +265,7 @@ def run_simplex(self) -> dict[Any, Any]: # If there are no more negative values in objective row if self.stop_iter: # Delete artificial variable columns and rows. Update attributes - self.tableau = self.change_stage(self.tableau, self.objectives) + self.tableau = self.change_stage(self.tableau) # Pivot again continue @@ -295,10 +288,8 @@ def interpret_tableau( >>> t.interpret_tableau(tableau, ["x1", "x2", "s1", "s2", "RHS"]) {'P': 5.0, 'x1': 1.0, 'x2': 1.0} """ - output_dict = {} - # P = RHS of final tableau - output_dict["P"] = abs(tableau[0, -1]) + output_dict = {"P": abs(tableau[0, -1])} for i in range(self.n_vars): # Gives ids of nonzero entries in the ith column @@ -316,7 +307,7 @@ def interpret_tableau( # Check for basic variables for title in col_titles: # Don't add RHS or slack variables to output dict - if title[0] not in "R-s": + if title[0] not in "R-s-a": output_dict.setdefault(title, 0) return output_dict From 6691506af5d2955b35b1b2937d36d59ebd2638b4 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu Date: Sat, 17 Jun 2023 19:58:09 +0100 Subject: [PATCH 09/13] removed README to add in separate PR --- linear_programming/README.md | 19 ------------------- 1 file changed, 19 deletions(-) delete mode 100644 linear_programming/README.md diff --git a/linear_programming/README.md b/linear_programming/README.md deleted file mode 100644 index 4f520e173493..000000000000 --- a/linear_programming/README.md +++ /dev/null @@ -1,19 +0,0 @@ -# Linear Programming -Sub-field of mathematical optimisation, where the best solution is found while satisfying linear functions. - -Taking this linear program as an example -``` - Maximise P = 8x1 + 6x2 - subject to: x1 + x2 <= 20 - x1 >= 5 - x2 = 6 - with: - x1 >= 0, x2 >= 0 -``` -x1 could be the number of blue cars and x2 the number of red cars produced. The company would be trying to maximise profit P. - -Resources: - -https://en.wikipedia.org/wiki/Linear_programming -https://brilliant.org/wiki/linear-programming/ -https://byjus.com/maths/linear-programming/ \ No newline at end of file From aa1414504b82d94e6525936556a6a67e74c5c374 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu <113149540+imengus@users.noreply.github.com> Date: Sat, 17 Jun 2023 21:18:04 +0100 Subject: [PATCH 10/13] Update linear_programming/simplex.py Co-authored-by: Tianyi Zheng --- linear_programming/simplex.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index f59433989377..65dd7ecfdc41 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -33,8 +33,7 @@ def __init__(self, tableau: np.ndarray, n_vars: int) -> None: raise ValueError("RHS must be > 0") self.tableau = tableau - self.n_rows = len(tableau[:, 0]) - self.n_cols = len(tableau[0]) + self.n_rows, _ = tableau.shape # Number of decision variables x1, x2, x3... self.n_vars = n_vars From 7b6f8983fc60749ccdce9e21bc2abc32f6f8373f Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu <113149540+imengus@users.noreply.github.com> Date: Sat, 17 Jun 2023 21:21:05 +0100 Subject: [PATCH 11/13] Update linear_programming/simplex.py Co-authored-by: Tianyi Zheng --- linear_programming/simplex.py | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index 65dd7ecfdc41..c99b7bc686d7 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -246,9 +246,7 @@ def run_simplex(self) -> dict[Any, Any]: {'P': 7.0, 'x1': 5.0, 'x2': 2.0} """ # Stop simplex algorithm from cycling. - iter_num = 0 - - while iter_num < 100: + for _ in range(100): # Completion of each stage removes an objective. If both stages # are complete, then no objectives are left if not self.objectives: @@ -265,12 +263,8 @@ def run_simplex(self) -> dict[Any, Any]: if self.stop_iter: # Delete artificial variable columns and rows. Update attributes self.tableau = self.change_stage(self.tableau) - - # Pivot again - continue - - self.tableau = self.pivot(self.tableau, row_idx, col_idx) - iter_num += 1 + else: + self.tableau = self.pivot(self.tableau, row_idx, col_idx) return {} def interpret_tableau( From 781711ce2b7ee96200287d0ddb0806aeec7245fb Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu Date: Sat, 17 Jun 2023 21:30:56 +0100 Subject: [PATCH 12/13] fix class docstring --- linear_programming/simplex.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index f59433989377..9674b32c5ceb 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -20,8 +20,7 @@ class Tableau: """Operate on simplex tableaus - >>> tableau = np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]) - >>> t = Tableau(tableau, 2) + >>> t = Tableau(np.array([[-1,-1,0,0,-1],[1,3,1,0,4],[3,1,0,1,4.]]), 2) Traceback (most recent call last): ... ValueError: RHS must be > 0 From 7593cbbb713a53e17261e17dab682824e4dac701 Mon Sep 17 00:00:00 2001 From: Ilkin Mengusoglu Date: Sat, 17 Jun 2023 21:37:40 +0100 Subject: [PATCH 13/13] add comments --- linear_programming/simplex.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/linear_programming/simplex.py b/linear_programming/simplex.py index 539fc23980d8..ba64add40b5f 100644 --- a/linear_programming/simplex.py +++ b/linear_programming/simplex.py @@ -27,8 +27,8 @@ class Tableau: """ def __init__(self, tableau: np.ndarray, n_vars: int) -> None: - rhs = tableau[:, -1] - if np.any(rhs, where=rhs < 0): + # Check if RHS is negative + if np.any(tableau[:, -1], where=tableau[:, -1] < 0): raise ValueError("RHS must be > 0") self.tableau = tableau @@ -296,6 +296,7 @@ def interpret_tableau( if n_nonzero == nonzero_val == 1: rhs_val = tableau[nonzero_rowidx, -1] output_dict[col_titles[i]] = rhs_val + # Check for basic variables for title in col_titles: # Don't add RHS or slack variables to output dict