@Ntonk

Почему зацикливается модифицированный симплекс-метод (реализация на Python с numpy)?

Пытаюсь реализовать модифицированный симплекс-метод в мультипликативной форме на Python с использованием numpy. Постоянно сталкиваюсь с проблемой того, что он у меня работает только на максимизацию, причем на маленьких матрицах - большие считает с ошибками, значения целевой функции сильно отличаются от тех, которыми они должны быть. На минимизации чаще всего зацикливается. Теория изложена, нарпимер, здесь. Вот пример моего кода для минимизации:

class Simplex(object):
    def __init__(self, A: np.ndarray, b: np.ndarray, c: np.ndarray):
        base_size = A.shape[0]  # Кол-во базисных векторов
        I = np.eye(base_size)  # Единичная матрица - начальный базис системы. 
        self.extended_matrix = np.concatenate((A, I), axis=1)  # Расширенная матрица системы
        self.free_size = A.shape[1]  # Кол-во свободных векторов
        self.b = b  # Правые части ограничений
        self.base_inv = I  # Начальная обратная матрица базиса равна его прямой матрице, т.к. является единичной
        self.c = np.concatenate((c, np.zeros(base_size)))  # Коэффициенты целевой ф-ции с учетом базисных переменных
        self.c_i = [i for i, coeff in enumerate(self.c)]  # Индексы всех переменных
        
    @property
    def c_f_indices(self):
        """
        Индексы свободных переменных. 
        """
        return self.c_i[:self.free_size]
    
    @property
    def c_T_B(self):
        """
        Коэффициенты целевой функции при базисных переменных.
        """
        c_B_indices = self.c_i[self.free_size:]  # Индексы базисных переменных.
        return self.c[c_B_indices]
    
    @property
    def c_T_f(self):
        """
        Коэффициенты целевой функции при свободных переменных.
        """
        return self.c[self.c_f_indices]
        
    @property
    def free(self):
        """
        Свободные векторы.
        """
        return self.extended_matrix[:, self.c_f_indices]
    
    @property
    def y_T(self):
        """
        Вектор относительных оценок (разрешающие множители по Канторовичу).
        """
        return self.c_T_B @ self.base_inv
    
    @property
    def deltas(self):
        """
        Текущие оценки векторов.
        """
        return (self.y_T @ self.free) - self.c_T_f 
    


    def _swap(self, exits: int, enters: int) -> None:
        """
        Включение/исключение в/из базис(-а) соответствующих векторов
        """
        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
            x_B = self.base_inv @ self.b  # Текущие значения базисных переменных
            enters_base = np.argmax(self.deltas)  # Нахождение вектора, включаемого в базис. В случае максимизации будет argmin
            
            # Нахождение вектора, исключаемого из базиса:
            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")
            
            # Нахождение матрицы соотношения E_r и вычисление новой обратной матрицы базиса системы:
            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
            
            # Включение/исключение в/из базис(-а) соответствующих векторов:
            self._swap(exits_base, enters_base)
            
        return self.c_T_B @ self.base_inv @ self.b # Итоговое значение целевой функции


Я пробовал, в том числе, сортировать индексы свободных переменных, но все равно приходил к зацикливанию. И, кстати, нашел чью-то похожую реализацию, которая тоже допускает грубые ошибки на больших матрицах. Две недели ломания головы ни к чему не привели. Где может крыться ошибка? Заранее благодарю.
  • Вопрос задан
  • 62 просмотра
Пригласить эксперта
Ваш ответ на вопрос

Войдите, чтобы написать ответ

Войти через центр авторизации
Похожие вопросы