企业服务网,网站建设和seo,wordpress可视化拖拽编辑,在网站做博客0.引言 上一篇博客介绍了使用Yalmip工具箱求解单阶段鲁棒优化的方法。这篇文章将和大家一起继续研究如何使用Yalmip工具箱求解两阶段鲁棒优化(默认看到这篇博客时已经有一定的基础了#xff0c;如果没有可以看看我专栏里的其他文章)。关于两阶段鲁棒优化与列与约束生成算法的原…0.引言 上一篇博客介绍了使用Yalmip工具箱求解单阶段鲁棒优化的方法。这篇文章将和大家一起继续研究如何使用Yalmip工具箱求解两阶段鲁棒优化(默认看到这篇博客时已经有一定的基础了如果没有可以看看我专栏里的其他文章)。关于两阶段鲁棒优化与列与约束生成算法的原理之前的博客已经详细地介绍过了这里就不再过多介绍主要是结合实例来讲解编程思路。这篇博客用到了两个算例1个是两阶段鲁棒优化问题和列与约束生成算法的开山鼻祖[1]另一个是电气专业中两阶段鲁棒优化问题最热门的文章之一[2]相信大家在网上见到过无数号称完美复现的代码但实际上大部分都是有问题的(包括我自己早期写的代码也是被网上的代码带歪了后面理解慢慢深入才发现问题所在)。 求解两阶段鲁棒优化问题一共有两个难点一是求解max-min或者min-max形式的子问题其实就是求解一个单阶段鲁棒优化上一篇博客我已经非常详细地介绍了求解方式借助Yalmip工具箱共有三种不同的方式可以解决。二是主问题和子问题的迭代求解也就是列与约束生成算法(CCG)的实现。很多代码在复现CCG算法时并没有向主问题同时添加列(变量)和约束这也是代码中最常见的问题。针对这两个难点我将用两个不同的算例详细地进行讲解。 此外文献[1]和[2]中都是采用了先将约束条件写成紧凑的矩阵形式然后再对子问题进行处理的方式很多朋友和我反映这部分太难处理了实际问题的约束建模过程中经常包括循环语句想要转成矩阵形式确实很不容易。这篇文章中我将分别采用两种不同的方式求解鲁棒优化。一是采用原始的约束条件省去将约束条件转为矩阵形式的步骤这种方式数学公式可能会更繁琐但比起矩阵形式的转换理解起来会更容易一些。二是采用矩阵形式进行编程在博客中我教大家一种非常简单就能将约束写为矩阵形式的方法文中只是介绍了如何使用之后也会单独写博客对此详细展开。 总之这篇博客干货满满可以认真通读一遍跟着博客中的思路亲自动手使用MatlabYalmip实现两阶段鲁棒优化的编程(博客中提到的所有例子我都提供了相应的代码)。相信大家理解后面对任何类型的两阶段鲁棒优化问题都能迅速使用类似的方法进行解决。 博客中主要包含8大内容 ①.拿到一个复杂的两阶段鲁棒优化问题的分析步骤和方法。 ②.采用Yalmip工具箱中的uncertain函数和鲁棒优化模块求解两阶段鲁棒优化的子问题。 ③.Yalmip工具箱中的鲁棒优化模块和常规的求解思路有什么异同。 ④.使用KKT条件求解两阶段鲁棒优化的子问题并使用CCG算法进行迭代求解。 ⑤.使用对偶变换求解两阶段鲁棒优化的子问题并使用CCG算法进行迭代求解。 ⑥.采用Yalmip工具箱的内置函数将线性约束写成紧凑矩阵形式的方法。 ⑦.矩阵形式的两阶段鲁棒优化问题如何快速写出子问题内层优化的KKT条件并使用CCG算法进行迭代求解。 ⑧.矩阵形式的两阶段鲁棒优化问题如何快速写出子问题内层优化的对偶问题并使用CCG算法进行迭代求解。 由于博客篇幅较长将分上下两篇发布其中上篇使用的是文献[1]中的算例包含上述①-⑤的内容。下篇使用文献[2]中的算例包含上述①、④、⑥-⑧的内容。 这篇博客是下篇的内容。
1.两阶段鲁棒优化基本形式 如文献[1]中所述标准的两阶段鲁棒优化问题的形式为 其中y为第一阶段决策变量u为不确定变量x为第二阶段决策变量。和分析单阶段鲁棒优化问题的五个特征一样拿到一个复杂的两阶段鲁棒优化问题先不用慌按照下面的步骤进行分析即可 1确定第一阶段决策变量有哪些将其与变量y对应。 2确定第二阶段决策变量有哪些将其与变量x对应。 3确定不确定变量有哪些将其与变量u对应。 4确定优化问题中不确定集合的形式并考虑是否可以直接使用Yalmip中的鲁棒优化模块进行求解。 5确定目标函数是否有仅包含第一阶段决策变量的项如果有的话可以单独拿出来。 6确定子问题的目标函数将其与鲁棒优化的标准形式相对应。 7确定约束条件考虑是否包含非线性约束是否需要线性化。 8求解max-min或者min-max类型的子问题。 9使用迭代方式将子问题产生的变量和约束不断添加到主问题中最终得到最优解。 下面分别以文献[1]和[2]中的优化问题进行讲解说明
2.两阶段鲁棒运输问题编程实战 更多内容请关注MatlabYalmip两阶段鲁棒优化通用编程指南(上)
鲁棒优化入门(6)—MatlabYalmip两阶段鲁棒优化通用编程指南(上)
3.微电网两阶段鲁棒优化调度编程实战 第二个算例来源于文献[2]其两阶段鲁棒优化问题的形式为 s.t. 我之前写过一篇博客解析这篇论文但是使用的是紧凑的矩阵形式。这次博客我先尝试使用一般形式的约束进行求解方便大家体会采用矩阵形式求解两阶段鲁棒优化的方便之处。 首先逐步进行分析 1确定第一阶段决策变量有哪些。 2确定第二阶段决策变量有哪些。 3确定不确定变量有哪些。 4确定优化问题中不确定集合的形式并考虑是否可以直接使用Yalmip中的鲁棒优化模块进行求解。 该优化问题的不确定集中为箱式不确定集但包含0-1变量不可以直接使用Yalmip中的鲁棒优化模块进行求解。 5确定目标函数是否有仅包含第一阶段决策变量的项如果有的话可以单独拿出来。 该优化问题第一阶段和第二阶段目标函数相同。 6确定子问题的目标函数将其与两阶段鲁棒优化的标准形式相对应。 7确定约束条件考虑是否包含非线性约束是否需要线性化。 子问题中的约束条件均为线性但不确定变量的约束条件包含0-1变量则含变量的优化问题并不满足强对偶定理KKT条件和对偶变换都无法适用。那是不是无法使用常规的方法求解单阶段鲁棒优化形式的子问题呢 当然不是如果出现这种情况我们就可以把不确定变量写在外层。那么对于子问题的内层优化问题来说变量是一个确定值即使其中包含0-1变量也只是一个常数。实际上子问题的内层优化并不包含0-1变量因此还是可以使用常规的KKT条件或强对偶定理进行求解。但是由于不确定变量中包含0-1变量直接通过uncertain函数直接使用Yalmip鲁棒优化模块的方法不可行。 分析完成后分别使用KKT条件和强对偶定理求解子问题具体如下
3.1 KKT条件求解子问题 为了方便求解我们首先把子问题的内层min优化问题写出来并将所有约束写成≤0的形式或0的形式 根据拉格朗日函数可以写出内层优化的KKT条件。 其中表示上三角全为1主对角线以下全为0的24阶方阵即 要想求解这个问题还需要写出引入16组0-1变量使用大M法将16组互补松弛条件转为线性约束具体如下 其中q1-q16都是0-1变量qi和λi的维度相同。 将内层优化的KKT方程组添加到外层优化中就可以将双层优化问题转为单层优化问题。经过上面的处理便顺利将max-min形式的子问题转为混合整数线性规划问题并可以使用Yalmip进行求解代码在压缩包中的Problem2文件夹中运行Problem2_subproblem_KKT1.m文件即可得到结果。 PS假设第一阶段变量的取值 运行结果如下 从运行结果可以看出使用KKT条件可以正常求出最优解但由于引入了非常多的拉格朗日乘子和相应的0-1变量对于初始的约束形式手动写KKT条件非常麻烦一旦出错也很难发现问题在哪。
3.2 将两阶段鲁棒优化问题写成矩阵形式 对于约束条件比较复杂的两阶段鲁棒优化问题先把目标函数和约束条件写成紧凑的矩阵形式再使用KKT条件或对偶变换求解会方便很多。文献[1]中将优化问题写成了紧凑的矩阵形式但存在一定问题我在此重新梳理一遍。 首先将优化问题写成紧凑的矩阵形式(等式约束都改写成不等式形式例如x0可以写成x≥0-x≥0) 其中第一阶段决策变量x是一个48维的列向量 第二阶段决策变量y是一个192维的列向量 c和y的维度一样但c是一个192维的常数列向量 那么如何确定矩阵G矩阵E矩阵M向量h的取值呢之前的博客中我采用了手动推导的方式(微电网两阶段鲁棒优化经济调度方法matlab代码_)。 这次我来教大家使用Yalmip工具箱中的函数depends、getbase、getbasematrix、see写出约束矩阵取值的方法函数的语法可以参考官方文档。 下面是一个简单的线性规划问题 我用这个例子来帮助大家体会depends、getbase、getbasematrix、see等函数的用法代码如下
clc
clear
close all
warning off
yalmip(clear)%% 决策变量
sdpvar x1 x2
%% 目标函数
obj x1 2*x2;z [-2*x1 3*x2 - 12 ;x1 x2 - 14;-3*x1 x2 3;3*x1 x2 - 30];x1_index depends(x1)
x2_index depends(x2)
M1 full(getbase(z))
M2 full(getbasematrix(z,depends(x1)))
M3 full(getbasematrix(z,depends(x2)))
see(z) 运行结果如下 观察可知depends函数返回的是变量在Yalmip工具箱中的编号getbase函数以稀疏矩阵的形式返回sdpvar变量中的常数项系数以及变量系数getbasematrix函数以稀疏矩阵的形式返回sdpvar变量中指定编号的变量系数see函数直接在命令行输出sdpvar变量中的常数项、变量系数以及用到的变量编号。 假设我们想把上面的优化问题写成紧凑的矩阵形式 求矩阵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 0disp(求解成功!!!);x value(x)
elsedisp([求解失败原因为,sol.info]);
end 运行结果如下 采用相同的方法将文献[1]所提确定性优化问题改写成矩阵形式求解的代码在压缩包中的Problem2文件夹中运行Problem2_matrix.m文件即可得到结果。 3.3 矩阵形式的CCG算法与KKT条件求解 原始子问题为 内层优化中变量x和不确定变量u都可以看作常数其拉格朗日函数可以写作 其中π和θ都是拉格朗日乘子。根据拉格朗日函数进一步写出KKT条件 这就是初始的KKT条件可以进一步化简。首先根据式(1)和(4)可以得到 然后根据式(1)和(3)可以得到 经过这样处理就可以消去变量θ减少决策变量的数目得到 其中包含两项互补松弛条件可以通过引入0-1变量使用大M法转换为线性约束 将KKT条件添加到子问题的外层优化中子问题便转为 子问题变成了一个混合整数线性规划可以采用求解器有效地进行求解。 再复习一下我们把文献[2]中所提两阶段鲁棒优化分解成一个主问题MP2_KKT和一个子问题SP2_KKT 主问题MP2_KKT 子问题SP2_KKT 使用KKT条件CCG算法求解该两阶段鲁棒优化问题的代码在压缩包中的Problem2文件夹中运行Problem2_KKT.m文件即可得到结果。
3.4 矩阵形式的CCG算法与对偶变换求解 原始子问题为 其内层优化min问题的对偶问题可以写做 将其与外层max问题合并得到 该问题的目标函数中包含的非线性项可以进行线性化也可以直接使用求解器进行求解。为了方便起见我直接使用了求解器进行求解。但是由于问题规模比较大求解时间会比较长(大概要3小时左右才能得到最优解)。 首先将不确定变量u写成 式中符号表示矩阵的Hadamard乘积(哈达玛积 Hadamard Product - 知乎)即对应元素相乘对应于matlab中的“.*”符号。在此基础上可以对原优化问题的目标函数进行化简 两阶段鲁棒优化问题可以分为主问题和子问题 对偶变换版本主问题MP2_dual 对偶变换版本子问题SP2_dual 使用CCG算法与对偶变换求解两阶段鲁棒优化的步骤概括如下 使用对偶变换CCG算法求解该两阶段鲁棒优化问题的代码在压缩包中的Problem2文件夹中运行Problem2_dual.m文件即可得到结果。 参考文献
[1]Zeng B, Zhao L. Solving two-stage robust optimization problems using a column-and-constraint generation method[J]. Operations Research Letters, 2013, 41(5): 457-461
[2]刘一欣,郭力,王成山.微电网两阶段鲁棒优化经济调度方法[J].中国电机工程学报,2018,38(14):4013-40224307.
PS 完整资料可以私信博主获取。
