运筹学——matlab实现对偶单纯形法

前言

之前我已经发布了matlab实现单纯形法的相关内容，想了解的小伙伴可以点击这个链接：
https://blog.csdn.net/abbcdc/article/details/111324457

代码

这个代码仍然比较简单，适合初学者阅读，如果想要更好的学习效果建议先看看上面链接中使用matlab实现单纯形法的内容，再来阅读这个代码会觉得更轻松~

function [x_opt,fx_opt,iter] = DSimplex_eye(A,b,c)
% 对偶单纯形法求解标准形线性规划问题: max cx s.t. Ax = b x>=0
% 输入参数: c为目标函数系数, A为约束方程组系数矩阵, b为约束方程组常数项
% 输出参数: x_opt为最求解，fx_opt为最优函数值，iter为迭代次数
format rat                    %元素使用分数表示
[m,n] = size(A);              %m约束条件个数, n决策变量数
v=nchoosek(1:n,m);
ind_B=[];
for i=1:size(v,1)
    if A(:,v(i,:))==eye(m)
        ind_B=v(i,:);
    end
end
ind_N = setdiff(1:n, ind_B);  %非基变量的索引,原理是返回在1:n中出现而不在ind_B即基变量索引中出现的元素，并从小到大排序
iter=0;
while true
    x0 = zeros(n,1);
    x0(ind_B) = b;                %初始基可行解
    cB = c(ind_B);                %计算cB
    if ~any(b < 0)                %此基可行解为最优解, any存在某个<0        
        x_opt = x0;
        fx_opt = c'*x_opt;
        return
    end
    index=find(b<0);
    for i = 1:numel(index)
        if all(A(index(i),:)>=0)
            x_opt=[];
            fx_opt = [];         %此时原问题无可行解，对偶问题存在无界解
            return
        end
    end
    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)
    % 换基
    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);
    iter=iter+1;
end
end

测试样例

A=[-1 -2 -1 1 0;
      -2 1  -3 0 1]; 
b=[-3 -4]'; 
c=[-2 -3 -4 0 0]';
[x_opt,fx_opt,iter] = DSimplex_eye(A,b,c);
disp('最优解为：');
disp(x_opt);
disp('最优函数值为：');
disp(fx_opt);
disp('迭代次数为：');
disp(iter);

结果

最优解为：
      11/5     
       2/5     
       0       
       0       
       0       

最优函数值为：
     -28/5     

迭代次数为：
       2       

觉得有用的小伙伴请多多分享，谢谢~