具有多个线性约束的二次规划
此示例显示了与默认的 quadprog 算法相比,'active-set' 'interior-point-convex' 算法在存在许多线性约束的情况下的表现如何。此外,'active-set' 算法中的拉格朗日乘数在非活动约束恰好为零,这在您寻找活动约束时非常有用。
问题描述
创建一个具有 N 变量和 10*N 线性不等式约束的伪随机随机二次问题。指定 N = 150。
rng default % For reproducibility N = 150; rng default A = randn([10*N,N]); b = 10*ones(size(A,1),1); f = sqrt(N)*rand(N,1); H = 18*eye(N) + randn(N); H = H + H';
检查得到的二次矩阵是否是凸的。
ee = min(eig(H))
ee = 3.6976
所有特征值都是正的,所以二次型 x'*H*x 是凸的。
不包括线性等式约束或边界。
Aeq = []; beq = []; lb = []; ub = [];
使用两种算法求解问题
设置选项以使用 quadprog 'active-set' 算法。该算法需要一个初始点。将初始点 x0 设置为长度为 N 的零向量。
opts = optimoptions('quadprog','Algorithm','active-set'); x0 = zeros(N,1);
对解进行计时。
tic [xa,fvala,eflaga,outputa,lambdaa] = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,opts);
Minimum found that satisfies the constraints. Optimization completed because the objective function is non-decreasing in feasible directions, to within the value of the optimality tolerance, and constraints are satisfied to within the value of the constraint tolerance. <stopping criteria details>
toc
Elapsed time is 0.042058 seconds.
将解时间与默认 'interior-point-convex' 算法的时间进行比较。
tic [xi,fvali,eflagi,outputi,lambdai] = quadprog(H,f,A,b);
Minimum found that satisfies the constraints. Optimization completed because the objective function is non-decreasing in feasible directions, to within the value of the optimality tolerance, and constraints are satisfied to within the value of the constraint tolerance. <stopping criteria details>
toc
Elapsed time is 2.305694 seconds.
对于具有许多线性约束的问题,'active-set' 算法的速度要快得多。
检查拉格朗日乘数
'active-set' 算法仅报告与线性约束矩阵相关的拉格朗日乘数结构中的几个非零项。
nnz(lambdaa.ineqlin)
ans = 14
相比之下,'interior-point-convex' 算法返回一个所有元素均为非零的拉格朗日乘数结构。
nnz(lambdai.ineqlin)
ans = 1500
几乎所有这些拉格朗日乘数的尺寸都小于 N*eps。
nnz(abs(lambdai.ineqlin) > N*eps)
ans = 20
换句话说,'active-set' 算法明确指示了拉格朗日乘数结构中的有效约束,而 'interior-point-convex' 算法则没有。