Yalmip使用教程(6)-将约束条件写成矩阵形式

        博客中所有内容均来源于自己学习过程中积累的经验以及对yalmip官方文档的翻译：https://yalmip.github.io/tutorials/

        这篇博客将详细介绍如何借助yalmip工具箱将约束条件写成矩阵形式。

1.相关函数介绍

1.1 depends和getvariables函数

        depends和getvariables函数都可用于求出sdpvar类型变量在Yalmip工具箱内部的索引(可以简单理解为该变量是第几个使用的变量，比如索引为3，表示为第3个使用的变量)，两个函数的使用语法分别为：

k = depends(x)

k = getvariables(x)

        其中，x为sdpvar类型变量，k为变量的索引。

        对于线性变量来说，两者的返回结果完全相同，例如下面的代码：

        例1：对线性变量使用depends和getvariables函数

yalmip('clear')
x = sdpvar(1);
y = sdpvar(2,3);
x_index1 = getvariables(x)
y_index1 = getvariables(y)
x_index2 = depends(x)
y_index2 = depends(y)

        运行结果为：

        两个函数的区别在于非线性的变量。对于非线性变量，depends函数只会返回其中涉及的线性变量，而getvariables函数则将直接返回非线性变量的索引。下面是一个例子：

        例2：对非线性变量使用depends和getvariables函数

yalmip('clear')
x = sdpvar(1);
z = x^2;
z_index1 = getvariables(z)
z_index2 = depends(z)

        运行结果为：

        上面的代码中涉及到1个线性变量x和一个非线性变量z，其中x的索引为1，z的索引为2。从结果中可以看到，getvariables函数将直接返回非线性变量z的索引，而depends函数只会输出非线性变量z中涉及到的线性变量x的索引。

1.2 getbase函数

        getbase函数将返回sdpvar类型变量中的full basis(我将其理解为变量的系数矩阵)，使用语法如下：

B = getbase(x);

        其中x表示sdpvar类型的变量，B表示所涉及基本变量(也就是Yalmip内部中分配了索引号的1维变量)的系数矩阵。光看这个函数的文字描述也是比较难理解它的作用，我们还是结合一个实例进行讲解。

例3：

yalmip('clear')
x = sdpvar(1);
y = sdpvar(1);
z = [1;2*x;3*y;4*x+5*y + 6*x^2];
full(getbase(z))

        运行结果为：

ans =
     1     0     0     0
     0     2     0     0
     0     0     3     0
     0     4     5     6

        首先我们可以知道，变量z为4×1的矩阵形式，其中共涉及3个基本变量，x，y和x²。那么我们可以把变量z用这三个基本变量进行表达，也就是：

        对于上面的例子来说：

        也就是系数矩阵为：

从分析结果可知，getbase函数所返回的就是所涉及基本变量的系数矩阵。

1.3 getbasematrix函数

        getbasematrix函数将返回sdpvar类型变量中指定基本变量的系数矩阵，使用语法如下：

B = getbasematrix(x,index)

        其中x表示sdpvar类型的变量，index表示指定的基本变量索引，B表示所涉及基本变量(也就是Yalmip内部中分配了索引号的1维变量)的系数矩阵。

        如果理解了getbase函数的用法，那么getbasematrix函数也就不难理解，无非是将返回所有涉及的基本变量系数矩阵变成了指定基本变量的系数矩阵。下面是一个示例：

        例4：

yalmip('clear')
x = sdpvar(1);
y = sdpvar(1);
z = [1;2*x;3*y;4*x+5*y + 6*x^2];
xB = full(getbasematrix(z,1))
yB = full(getbasematrix(z,2))

        运行结果如下：

1.4 recover函数

        recover函数通过索引来创建变量，标准语法如下：

x = recover(index)

        其中index表示变量索引，x表示所创建的变量。

        从本质上来讲，recover和getvariables函数是一个逆运算，即recover函数的输入是getvariables函数的输出，recover函数的输出是getvariables函数的输入。下面的代码说明了这一点：

yalmip('clear')
x = sdpvar(1);
assign(x,2);
x_index = getvariables(x);
x1 = recover(x_index);
x_value = value(x)
x1_value = value(x1)

        运行结果如下：

        从结果上来看，我们首先使用getvariables函数获取了变量x的索引，再利用recover函数，通过索引得到变量x1。由于我们提前使用了assign函数给变量x赋值为2，因此变量x的取值为2。但我们并没有给变量x1赋初值，但结果显示变量x1的取值也为2，说明变量x和变量x1本质上是同一个变量。也就是说通过recover函数和变量的索引号得到的新变量，和原变量是完全等价的。

1.5 see函数

        see函数在命令行返回sdpvar类型变量中所涉及基本变量的索引、常数矩阵和基本变量的系数矩阵，语法如下：

see(x)

        下面是一个例子：

        例5：

yalmip('clear')
x = sdpvar(1);
y = sdpvar(1);
z = [1;2*x;3*y;4*x+5*y + 6*x^2];
see(z)

        运行结果如下：

该结果和getvariables、getbase等函数返回结果的含义相同，此处不再赘述。

2.将优化问题改写为矩阵形式

        我们知道，在线性规划中一个变量对应优化问题中的一列，一个约束条件对应优化问题中的一行。对于有N个变量和M条约束的线性规划问题，那约束矩阵就有M行N列。我们可以将每一条约束都改写成≤0或≥0的形式，并使用一个中间变量表示这条约束，再通过上述函数确定中间变量所涉及的基本变量索引与系数矩阵，就可以间接表示出约束矩阵，进一步将优化问题改写为矩阵形式。

        下面通过1个实际例子进行说明。

        例6：

        假设我们想把上面的优化问题写成紧凑的矩阵形式：

        求矩阵A、b、c，然后使用矩阵形式求解优化问题的代码如下：

clc
clear
close all
warning off
yalmip('clear')
 
%% 决策变量
sdpvar x1 x2
 
%% 求系数矩阵
obj0 = x1 + 2*x2;
 
z = [-2*x1 + 3*x2 - 12 ;
    x1 + x2 - 14;
    -3*x1 + x2 + 3;
    3*x1 + x2 - 30];
 
M1 = full(getbase(z));
M2 = full(getbase(obj0));
index = depends([x1 x2]);
A = M1(:,index + 1);
b = -M1(:,1);
c = M2(index + 1)';
 
%% 矩阵形式的目标函数
x = [x1;x2];
obj = c'*x;
C = [A*x <= b , x >= 0];
 
%% 求解优化问题
ops = sdpsettings('verbose', 3, 'solver', 'gurobi');
sol = optimize(C , -obj ,ops);
 
%% 判断求解是否成功
if sol.problem == 0
    disp('求解成功!!!');
    x = value(x)
else
    disp(['求解失败，原因为',sol.info]);
end

        运行结果如下：

        上面的例子是一个简单的线性规划问题，直接写出约束矩阵也很容易，采用Yalmip函数反而多此一举。但是实际问题通常涉及0-1变量，多重下标，需要用循环语句表达约束等等，手动写系数矩阵非常麻烦且容易出错，借助Yalmip的相关函数将会简单很多。

3.测试题

3.1测试1

        将下列优化问题建模为混合整数线性规划问题，并利用Yalmip函数改写为矩阵形式：

3.2测试2

        其中数据为：

        请将上述优化问题改写为下面的矩阵形式，并使用矩阵形式进行求解：

3.3测试3

        将旅行商问题(TravelingSalesmanProblem，TSP)的0-1规划模型改写成矩阵形式并进行求解：

3.4测试题参考答案

        第六章测试题的参考答案可以从下面的链接中获取：