【运筹学】学习笔记--单纯形法Python实现（第2版）

OR研究生

于 2024-07-28 01:20:17 发布

阅读量209

点赞数 8

文章标签： python 学习笔记算法数据分析数学建模

本文链接：https://blog.csdn.net/TuTuTuhong/article/details/140744169

版权

1 利用Python单纯形法

1.1 创建单纯形表

import numpy as np
import sys
# 本例中需要采用标准型LP问题（添加人工变量、目标函数最小值、变量取值非负）
# 单纯形表前m行n列为系数表，m+1行为目标函数系数，n+1行为约束值，-(m+1,n+1)为目标函数值；
# 注：此案例中默认添加的人工变量系数矩阵为初始基矩阵；目前没有出基变量的索引（相当于没有解的索引）；
def create_tableau(A,b,c,num_constraints,num_variables):
    # np.hstack水平组合向量或矩阵；array.reshape(-1,1)将b转化为n行1列的矩阵；
    tableau=np.hstack((A, b.reshape(-1,1) ))
    # 需要增加判断，是求最小值还是最大值；
    c_row=np.hstack((c, np.zeros(1)))
    tableau=np.vstack((tableau,c_row))
    print('Iteration 0：\n',tableau)
    num_constraints=len(A)
    basic_vars = np.arange(num_variables, num_variables + num_constraints)
    non_basic_vars = np.arange(num_variables)
    print("Basic Variables:", basic_vars)
    print("Non-Basic Variables:", non_basic_vars)
    return tableau,basic_vars,non_basic_vars

1.2 定义换基过程

# 首先判断目标函数系数向量中非基变量的值是否非负，若全是非负数，则达到最优；若存在负数，则取最小值对应的变量作为入基变量；
# 原理是利用目标函数系数向量*变化量向量（改进方向/非基变量列），观察是否减小；
# 由于目标函数系数向量跟着变化（单纯形字典），用非基变量表示目标函数，而换基时是固定一个非基变量，变化一个非基变量，因此可以只看目标函数系数即可
def pivot(tableau,basic_vars,non_basic_vars):
    tableau_N=tableau[-1,non_basic_vars]
    arg_min_N=np.argmin(tableau_N)
    entering_vari_col=non_basic_vars[arg_min_N]
    if tableau[-1,entering_vari_col]>=0:
        print('当前解为最优解') 
        return False
    else:
        print(f'入基变量为x{entering_vari_col+1};',end='')
        step=np.where(tableau[:-1,entering_vari_col]>0)[0]
        if len(step)==0:
            print('此LP问题无界')
            sys.exit()
        else:
            #最小比值法确定步长和出基变量；存在疑问：基变量的变化量会为0吗？如果变化为0，则约束没变
            ratios=tableau[:-1,-1]/tableau[:-1,entering_vari_col]# b列是基变量的取值，始终非负
            ratios=np.where(ratios<0,np.inf,ratios)
            leaving_vari_row=np.argmin(ratios) #需要取正值最小值的索引
            print(f'出基变量为x{basic_vars[leaving_vari_row]+1};')
            tableau[leaving_vari_row,:]/=tableau[leaving_vari_row,entering_vari_col]
            for i in range(len(tableau)):
                if i!=leaving_vari_row:
                    tableau[i,:]-=tableau[i,entering_vari_col]*tableau[leaving_vari_row,:]
            print(tableau)
            convert=basic_vars[leaving_vari_row]
            basic_vars[leaving_vari_row]=non_basic_vars[arg_min_N]
            non_basic_vars[arg_min_N]=convert
            print("Basic Variables:", basic_vars)
            print("Non-Basic Variables:", non_basic_vars)
            return True

1.3 提取最优解和解的情况

def optimal(tableau,basic_vars,non_basic_vars):
    optimal_value = -tableau[-1, -1]  # tableau[-1, -1]是目标函数最小值的负数
    optimal_x=np.zeros(len(tableau[0])-1)
    for i,arg in enumerate(non_basic_vars):
        optimal_x[arg]=0
    for i,arg in enumerate(basic_vars):
        optimal_x[i]=tableau[arg,-1]
    return optimal_value,optimal_x

2 输入参数等

c = np.array([-12, -9, 0, 0, 0, 0])  # 目标函数的系数
A = np.array([[1, 0, 1, 0, 0, 0],
              [0, 1, 0, 1, 0, 0],
              [1, 1, 0, 0, 1, 0],
              [4, 2, 0, 0, 0, 1]])  # 约束矩阵
b = np.array([10, 15, 17.5, 48])  # 约束的右端项
num_variables=2 #需要手动输入下
num_constraints=len(A)
tableau,basic_vars,non_basic_vars=create_tableau(A,b,c,num_constraints,num_variables)
iteration=1
while True:
    print(f"\nIteration {iteration}:",end='')
    iteration+=1
    if pivot(tableau,basic_vars,non_basic_vars)==False:
        break

opt_value,solution=optimal(tableau,basic_vars,non_basic_vars)
print("\n最优解如下：")
print("目标函数最小值:", opt_value)
print("解的情况：",end='')
for i,val in enumerate(solution):
    print(f'x{i+1}={val}',end='; ')