Interior Point Method
主双対内点法:中心パスに沿った最適解へのアプローチ
1. イントロダクション
シンプレックス法が領域の境界(頂点)を探索するのに対し、内点法(Interior Point Method)は領域の内部を通り、最適解へ至る「中心パス(Central Path)」を追跡するアルゴリズムです。
大規模な線形計画問題においてはシンプレックス法を凌駕する性能を示すことから、現代の数理最適化ソルバーの多くで採用されています。本稿では、その中でも標準的な「主双対内点法(Primal-Dual Path Following Method)」について、数学的導出から具体的な数値計算の手順までを詳説します。
2. 理論的導出 (Derivation)
2.1 主問題と双対問題
標準形の線形計画問題とその双対問題を考えます。
Min $\boldsymbol{c}^T \boldsymbol{x}$
s.t. $\boldsymbol{Ax} = \boldsymbol{b}, \quad \boldsymbol{x} \ge 0$
Max $\boldsymbol{b}^T \boldsymbol{y}$
s.t. $\boldsymbol{A}^T \boldsymbol{y} + \boldsymbol{s} = \boldsymbol{c}, \quad \boldsymbol{s} \ge 0$
2.2 最適性条件 (KKT条件)
これらが最適解であるための必要十分条件は、以下の KKT (Karush-Kuhn-Tucker) 条件で記述されます。
- $\boldsymbol{A x} - \boldsymbol{b} = \boldsymbol{0}$ (主実行可能条件)
- $\boldsymbol{A}^T \boldsymbol{y} + \boldsymbol{s} - \boldsymbol{c} = \boldsymbol{0}$ (双対実行可能条件)
- $\boldsymbol{X S e} = \boldsymbol{0}$ (相補性条件:$x_i s_i = 0$)
- $\boldsymbol{x} \ge 0, \ \boldsymbol{s} \ge 0$
※ $\boldsymbol{X} = \text{diag}(x_1 \dots x_n)$, $\boldsymbol{S} = \text{diag}(s_1 \dots s_n)$, $\boldsymbol{e} = (1 \dots 1)^T$
2.3 バリアパラメータと中心パス
内点法では、相補性条件 $x_i s_i = 0$ を直接満たすのではなく、正のパラメータ $\mu$(バリアパラメータ)を用いて緩和します。
\[ x_i s_i = \mu \quad (\mu > 0) \]この方程式系の解 $(\boldsymbol{x}(\mu), \boldsymbol{y}(\mu), \boldsymbol{s}(\mu))$ が描く滑らかな曲線を中心パス(Central Path)と呼びます。アルゴリズムは、ニュートン法を用いてこのパスを追跡しながら、徐々に $\mu \to 0$ へと収束させます。
2.4 ニュートン方程式
現在地 $(\boldsymbol{x, y, s})$ からの修正量 $(\Delta \boldsymbol{x}, \Delta \boldsymbol{y}, \Delta \boldsymbol{s})$ を求めるため、KKT条件を線形近似して以下の巨大な連立一次方程式を解きます。
3. 具体的な計算手順 (Step-by-Step Calculation)
理論を理解するため、シンプレックス法のページと同じ問題を内点法で解いてみます。
Min $-4x_1 - 3x_2$
s.t. $x_1 + x_2 \le 40, \quad 2x_1 + x_2 \le 60$
スラック変数 $x_3, x_4$ を導入して等式標準形にします。
\[ \boldsymbol{A} = \begin{pmatrix} 1 & 1 & 1 & 0 \\ 2 & 1 & 0 & 1 \end{pmatrix}, \quad \boldsymbol{b} = \begin{pmatrix} 40 \\ 60 \end{pmatrix}, \quad \boldsymbol{c} = \begin{pmatrix} -4 \\ -3 \\ 0 \\ 0 \end{pmatrix} \]制約を満たす厳密に正の点(内点)を選びます。
- 主変数 $\boldsymbol{x}^{(0)}$: $(10, 10, 20, 30)^T$
- 双対変数 $\boldsymbol{y}^{(0)}$: $(-2, -2)^T$ (仮定)
- 双対スラック $\boldsymbol{s}^{(0)}$: $\boldsymbol{c} - \boldsymbol{A}^T \boldsymbol{y} = (2, 1, 2, 2)^T$ (すべて正なのでOK)
1. ギャップと $\mu$ の計算
現在の相補性ギャップ $\boldsymbol{x}^T \boldsymbol{s} = 10(2) + 10(1) + 20(2) + 30(2) = 130$
現在の $\mu = 130 / 4 = 32.5$
ターゲット $\mu_{target} = 0.5 \times 32.5 = 16.25$
2. ニュートン方程式の構築 ($10 \times 10$)
以下の連立方程式を解きます。
3. 解と更新
方程式を解くと以下の更新ベクトルが得られます。
$\Delta \boldsymbol{x} \approx (4.14, 1.76, -5.90, -10.04)^T$
これを用いて更新した次の点 $\boldsymbol{x}^{(1)}$ は、$(14.14, 11.76, \dots)$ となり、初期値 $(10, 10)$ から最適解 $(20, 20)$ の方向へ前進していることが確認できます。
このプロセスを数回繰り返すと $\mu \to 0$ となり、解は理論上の最適解 $(20, 20)$ へと滑らかに収束します。
4. シンプレックス法との比較
最後に、両アルゴリズムの挙動を比較します。
Fig 1. 探索ルートの違い(カクカク進むシンプレックス法 vs 滑らかに進む内点法)
| 特徴 | シンプレックス法 | 内点法 |
|---|---|---|
| 移動ルート | 領域の境界(頂点)を移動 | 領域の内部(中心パス)を移動 |
| 計算の核心 | ピボット操作(基底の入替) | ニュートン法(巨大な連立一次方程式) |
| 得意な問題 | 中規模、厳密解が必要な場合 | 超大規模問題(変数が数万以上) |