Yalmip使用学习

yalmip学习

1. yalmip简介

1.1 什么是yalmip

yalmip是由Lofberg开发的一种免费的优化求解工具,其最大特色在于集成许多外部的最优化求解器,形成一种统一的建模求解语言,提供了Matlab的调用API,减少学习者学习成本。

1.2 yalmip安装方式

这里以MATLAB的安装方式为例,在官网上下载最新包,将其解压至matlabtoolbox文件夹下(当然也放置在其他文件夹),打开matlab软件添加Path路径即可。最后键入which sdpvar命令,显示sdpvar路径则安装成功。

2.yalmip求解优化问题的四部曲

2.1 创建决策变量

yalmip一共有三种方式创建决策变量,分别为:

  1. sdpvar-创建实数型决策变量
  2. intbar-创建整数型决策变量
  3. binvar-创建0/1型决策变量

不过值得注意的是,在创建n*n的决策变量时,yalmip默认是对称方阵,所以要创建非对称方针时,需要这样写:

xxxvar(n,n,'full')

2.2 添加约束条件

比起matlab自带的各种优化函数所要写明的约束条件,yalmip的约束条件写起来是非常舒适直观的。

比如要写入0<=x1+x2+x3<=1

那么可以这样写:

% 创建决策变量
x = sdpvar(1,3);
% 添加约束条件
C = [0<=x(1)+x(2)+x(3)<=1];

是不是非常爽呢。这才是人类语言(和我初见python的感觉差不多)

2.3 参数配置

关于参数设置,我们大多数是用来设置求解器solver的,当然还有其它的选项,可以通过doc sdpsettings查看。

2.4 求解问题

最后就是求解问题了。

首先要明确求解目标zyalmip默认是求解最小值问题,所以遇到求解最大值的问题,只需要在原问题的基础上添加一个负号即可。

求解调用格式:

optimize(target,constraints,opstions)

2.5 几个常用的其它指令

  1. check:可以检查约束条件是否被满足(检查约束条件的余值)
  2. value:可以查看变量或表达式的值
  3. assign: 可以给变量赋值,这个命令调试时很重要

3.举两个栗子

3.1 简单例子

如题

