拉格朗日松弛简介
拉格朗日松弛（Lagrangian Relaxation）是一种约束优化问题里处理约束的思想。其将约束分为简单约束和困难约束，通过一个拉格朗日乘子将困难约束罚至目标函数上。拉格朗日松弛是1971年Held和Krap在研究旅行商问题（travelling salesman problem）时提出的。
假设我们考虑问题(不一定是凸问题)
其中(X)是(mathbb{R}^n)中的一个有限集。有时候，(g_i(m{x}) leq 0;i = 1,cdots,m)使求解非常困难，称为困难约束，而(m{x} in X)是容易求解的简单约束。拉格朗日松弛就是利用引入的拉格朗日乘子(lambda_i ge 0)将困难约束作为惩罚项写到目标函数里:
并将其称作拉格朗日函数，即(L(lambda) = minlimits_{m{x}} {f(m{x})+Sigma_{j=1}^m lambda_ig_i(m{x})|m{x} in X})。注意拉格朗日函数是关于(lambda)的分片线性函数，因为它是一系列仿射函数的逐点最小函数，所以是凹函数。我们称
为拉格朗日对偶(Lagrangian dual)，有些地方也称为拉格朗日乘子问题(Lagrangian multiplier problem)。困难约束可以是等式约束也可以是不等式约束。上述做法降低了求解的复杂程度，并使我们得到一个原问题最优解的下界。有时这个下界可以非常接近甚至等于最优解。
写在中间
拉格朗日松弛（Lagrangian Relaxation）是我们凸优化课程读书报告的主题，读书报告全文可见 LagrangianRelaxation.pdf。
读书报告的内容包括：
- 拉格朗日松弛的概念
- 困难约束为线性等式时
- 困难约束为线性不等式时
- 求解拉格朗日对偶：次梯度法（Subgradient）
- 例子：使用拉格朗日松弛求解约束最短路径问题
- 约束最短路径问题
- 一个具体的例子
- (L(lambda))中(lambda)取值的影响
- 求解(L^*)
报告里的内容就不在博客里再写一遍了（太麻烦了），大家可以直接去看全文 LagrangianRelaxation.pdf。这里只说一下关于程序实现的部分。
使用拉格朗日松弛求解约束优化问题的程序
强烈建议大家通读过全文之后再来看程序。
编程使用Matlab。这里求解凸优化问题使用了CVX工具包，安装和教程可见：CVX: Matlab Software for Disciplined Convex Programming。
程序1：指定乘子(lambda)，求解松弛后的问题
% Constrainted Shortest Path
% 在指定拉格朗日乘子lambda的前提下，
% 使用Lagrangian Relaxation 对 Constrainted Shortest Path 问题求解
clc
clear all
nInf = 10000;
% parametters
Cost = [nInf, 1, 10, nInf, nInf, nInf;
nInf, nInf, nInf, 1, 2, nInf;
nInf, 1, nInf, 5, 12, nInf;
nInf, nInf, nInf, nInf, 10, 1;
nInf, nInf, nInf, nInf, nInf, 2;
nInf, nInf, nInf, nInf, nInf, nInf];
Time = [0, 10, 3, nInf, nInf, nInf;
nInf, 0, nInf, 1, 3, nInf;
nInf, 2, 0, 7, 3, nInf;
nInf, nInf, nInf, 0,1, 7;
nInf, nInf, nInf, nInf, 0, 2;
nInf, nInf, nInf, nInf, nInf, 0];
n = 6; % number of nodes
T = 10; % max total time
lambda = 1; % lagrangian multiplier
cvx_begin % quiet
variable X(n, n)
minimize( sum(sum( Cost .* X) ) + lambda * (sum(sum(Time .* X)) - T) )
subject to
sum(X(1, :)) - sum(X(:, 1)) == 1;
for i = 2 : n-1
sum(X(i, :)) - sum(X(:, i)) == 0;
end
sum(X(n, :)) - sum(X(:, n)) == -1;
0 <= X <= 1; % linear relaxation for the integer constraints
cvx_end
程序2：使用次梯度法找到最佳的(lambda)值（求解对偶问题）
注意cvx_solver Mosek。为了获得更高的求解精度，这里使用了一个求解优化问题的商业求解器：Mosek。其和CVX配合使用的方法可见 Using MOSEK with CVX。
% Constrainted Shortest Path
% 使用 subgradient 方法找到最佳的 lambda 值
clc
clear all
nInf = 10000;
Cost = [nInf, 1, 10, nInf, nInf, nInf;
nInf, nInf, nInf, 1, 2, nInf;
nInf, 1, nInf, 5, 12, nInf;
nInf, nInf, nInf, nInf, 10, 1;
nInf, nInf, nInf, nInf, nInf, 2;
nInf, nInf, nInf, nInf, nInf, nInf];
Time = [0, 10, 3, nInf, nInf, nInf;
nInf, 0, nInf, 1, 3, nInf;
nInf, 2, 0, 7, 3, nInf;
nInf, nInf, nInf, 0,1, 7;
nInf, nInf, nInf, nInf, 0, 2;
nInf, nInf, nInf, nInf, nInf, 0];
n = 6;
T = 10;
% formulate the cost constraints as Ax = b
A = zeros(n, n * n);
for i = 1 : n
for j = 0 : n - 1
if j == i - 1
A(i, j*n + 1 : j*n + n) = - ones(1, n);
A(i, j*n + i) = 0;
else
A(i, j * n + i) = 1;
end
end
end
b = zeros(n, 1);
b(1, 1) = 1;
b(n, 1) = -1;
% subgradient iteration
k = 0;
lambda = 1;
while (1)
[x, L] = cvx_optimize(Cost, Time, n, T, A, b, lambda);
g = sum(sum(Time(:) .* x)) - T;
fprintf('k: %i, lambda: %f, L: %f, g: %f
', k, lambda, L, g);
if g <= 0 && lambda * g == 0
break;
end
theta = 0.1 * (20 - L) / (norm(g) ^ 2);
lambda = max(0, lambda + theta' * g);
k = k + 1;
end
function [x, fval] = cvx_optimize(Cost, Time, n, T, A, b, lambda)
cvx_solver Mosek
cvx_begin quiet
variable x(n * n)
minimize( sum(sum( Cost(:) .* x) ) + lambda * (sum(sum(Time(:) .* x)) - T) )
subject to
A * x <= b;
0 <= x <= 1
cvx_end
fval = sum(sum( Cost(:) .* x) ) + lambda * (sum(sum(Time(:) .* x)) - T);
end