Блог

Я пытаюсь реализовать пересмотренный алгоритм симплексного метода (RSM) с использованием Python и numpy. Я застрял в том, что он либо работает только над максимизацией (корректно на крошечных матрицах, таких как 2×4, 5×5 и т. Д., И неправильно на больших), Либо вводит бесконечный цикл в большинстве случаев минимизации. Приведенный ниже код демонстрирует мою попытку реализовать минимизацию:

     def __init__(self, A: np.ndarray, b: np.ndarray, c: np.ndarray):
        base_size = A.shape[0]  # Number of basis vectors
        I = np.eye(base_size)  # Identity matrix, an initial basis 
        self.extended_matrix = np.concatenate((A, I), axis=1)  # Extended matrix of the system
        self.free_size = A.shape[1]  # Number of free vectors
        self.b = b  # Right parts of the constraints
        self.base_inv = I  # Initial inverse basis matrix
        self.c = np.concatenate((c, np.zeros(base_size)))  # Objective function quotients including those related to the basis variables
        self.c_i = [i for i, coeff in enumerate(self.c)]  # Indices for all variables
        
    @property
    def c_f_indices(self):
        """
        Indices of the free variables.
        """
        return self.c_i[:self.free_size]
    
    @property
    def c_T_B(self):
        """
        Objective function quotients related to the basis variables.
        """
        c_B_indices = self.c_i[self.free_size:]  # Basis variables indices.
        return self.c[c_B_indices]
    
    @property
    def c_T_f(self):
        """
        Objective function quotients related to the free variables.
        """
        return self.c[self.c_f_indices]
        
    @property
    def free(self):
        """
        Free vectors.
        """
        return self.extended_matrix[:, self.c_f_indices]
    
    @property
    def y_T(self):
        """
        Lagrange multipliers.
        """
        return self.c_T_B @ self.base_inv
    
    @property
    def deltas(self):
        """
        Net evaluations. 
        """
        return (self.y_T @ self.free) - self.c_T_f 
    


    def _swap(self, exits: int, enters: int) -> None:
        """
        In-/excluding respective vectors into/from the basis.
        """
        self.c_i[enters], self.c_i[exits   self.free_size] = self.c_i[exits   self.free_size], self.c_i[enters]
    
    def optimize(self):
        while any([delta > 0 for delta in self.deltas]): # < 0 in case of maximization
            x_B = self.base_inv @ self.b  # Current basis variables
            enters_base = np.argmax(self.deltas)  # Vector to enter the basis; argmin in case of maximization
            
            # Finding the vector to leave the basis:
            alpha = self.base_inv @ self.free[:, enters_base]

            try:
                exits_base = np.argmin([b/a if a > 0 else np.inf for b, a in zip(x_B, alpha)])
                assert alpha[exits_base] != 0
            except (ValueError, AssertionError):
                raise Exception("No solutions")
            
            # Finding the E_r matrix, which current inverse basis will be left-multiplied with in order to achieve the new inverse basis:
            zetas = [-alpha[j] / alpha[exits_base] if j != exits_base else 1/alpha[exits_base] for j, a_j in enumerate(alpha)]
            E = np.eye(self.free.shape[0])
            E[:, exits_base] = zetas
            self.base_inv = E @ self.base_inv
            
            # In-/excluding respective vectors into/from the basis:
            self._swap(exits_base, enters_base)
            
        return self.c_T_B @ self.base_inv @ self.b # Final objective function value
 

Я также пытался сортировать c_f_indices, но все равно получаю бесконечный цикл. Аналогичная реализация RSM также дает неправильные результаты на больших матрицах (например, 16×8, например) и отлично работает на крошечных.

Где корень проблемы?

Как уже упоминалось Эрвином Калвелагеном, обычное сводное правило Данцига может приводить к циклам и остановкам, когда в течение длительных периодов времени не происходит улучшения целевого значения.

В общем, это явление известно как вырождение LP. Ограниченная числовая точность и ошибки округления могут способствовать вырождению проблемы. Вот почему он обычно более распространен в больших LP.

Есть несколько способов борьбы с этой проблемой. Как упоминал Эрвин, чаще всего используется pertubation . Однако, если вы делаете это как учебный проект, я бы предложил использовать решение, в котором вы используете более усовершенствованное правило поворота, такое как правило Заде или Правило Каннингема, где вы просто сохраняете таблицу или переменную, чтобы убедиться, что вы выбираете другую переменную для ввода в базу после циклического выполнения.