$$ \begin{equation} z = max(\frac{x_1+2x_2}{2x_1+x_2})\\ \left\{ \begin{array}{c} x_1 + x_2 \ge 2 \\ x_2 - x_1 \le 1 \\ x_1 \le 1 \end{array} \right. \end{equation} $$

代码如下,附有详细解释,就不说明了:

% 清除工作区
clear;clc;close all;
% 创建决策变量
x = sdpvar(1,2);
% 添加约束条件
C = [
    x(1) + x(2)  >= 2
    x(2)-x(1) <=1
    x(1)<=1
    ];
% 配置 找不到lpsolve求解器请检查是否安装,或者直接不适用'solver'字段。
% ops = sdpsettings('verbose',1);
ops = sdpsettings('verbose',1,'solver','lpsolve');
% 目标函数
z = -(x(1)+2*x(2))/(2*x(1)+x(2)); % 注意这是求解最大值
% 求解
result = optimize(C,z,ops);
if result.problem == 0 % problem =0 代表求解成功
    value(x)
    -value(z)   % 反转
else
    disp('求解出错');
end

求解结果:

ans =

   0.500000999998279   1.500000333331623


ans =

   1.399999360001426

3.2 解决经典的TSP问题

关于TSP的理论,这里我就不详细介绍了,百度有很多。在遇到yalmip之前,我学习的求解TSP的第一解法就是利用lingo来求解,后来学习了几种智能算法,如遗传算法,模拟退火,蚁群算法等等都可以解决这个问题。现在,学习了yalmip之后,我们可以完全抛弃lingo那种简陋的ide。废话不多说,先贴上约束条件:

$$ \begin{align*} &\min Z = \sum_{i=1}^{n}\sum_{j=1}^{n} d_{ij}x_{ij} \\ s.t. &\left\{ \begin{array}{c} \sum_{i=1,i\neq j}^{n}x_{ij} = 1,\qquad j = 1,\cdots,n \\ \sum_{j=1,j\neq i}^{n}x_{ij} = 1,\qquad i = 1,\cdots,n\\ u_i-u_j + nx_{ij} \le n-1,\qquad 1< i\neq j \le n \\ x_{ij} = 0 \text{或} 1,\qquad u,j=1,\cdots,n \\ u_i\text{为实数},\qquad i=1,\cdots,n \end{array} \right. \end{align*} $$

再来看看代码:

% 利用yamlip求解TSP问题
clear;clc;close all;
d = load('tsp_dist_matrix.txt')';
n = size(d,1);
% 决策变量
x = binvar(n,n,'full');
u = sdpvar(1,n);
% 目标
z = sum(sum(d.*x));
% 约束添加
C = [];
for j = 1:n
    s = sum(x(:,j))-x(j,j);
    C = [C,   s  == 1];
end
for i = 1:n
    s = sum(x(i,:)) - x(i,i);
    C = [C, s  == 1];
end
for i = 2:n
    for j = 2:n
        if i~=j
            C = [C,u(i)-u(j) + n*x(i,j)<=n-1];
        end
    end
end
% 参数设置
ops = sdpsettings('verbose',0);
% 求解
result= optimize(C,z,ops);
if result.problem == 0
    value(x)
    value(z)
else
    disp('求解过程中出错');
end

这里用到的tsp_dist_matrix.txt如下:

0 7 4 5 8 6 12 13 11 18
7 0  3 10 9 14 5 14 17 17
4 3 0 5 9 10 21 8 27 12
5 10 5 0 14 9 10 9 23 16
8 9 9 14 0 7 8 7 20 19
6 14 10 9 7 0 13 5 25 13
12 5 21 10 8 13 0 23 21 18
13 14 8 9 7 5 23 0 18 12
11 17 27 23 20 25 21 18 0 16
18 17 12 16 19 13 18 12 16 0

最后来看看结果吧:

>> value(x)

ans =

   NaN     0     0     0     0     0     0     0     1     0
     0   NaN     1     0     0     0     0     0     0     0
     0     0   NaN     1     0     0     0     0     0     0
     1     0     0   NaN     0     0     0     0     0     0
     0     0     0     0   NaN     0     1     0     0     0
     0     0     0     0     1   NaN     0     0     0     0
     0     1     0     0     0     0   NaN     0     0     0
     0     0     0     0     0     1     0   NaN     0     0
     0     0     0     0     0     0     0     0   NaN     1
     0     0     0     0     0     0     0     1     0   NaN

>> value(z)

ans =

    77

最后,发现了yalmip的一个bug,在书写yalmip的约束条件时,如下:
x(1) + x(2)-2 >= 0
注意x(2)-2这个-2是紧贴x(2)的。这样求解释正确的,如果换成
x(1) + x(2) -2 >= 0
x(2)-2这个-2不是紧贴x(2)的。这样求解则会出现问题。

建议把所有常数写在等式的右边

Last modification:October 20th, 2019 at 12:41 pm
如果觉得我的文章对你有用,请随意赞赏

Leave a Comment

54 comments

  1. lili

    博主,为什么第二个例子模型里的非回路约束和代码里的不一致呢, 模型里的约束是怎么转变成C = [C,u(i)-u(j) + n*x(i,j)<=n-1],u(i) u(j)是什么意思呢?谢谢

    1. Raven
      @lili

      这里应该是我的latex公式写错了,已改正,代码是没问题的。约束条件的推导,建议你百度下资料,应该挺多讲这个的。

      1. lili
        @Raven

        谢谢您,明白了

  2. wenshijie
    该评论仅登录用户及评论双方可见
    1. Raven
      @wenshijie

      你给optimize函数的参数应该有问题。

      1. wenshijie
        @Raven
        该评论仅登录用户及评论双方可见
  3. Xiao

    博主,问一下这个怎么求解线性矩阵不等式?找了很久没找到有关的教程

    1. Raven
      @Xiao

      举个例子呢。

  4. xlchen

    博主,请收下小弟的膝盖

  5. Lee

    博主,我的代码已经再次修改了,但是现在错误显示却仍是下面这样,我能把代码发到你邮箱,拜托您给我看看源代码,试试看在你的电脑上运行会是一样的问题吗?我现在太需要别人的建议与意见了,十分感谢

    1. Raven
      @Lee

      不好意思,最近一星期都在外地,得回家才有电脑调试了。

  6. Lee

    博主,我的代码已经修改的没有问题了,但是现在错误显示却是下面这样,能拜托您给我看看源代码然后指出错误吗

    Starting YALMIP integer branch & bound.Lower solver : fmincon-standardUpper solver : rounderMax time : 3600Max iterations : Inf

    Warning : The continuous relaxation may be nonconvex. This means
    that the branching process is not guaranteed to find a
    globally optimal solution, since the lower bound can be

    Hence, do not trust the bound or the gap...
    Node Upper Gap(%) Lower Open Elapsed time
    1 : Inf Inf NaN 0 2.1 Infeasible problem1 Finishing. Cost: Infwent wrong!
    1. Raven
      @Lee

      最近两天在外,不方便调试,不过你给的日志信息也指出问题了,The continuous relaxation may be nonconvex。应该是非凸问题无法找到全局最优。不过我也不是数学专业,不是很懂。

  7. Lee

    博主你好,我有个问题想请教一下,m(i,j)=(m(u,i)-r(i))*x(i,j)这样的一个约束,其中u,i,j均属于一个点集U,x定义为一个0-1变量,我在MATLAB中利用yalmip该怎么编程呢?

    1. Raven
      @Lee

      01变量用binvar即可。至于你说的u,i,j是一个点集,m(i,j)=(m(u,i)-r(i))*x(i,j)是你的约束条件是吧,那可以直接用循环去遍历整个点集然后加限制即可。

  8. liyaohong

    博主你好,yalmip能不能求解if else类的优化问题呢?

    1. Raven
      @liyaohong

      这样说其实我不太明白,有具体案例吗?

  9. liuliu

    就直接显示:求解过程中出错。没有其他报错~

    1. Raven
      @liuliu

      电脑有安装其它优化器吗,如lpsolver或者gurobi。Matlab自带的优化器估计是解不了的。

      1. liuliu
        @Raven

        楼主,您有没有求解过参数需要排序的01优化问题?

        1. Raven
          @liuliu

          参数排序?暂时还没做过。

      2. liuliu
        @Raven

        我只安装了yalmip

        1. Raven
          @liuliu

          对于大规模优化问题,matlab自带的优化器还不够,建议安装gurobi优化器。

          1. liuliu
            @Raven

            那您第二个例子,是只需要安装yalmip就可以求解出来嘛?还是说需要再用其他求解器,我目前是用MATLAB2017a版本,也只安装了yalmip

            1. Raven
              @liuliu

              经过测试,matlab自带的优化器也可解。我的matlab版本是2016b

              1. liuliu
                @Raven

                根据您的建议安装gurobi优化器后,现在可以运行了,还是非常感谢楼主~

            2. Raven
              @liuliu

              我当初是安装了gurobi后做的,matlab自带的我改没测试,晚点回去试试在回复你吧。

              1. liuliu
                @Raven

                好的,谢谢您了,方便加您个微信或qq嘛?

                1. Raven
                  @liuliu

                  问题不多的话就邮件沟通吧。

                  1. liuliu
                    @Raven

                    好的,麻烦您了~

  10. liuliu

    楼主,为啥第二个例子运行不出来结果,完全按照你的程序来的?

    1. Raven
      @liuliu

      把报错贴出来吧。