matlab setdiff函数_【运筹学】对偶单纯形法的Matlab实现

前两篇讨论了单纯形法的Matlab实现：

张敬信：【运筹学】单纯形法求解线性规划问题的Matlab实现​zhuanlan.zhihu.com

张敬信：【运筹学】单纯形法之大M法和两阶段法​zhuanlan.zhihu.com

本篇讨论对偶单纯形法的Matlab实现。

一. 对偶单纯形法的原理和步骤

原始问题的最优基也是对偶问题的可行基。换句话说，当原始问题的基 B 既是原始可行基又是对偶可行基时，B 成为最优基。

对偶单纯形法思想：

对偶单纯形法的步骤：

二. 对偶单纯形法的 Matlab 实现

算法步骤与单纯形法基本是类似的，除了进出基的判断规则有变化，所以，对 SimplexMax 函数稍加修改即可。

优点：

1. 代码采用向量化Matlab语言实现，非常简洁；
2. 可以把整个单纯形表直接保存下来，再输出到Excel非常方便。

缺点（欢迎指正）：

1. 各种特殊解的情况，暂时没有讨论全面；
2. 测

function [x,z,ST,res_case] = DualSimplexMax(c,A,b,ind_B)
% 对偶单纯形法求解标准形线性规划问题: max cx s.t. Ax <= b x>=0
% 输入参数: c为目标函数系数, A为约束方程组系数矩阵, b为约束方程组常数项, ind_B为基变量索引
% 输出参数: x最优解, z最优目标函数值, ST存储单纯形表数据, res_case=0表示有最优解

[m,n] = size(A);              %m约束条件个数, n决策变量数
ind_N = setdiff(1:n, ind_B);  %非基变量的索引
ST = [];
format rat
% 循环求解
while true
    x0 = zeros(n,1);
    x0(ind_B) = b;                %初始基可行解
    cB = c(ind_B);                %计算cB
    Sigma = zeros(1,n);
    Sigma(ind_N) = c(ind_N) - cB*A(:,ind_N);   %计算检验数
    [~,q] = min(b);               %选出最小的b
    r = ind_B(q);                 %确定出基变量索引r
    Theta = Sigma ./ A(q,:);      %计算θ
    Theta(Theta<=0) = 10000;
    [~,s] = min(Theta);           %确定进基变量索引s, 主元为A(q,s)
    vals = [cB',ind_B',b,A];
    vals = [vals; NaN, NaN, NaN, Sigma];
    ST = [ST; vals];
    if ~any(b < 0)               %此基可行解为最优解, any存在某个>0        
        x = x0;
        z = c * x;
        res_case = 0;
        return
    end
    % 换基
    ind_B(ind_B == r) = s;      %新的基变量索引
    ind_N = setdiff(1:n, ind_B); %非基变量索引
    % 更新A和b
    A(:,ind_N) = A(:,ind_B)  A(:,ind_N);
    b = A(:,ind_B)  b;
    A(:,ind_B) = eye(m,m);
end

测试实例：

改写为目标最大化，并化为标准型：

A = [-2 -4 0 1 0; -2 0 -5 0 1];
b = [-2; -3];
c = [-12 -16 -15 0 0];
ind = [4 5];
[x, z, ST, ca] = DualSimplexMax(c, A, b, ind)

运行结果：

也可以进一步将对偶单纯形表 ST 写入Excel（略）。

主要参考文献：

胡运权，运筹学原理及应用（第六版），及其课件

对偶单纯形法详解，百度文库

[图文]对偶单纯形法详解 - 百度文库​wenku.baidu.com

————————————————————————————————

原创作品，转载请注明，版权所有，禁止盗用。