簡単なSOCP例

このデモンストレーションは、nAG Library for Pythonに同梱されている例のセットから取られています。コマンドシェルで以下のコマンドを実行すると、すべての例を表示できます。

python -m naginterfaces.library.examples --locate

nAGのhandle_solve_socp_ipm関数（e04ptとしても知られる）は、大規模な二次錐計画問題（SOCP）のためのnAG最適化モデリングスイートのソルバーです。これは内点法（IPM）に基づいています。

$$\begin{array}{ll}\underset{x\in {\mathbb{R}}^n}{minimize}& c^{T}x\\ \text{subject to}& l_A\le Ax\le u_A\text{,}\\ & l_x\le x\le u_x\text{,}\\ & x\in \mathcal{K}\text{,}\end{array}$$

ここで、$\mathcal{K}={\mathcal{K}}^{n_1}\times \cdots \times {\mathcal{K}}^{n_r}\times {\mathbb{R}}^{n_l}$は二次（二次錐型）錐と$n_l$次元実空間のデカルト積であり、$n=\sum _{i=1}^r{n}_i+n_l$は決定変数の数です。ここで$c$、$x$、$l_x$、$u_x$は$n$次元ベクトルです。

$A$は$m$行$n$列の疎行列で、$l_A$と$u_A$は$m$次元ベクトルです。$x\in \mathcal{K}$は変数のサブセットを二次錐に分割し、各${\mathcal{K}}^{n_i}$は二次錐または回転二次錐のいずれかであることに注意してください。これらは以下のように定義されます：

• 二次錐：

$${\mathcal{K}}_q^{n_i}:=\left\{z=\left\{z_1,z_2,\dots ,z_{n_i}\right\}\in {\mathbb{R}}^{n_i}:z_1^2\ge \sum _{j=2}^{n_i}{z}_j^2,z_1\ge 0\right\}\text{.}$$

• 回転二次錐：

$${\mathcal{K}}_r^{n_i}:=\left\{z=\left\{z_1,z_2,\dots ,z_{n_i}\right\}\in {\mathbb{R}}^{n_i}:2z_1{z}_2\ge \sum _{j=3}^{n_i}{z}_j^2,z_1\ge 0,z_2\ge 0\right\}\text{.}$$

このルーチンの完全な説明については、nAGライブラリマニュアルのe04ptcを参照してください。

この例は、handle_set_group関数のドキュメントから派生したもので、以下のSOCP問題を解きます：

minimize $$10.0x_1+20.0x_2+x_3$$

from naginterfaces.base import utils
from naginterfaces.library import opt

# 問題のサイズ:
n = 3

# 問題ハンドルを作成する:
handle = opt.handle_init(nvar=n)

# 目的関数を設定する
opt.handle_set_linobj(handle, cvec=[10.0, 20.0, 1.0])

以下の制約条件に従う $$\begin{array}{rllll}-2.0& \le & x_1& \le & 2.0\\ -2.0& \le & x_2& \le & 2.0\end{array}$$

# 箱型制約を設定する
opt.handle_set_simplebounds(
    handle,
    bl=[-2.0, -2.0, -1.e20],
    bu=[2.0, 2.0, 1.e20]
)

一般線形制約

# 線形制約を設定する
_ = opt.handle_set_linconstr(
    handle,
    bl=[-1.e20, 1.0],
    bu=[1.5, 1.e20],
    irowb=[1, 1, 1, 2, 2, 2],
    icolb=[1, 2, 3, 1, 2, 3],
    b=[-0.1, -0.1, 1.0, -0.06, 1.0, 1.0]
    )

そして円錐制約

$$\left(x_3,x_1,x_2\right)\in {\mathcal{K}}_q^3\text{.}$$

# 円錐制約を設定
_ = opt.handle_set_group(
    handle,
    gtype='Q',
    group=[ 3,1, 2],
    idgroup=0
)

アルゴリズムのオプションをいくつか設定します

# アルゴリズムのオプションをいくつか設定します。
for option in [
        'Print Options = NO',
        'Print Level = 1'
]:
    opt.handle_opt_set(handle, option)

# 省略された反復出力のために明示的なI/Oマネージャーを使用する：
iom = utils.FileObjManager(locus_in_output=False)

最後に、ソルバーを呼び出します

# SOCP内点法ソルバーを呼び出す
result = opt.handle_solve_socp_ipm(handle, io_manager=iom)

 E04PT, Interior point method for SOCP problems
 Status: converged, an optimal solution found
 Final primal objective value -1.951817E+01
 Final dual objective value   -1.951817E+01

最適解は以下の通りです

result.x

array([-1.26819151, -0.4084294 ,  1.3323379 ])

そして目的関数の値は

result.rinfo[0]

-19.51816515094211

最後に、ハンドルを破棄して後片付けをします

# ハンドルを破棄する:
opt.handle_free(handle)