二次計画問題から 2 次錐問題への変換

この例では、半正定値二次計画問題を coneprog ソルバーで使用される 2 次錐形式に変換する方法を示します。二次計画問題の形式は次のとおりです。

$$\underset{x}{\min }\frac{1}{2}{x}^{T}Hx+f^{T}x$$,

範囲と線形制約の影響を受ける可能性があります。coneprog は次の形式の問題を解きます。

$$\underset{x}{\min }{f}_{sc}^{T}x$$

条件

$$\Vert A_{sc}x-b_{sc}\Vert \le d_{sc}^{T}x-\gamma$$,

範囲と線形制約の影響を受ける可能性があります。

二次計画を coneprog 形式に変換するには、まず、行列 $$H$$ の行列の平方根を計算します。$$H$$ が対称半正定値行列であるとすると、コマンド

A = sqrtm(H);

は、A'*A = A*A = H となる半正定値行列 A を返します。したがって次のようになります。

$$x^{T}Hx=x^T{A}^{T}Ax=(Ax{)}^{T}Ax=\Vert Ax{\Vert }^2$$.

二次計画の形式を次のように変更します。

$$\underset{x}{\min }\frac{1}{2}{x}^{T}Hx+f^{T}x=-\frac{1}{2}+\underset{x,t}{\min }[(t+1/2)+f^{T}x]$$

ここで、$$t$$ は次の制約を満たします。

$$t+1/2\ge \frac{1}{2}{x}^{T}Hx$$.

制御変数 $$x$$ を $$u$$ に拡張します。これには、$$t$$ が最後の要素として含まれています。

$\mathit{u}=\left[\begin{array}{c}\mathit{x}\\ \mathit{t}\end{array}\right]$.

2 次錐制約行列とベクトルを次のように拡張します。

$${\mathit{A}}_{\mathrm{sc}}=\left[\begin{array}{cc}\mathit{A}& 0\\ 0& 1\end{array}\right]$$

$${\mathit{b}}_{\mathrm{sc}}=\left[\begin{array}{c}0\\ ⋮\\ 0\end{array}\right]$$

$${\mathit{d}}_{\mathrm{sc}}=\left[\begin{array}{c}0\\ ⋮\\ 0\\ 1\end{array}\right]$$

$\gamma =-1$.

係数ベクトル $$f$$ も拡張します。

${\mathit{f}}_{\mathrm{sc}}=\left[\begin{array}{c}\mathit{f}\\ 1\end{array}\right]$.

新しい変数に関して、二次計画問題は次のようになります。

$$\underset{u}{\min }\frac{1}{2}{u}^T{A}_{sc}^{2}u+f_{sc}^{T}u=-1/2+\underset{u}{\min }[1/2+f_{sc}^{T}u]$$

ここで、

$$u(end)+1/2\ge \frac{1}{2}{u}^T{A}_{sc}u$$.

この 2 次制約は、$$A_{sc}$$、$$d_{sc}$$、および $$\gamma$$ の以前の定義を使用する次の計算を通して錐制約になります。

$$\frac{1}{2}\Vert Hx{\Vert }^2=\frac{1}{2}{u}^T{A}_{sc}^{2}u\le t+\frac{1}{2}$$

$$\Vert Hx{\Vert }^2\le 2t+1$$

$$\Vert Hx{\Vert }^2+t^2\le t^2+2t+1=(t+1{)}^2$$

$$\sqrt{\Vert Hx{\Vert }^2+t^2}\le |t+1|=\pm (t+1)$$

ただし、$$\Vert Hx{\Vert }^2+t^2=\Vert A_{sc}u{\Vert }^2$$ です。したがって、$$\gamma =-1$$ および $$d_{sc}$$ の定義より、不等式は次のようになります。

$$\Vert A_{sc}u\Vert \le \pm (t-\gamma )$$

$$\Vert A_{sc}u\Vert \le \pm (d_{sc}^{T}u-\gamma )$$.

二次計画は、対応する錐計画と同じ解をもちます。唯一の違いは、錐計画に追加された項 $$-1/2$$ です。

数値例

quadprog ドキュメンテーションでこの例が紹介されています。

H = [1,-1,1
    -1,2,-2
    1,-2,4];
f = [-7;-12;-15];
lb = zeros(3,1);
ub = ones(size(lb));
Aineq = [1,1,1];
bineq = 3;
[xqp fqp] = quadprog(H,f,Aineq,bineq,[],[],lb,ub)

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.

xqp = 3×1

    1.0000
    1.0000
    1.0000

fqp = -32.5000

この例の冒頭の説明を参照しながら、2 次錐制約変数を指定してから、関数 coneprog を呼び出します。

Asc = sqrtm(H);
Asc((end+1),(end+1)) = 1;
d = [zeros(size(f(:)));1];
gamma = -1;
b = zeros(size(d));
qp = secondordercone(Asc,b,d,gamma);
Aq = Aineq;
Aq(:,(end+1)) = 0;
lb(end+1) = -Inf;
ub(end+1) = Inf;
[u,fval,eflag] = coneprog([f(:);1],qp,Aq,bineq,[],[],lb,ub)

Optimal solution found.

u = 4×1

    1.0000
    1.0000
    1.0000
    1.0000

fval = -33.0000

eflag = 1

錐の解 u の最初の 3 つの要素は、表示精度に対する二次計画法の解 xqp の要素と同じです。

disp([xqp,u(1:(end-1))])

    1.0000    1.0000
    1.0000    1.0000
    1.0000    1.0000

返される二次関数値 fqp は、$$2t+1$$ が正の場合は返される錐の値 - 1/2 に、$$2t+1$$ が負の場合は返される錐の値 + 1/2 になります。

disp([fqp-sign(2*u(end)+1)*1/2 fval])

  -33.0000  -33.0000

参考

quadprog | coneprog | secondordercone

関連するトピック

• 2 次制約から 2 次錐制約への変換
• 錐計画法を使用した区分線形ばね質量系のエネルギーの最小化 (ソルバーベース)