MATLAB使用蒙特卡洛算法求解线性整数规划和0-1规划

蒙特卡洛方法是一种利用计算机的随机数理论模拟实际的情况的一种方法。今天主要是以实例讲解蒙特卡洛方法的MATLAB编程实现求解线性整数规划和0-1规划。

实例1

首先使用intlinprog线性整数规划求解函数对该线性规划进行求解，该函数的语法如下：

    x = intlinprog(f,intcon,A,b)
    x = intlinprog(f,intcon,A,b,Aeq,beq)
    x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)
    x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub,x0)
    x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub,x0,options)
    x = intlinprog(problem)
    [x,fval,exitflag,output] = intlinprog(___)
    %intlinprog函数标准型
    %x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)
    %f——系数阵
    %intcon——变量个数
    %A,b:不等式限制条件Ax<=b中的A和b
    %Aeq,beq:等式限制条件中的Aeq*x=beq中的Aeq和beq
    %lb,ub:自变量的最小值和最大值
    用［x , Fval］代替上述各命令行中左边的x，则可得到在最优解x处的函数值Fval
    [x,fval]=intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)

程序

clc;
clear all;
f=[-5 -8];
A=[1 1;5 9];
b=[6 45];
lb=zeros(2,1);
intcon=[1 2];
[x,fval]=intlinprog(f,intcon,A,b,[],[],lb,[]);
fprintf('max f(x) 在x1 = %f x2 = %f  处取得最大值：%f
',x(1),x(2),-fval);

运行结果

LP:                Optimal objective value is -41.250000.                                           


Heuristics:        Found 1 solution using ZI round.                                                 
                   Upper bound is -39.000000.                                                       
                   Relative gap is 2.50%.                                                          


Cut Generation:    Applied 1 Gomory cut.                                                            
                   Lower bound is -40.000000.                                                       
                   Relative gap is 0.00%.                                                          




Optimal solution found.


Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal
value, options.AbsoluteGapTolerance = 0 (the default value). The intcon variables are
integer within tolerance, options.IntegerTolerance = 1e-05 (the default value).


max f(x) 在x1 = 0.000000 x2 = 5.000000  处取得最大值：40.000000

蒙特卡洛求解线性整数规划程序

主程序(使用floor向下取整函数)

floor函数

floor - 朝负无穷大四舍五入
    此 MATLAB 函数 将 X 的每个元素四舍五入到小于或等于该元素的最接近整数。
    Y = floor(X)
    Y = floor(t)
    Y = floor(t,unit)

clc;
clear all;
rand('state',sum(clock));%初始化随机数发生器
f0=-inf;
x0 = [];
num = 1e7;
tic%计时开始
for i=1:num
    x=0 + 6*rand(2,1);%随机产生初始解
    x1 = floor(x);%向下取整函数
    [f,g]=mengte3(x1);%调用自定义函数计算
    if sum(g<=0)==2
        if f0<=f  %求最大值 如果当前值更优，则更新值
            x0=x1;
            f0=f;
        end
    end
end
toc%计时结束
fprintf('max f(x) 在x1 = %f x2 = %f  处取得最大值：%f
',x0(1),x0(2),f0);

自定义函数mengte3.m

function [f,g]=mengte3(x)
%% f是目标函数  g(x)<=0
f=5*x(1)+8*x(2);
x1 = x(1);
x2 = x(2);
g=[x1+x2-6;
  5*x1+9*x2-45];
end

运行结果

历时 10.114315 秒。
max f(x) 在x1 = 0.000000 x2 = 5.000000  处取得最大值：40.000000

使用ceil向上取整函数求解

ceil - 朝正无穷大四舍五入
    此 MATLAB 函数 将 X 的每个元素四舍五入到大于或等于该元素的最接近整数。
    Y = ceil(X)
    Y = ceil(t)
    Y = ceil(t,unit)

主程序

clc;
clear all;
rand('state',sum(clock));%初始化随机数发生器
f0=-inf;
x0 = [];
num = 1e7;
tic%计时开始
for i=1:num
    x=0 + 6*rand(2,1);%随机产生初始解
    x1 = ceil(x);
    [f,g]=mengte3(x1);%调用自定义函数计算
    if sum(g<=0)==2
        if f0<=f  %求最大值 如果当前值更优，则更新值
            x0=x1;
            f0=f;
        end
    end
end
toc%计时结束
fprintf('max f(x) 在x1 = %f x2 = %f  处取得最大值：%f
',x0(1),x0(2),f0);

运行结果

历时 11.714213 秒。
max f(x) 在x1 = 3.000000 x2 = 3.000000  处取得最大值：39.000000

实例2

首先使用intlinprog线性整数规划求解函数对该线性规划进行求解：

程序

clc;
clear all;
f=[-6 -2 -3 -5];
A=[-3 5 -1 -6;2 1 1 -1;1 2 4 5];
b=[-4 3 10]';
intcon=[1 2 3 4];
lb=zeros(4,1);
ub=ones(4,1);
[x,fval]=intlinprog(f,intcon,A,b,[],[],lb,ub);
fprintf('max f(x) 在x1 = %f x2 = %f x3 = %f x4 = %f  处取得最大值：%f
',x(1),x(2),x(3),x(4),-fval);

运行结果

LP:                Optimal objective value is -14.500000.                                           


Heuristics:        Found 1 solution using ZI round.                                                 
                   Upper bound is -13.000000.                                                       
                   Relative gap is 0.00%.                                                          




Optimal solution found.


Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal
value, options.AbsoluteGapTolerance = 0 (the default value). The intcon variables are
integer within tolerance, options.IntegerTolerance = 1e-05 (the default value).


max f(x) 在x1 = 1.000000 x2 = -0.000000 x3 = 1.000000 x4 = 1.000000  处取得最大值：14.000000
>>

蒙特卡洛求解0-1规划程序

主程序（四舍五入round函数）

clc;
clear all;
rand('state',sum(clock));%初始化随机数发生器
f0=-inf;
x0 = [];
num = 1e7;
tic%计时开始
for i=1:num
    x=0 + 1*rand(4,1);%随机产生初始解
    x1 = round(x);%利用四舍五入round函数0-1变量
    [f,g]=mengte4(x1);%调用自定义函数计算
    if sum(g<=0)==3
        if f0<=f  %求最大值 如果当前值更优，则更新值
            x0=x1;
            f0=f;
        end
    end
end
toc%计时结束
fprintf('max f(x) 在x1 = %f x2 = %f  x3 = %f x4 = %f  处取得最大值：%f
',x0(1),x0(2),x0(3),x0(4),f0);

自定义函数mengte4.m

function [f,g]=mengte4(x)
%% f是目标函数  g(x)<=0
f=6*x(1)+2*x(2)+3*x(3)+5*x(4);
x1 = x(1);
x2 = x(2);
x3 = x(3);
x4 = x(4);


g=[-3*x1+5*x2-x3-6*x4+4;
  2*x1+x2+x3-x4-3;
  x1+2*x2+4*x3+5*x4-10];
end

运行结果

历时 11.086100 秒。
max f(x) 在x1 = 1.000000 x2 = 0.000000  x3 = 1.000000 x4 = 1.000000  处取得最大值：14.0

本文内容来源于网络，仅供参考学习，如内容、图片有任何版权问题，请联系处理，24小时内删除。

作 者 | 郭志龙
编 辑 | 郭志龙
校 对 | 郭志龙