Research Reports on Mathematical and Computing Sciences Series B : Operations Research

Department of Mathematical and Computing Sciences Tokyo Institute of Technology 2-12-1 Oh-Okayama, Meguro-ku, Tokyo 152, Japan

Semidefinite Programming Relaxation and Global Optimization^{†} (in Japanese)

Masakazu Kojima^{‡}

B-342 July 1998

Abstract. We consider a nonconvex quadratic program having a linear objective function and finitely many quadratic inequalty constraints. We present implementable variants of the SSDP (Successive Semidefinite Programming) Relaxation Method and the SSILP (Successive Semi-Infinite LP) Relaxation Method recently proposed by Kojima and Tunçel [8, 9]. Each iteration of the variant of the SSDP Relaxation Method is carried out by solving a finite number of semidefinite programs, while each iteration of the SSILP Relaxation Method is carried out by solving a finite number of linear programs. Given any ϵ > 0, both variants compute an upper bound of the optimal objective value within ϵ in a finite number of iterations.

Key words.

Global Optimization, Nonconvex Quadratic Program, Relaxation, Linear Programming, Semidefinite Programming, Linear Matrix Inequality, Interior-Point Methods

† This article is written for a talk given in the 10th RAMP Symposium held in Kyoto, September 1998.

‡ Department of Mathematical and Computing Sciences, Tokyo Institute of Technology, 2-12-1 Oh-Okayama, Meguro-ku, Tokyo 152, Japan.

半正定値計画緩和と大域的最適化^{1}

小島政和

東京工業大学情報理工学研究科数理・計算科学専攻〒 152-8552 東京都目黒区大岡山 2-12-1

                          kojima@is.titech.ac.jp

1 はじめに

ここ１５年ほどの間に，線形計画問題を含む凸最適化分野では非常に大きな発展があった．その発展は１９８４年の Karmarkar 法 [7] を原点とした内点法に支えられている．特に，主・双対内点法（[15] 等参照）は超大規模な線形計画問題を解く実用的な計算手法として定着し，半正定値計画問題を含む凸計画問題へも拡張され，さらに，それらが実用化されつつある．理論面では Nesterov-Nemirovskii [12] が提案した self-concordance 理論がその中核をなしていると言ってよい．この分野での発展は凸性を仮定しない非凸最適化（大域的最適化）分野にも大きな影響を与えている．特に，半正定値計画緩和の 0-1 整数計画問題への適用は理論家の間で注目をあびている（[5, 6, 10, 11] 等参照）．しかしながら，凸最適化問題と比較して，理論面のからも手法の面からも非凸最適化問題は極めて難しい．凸最適化分野での双対理論， self-concordance 理論，主・双対内点法のように，非凸最適化分野全体を統一的に扱う理論や手法も見あたらない．本稿では， “非凸最適化問題に対する大域的収束理論” をめざして，著者等が最近の論文 [8, 9] で提案した逐次半正定値計画緩和法（SSDP Relaxation Method）および逐次半無限線形計画緩和法（SSILP Relaxation Method）を紹介する．この２つの手法は Lovász-Schrijver [10] による 0-1 整数計画問題に対する lift-and-project procedure の一般の非凸型２次計画問題への拡張になっている．

一般的な数理計画問題から始めよう．

\begin{matrix}目的 : & c^{T} x → 最大化 \\ 条件 : & x ∈ F \end{matrix} \}  (1)

ただし， c は n-次元 Euclid 空間 R^{n} の非ゼロ定数（縦）ベクトル， x ∈ R^{n} は変数ベクトル， F ⊆ R^{n} は許容領域， T はベクトルまたは行列の転置を表わす．この問題では線形の目的関数 c^{T} x が設定されているが，非線形の目的関数 f(x) の場合には，条件に x_{n+1} − f(x) ≤ 0 を加えることによって， f(x) を線形の目的関数 x_{n+1} に置き換えることが出来る．したがって， (1) は n-次元 Euclid 空間で定義された最も一般的な数理計画問題と考えてよい．問題 (1) に対して次の条件を仮定する．

条件 1.1. 許容領域 F は有界閉集合（空集合であってもよい）．この条件は，数学的には対象とする問題のクラスを狭めることになるが，実用的には十分に許容できる条件であり，また，一般的な数理計画問題 (1) の大域的最大解を計算するための必要最小限の仮定と考えてもよいであろう（もちろん，これに異論のある人もいるかもしれない）．この条件の下では，許容領域 F が空集合で無ければ，問題 (1) の大域的最大解，すなわち， c^{T} x^{∗} ≥ c^{T} x (∀x ∈ F) なる x^{∗} ∈ F が存在する．ここでの目標は

1 1998 年 9 月に京都で行われる第 10 回 RAMP シンポジウムの予稿集に掲載予定

• F = ∅ であることを示すか，あるいは，

• 問題 (1) の大域的最大解 x^{∗} を計算せよ

となる．

一般に，上記の目標を “厳密に達成する” のは極めて困難であろう．現実的には，さらに仮定を置き，問題のクラスを狭め，かつ，目標を “近似的に達成する” のが精一杯かもしれない．本稿では，論文 [8, 9] に基づき，半正定値計画緩和および半無限線形計画緩和がこの目標にどのように関わり，有用であるかを解説する．半正定値計画緩和は Lovász-Schrijver [10] によって 0-1 整数計画問題に適用された緩和であり，その後，非凸型２次計画問題にも拡張されている（[3, 5, 6, 11, 14] 等参照）．論文 [10] では，半正定値計画緩和をさらに緩和した半無限線形計画緩和も提案されている（ [1, 13] 等参照）．また，これらの手法は “lift-and-project procedure” とも呼ばれている．

第２節では問題 (1) の緩和について述べる．第３節では問題 (1) を非凸型２次計画問題に特殊化する．第４節では非凸型２次計画問題に対する逐次半正定値計画緩和法（SSDP Relaxation Method）について説明する．第５節ではこの方法の収束および近似に関する結果を与える．第６節では逐次半正定値計画緩和法をさらに緩和した逐次半無限線形計画緩和法（SSILP Relaxation Method）とその収束および近似に関する結果について述べる．

2 緩和問題

問題 (1) の許容領域 F を凸集合で緩和した問題

\begin{matrix}目的 : & c^{T} x → 最大化 \\ 条件 : & x ∈ C \end{matrix} \}  (2)

を導入する．ただし， C は F を含む任意の有界閉凸集合とする．このとき，以下の性質が成り立つことが知られている．

(i) 緩和問題 (2) の最大値 ζ は問題 (1) の最大値 ζ^{∗} の上界．

(ii) (2) の最大解 x^{∗} が (1) の許容領域 F に含まれれば， x^{∗} は (1) の最大解．

さらに， (1) と (2) の共通の目的関数 c^{T} x が線形であることを考慮にいれると， C を F の凸包 c.hull(F) （F を含む最小の凸集合）にとった場合には以下の性質が導かれる．

(i)’ ζ^{∗} = ζ． (ii)’ 緩和問題 (2) の最大解の集合は有界閉凸集合となり，そのすべての端点は F に含ま

れる．さらに，（R^{n} の Lebesgue 測度の意味で）ほとんどすべての c ∈ R^{n} に対して， (2) は C = c.hull(F) の端点で一意的な最大解 x′ を持つ．この場合， x′ は (1) の最大解となる．

ゆえに， c.hull(F) を精度よく近似する閉凸集合 C を作って，緩和問題 (2) を解くことにより，元問題 (1) の近似最大値のみならず近似最大解を得ることが期待できる．実際，上述した緩和手法は組合せ最適化や非凸最適化の分野で分枝限定法等と組み合わされて頻繁に使われている．

3 一般の数理計画問題から非凸型２次計画問題へ

次節では問題 (1) の許容領域 F の凸包 c.hull(F) の近似について議論するが，その準備として， F および C を（無限個の）２次不等式で表す． F および C の代数的な表現に２次不等式を使う理由は幾つかある．

• 線形不等式についで扱いやすい不等式．

• 表現力に富んでいる．実際，理論的には問題 (1) の有界で閉な許容領域 F は無限個の２次不等式を使えば必ず表現できる（このことは F の補集合が開集合になることより分かる）．また，任意の多項式からなる不等式は，いくつかの２次不等式に置き換えることができる．例えば， x_{1}x_{2}^{2} − 2x_{1}^{2}x_{2} − 4 ≤ 0 は

x_{1}x_{2} − x_{3} ≤ 0, −x_{1}x_{2} + x_{3} ≤ 0, x_{2}x_{3} − 2x_{1}x_{3} − 4 ≤ 0

と等価である．

• 緩和問題 (2) の許容領域 C が凸２次不等式で記述されている場合には，内点法を適用して解くことができる（[4] 等参照）．

• 元の問題 (1) の許容領域 F が２次不等式で記述されてる場合には，上述した半正定値計画緩和を適用することができる．

• 許容領域 F が２次不等式で記述される問題 (1) は応用上重要なさまざまな数理計画問題を含む．非凸型２次計画問題の他， 0-1 （混合）整数線形（および２次）計画問題，線形相補性問題，双線形行列不等式等．

S^{n} を n × n 対称行列の集合， S_{+}^{n} ⊆ S^{n} を n × n 半正定値対称行列の集合とする． γ ∈ R, q ∈ R^{n}, Q ∈ S^{n} に対して， γ を定数項， 2q^{T} x を線形項， x^{T} Qx を２次項とする２次式を qf(·; γ, q, Q) で表す．すなわち，

qf(x; γ, q, Q) ≡ γ + 2q^{T} x + x^{T} Qx (∀x ∈ R^{n})

このとき，２次関数の集合 Q および凸２次関数の集合 Q_{+} は

Q ≡ {qf(·; γ, q, Q) : γ ∈ R, q ∈ R^{n}, Q ∈ S^{n}},
Q_{+} ≡ {qf(·; γ, q, Q) : γ ∈ R, q ∈ R^{n}, Q ∈ S_{+}^{n}}

と書ける．簡単のために， qf(·; γ, q, Q) ∈ Q の代わりに， p(·) ∈ Q と書くこともある．以下では，

条件 3.1. 既知な有界閉凸集合 C_{0} と P_{F} ⊆ Q が存在して，問題 (1) の許容領域は F ≡ {x ∈ C_{0} : p(x) ≤ 0 (∀p(·) ∈ P_{F})} と表現されている．

を仮定する． C_{0} は c.hull(F) の初期近似の役割を果たす．実際に計算する上では P_{F} が有限個の２次式から構成されている必要があるが，理論展開上（定理 5.1 および定理 5.2）， P_{F} は無限個の２次式からなっていても差し支えない．また，緩和問題 (2) の凸許容領域 C も，２次不等式を用いて構成する．

次節で使われる記号を導入しておく．

\overline{D} ≡ 長さ 1 の方向ベクトルの集合≡ \{d ∈ R^{n} : ∥d∥ = 1\},
Q • X ≡ Q, X ∈ S^{n} の内積≡ ∑_{i=1}^{n} ∑_{j=1}^{n} Q_{ij}X_{ij},
α(C, d) ≡ max\{d^{T} x : x ∈ C\} (∀d ∈ \overline{D}),
ℓsf(x; C, d) ≡ d^{T} x − α(C, d) (∀d ∈ \overline{D}, x ∈ R^{n}),
P^{L}(C, D) ≡ \{ℓsf(·; C, d) : d ∈ D\} (∀D ⊆ \overline{D}),
r2sf(x; C, d_{1}, d_{2}) ≡ −(d_{1}^{T} x − α(C, d_{1}))(d_{2}^{T} x − α(C, d_{2})) (∀d_{1}, d_{2} ∈ \overline{D}),
P^{2}(C; D_{1}, D_{2}) ≡ \{r2sf(·; C, d_{1}, d_{2}) : d_{1} ∈ D_{1}, d_{2} ∈ D_{2}\} (∀D_{1}, D_{2} ⊆ \overline{D})

と定義する．ただし， C は R^{n} の有界閉凸集合． ℓsf(·; C, d) と r2sf(·; C, d_{1}, d_{2}) は，それぞれ， C に対する “linear supporting function” と “rank-2 quadratic supporting function” になっている．

4 逐次半正定値計画緩和法

アルゴリズム 4.1. ステップ０： k = 0 とする． n 本の線形独立なベクトル v_{1}, v_{2}, …, v_{n} を選び， D_{1} ≡

\{v_{1}, v_{2}, …, v_{n}, −v_{1}, −v_{2}, …, −v_{n}\} ⊆ \overline{D} と置く．（例えば，各 v_{i} を第 i 座標単位ベ

クトル e_{i} にとればよい）．

ステップ１：空でない D_{2} ⊆ \overline{D} を選ぶ．

ステップ２： α(C_{k}, d) ≡ sup{d^{T} x : x ∈ C_{k}} (d ∈ D_{1} ∪ D_{2}) を計算し，

P_{k} = P^{L}(C_{k}, D_{1}) ⋃ P^{2}(C_{k}, D_{1}, D_{2})

と置く．ステップ３：

C_{k+1} = F̂(C_{0}, P_{F} ∪ P_{k}) （半正定値計画緩和）
≡ \{x ∈ C_{0} : \begin{matrix}∃X ∈ S^{n}; ( \begin{matrix}1 & x^{T} \\ x & X \end{matrix} ) ∈ S_{+}^{1+n}, \\ γ + 2q^{T} x + Q • X ≤ 0 (∀qf(·; γ, q, Q) ∈ P_{F} ∪ P_{k}) \end{matrix} \} .

ステップ４： k = k + 1 として，ステップ 2 へもどる．

このアルゴリズムは

(a) F ⊆ C_{k+1} ⊆ C_{k} ⊆ C_{0} (∀k = 0, 1, 2, …)

を満たす有界閉凸集合の列 {C_{k} (k = 0, 1, 2, …)} を生成する． {C_{k} (k = 0, 1, 2, …)} が c.hull(F) に収束する，あるいは， c.hull(F) を近似するためにはステップ１で特別な D_{2} を選択する必要がある．これについては次節で述べる．

アルゴリズム 4.1 をコンピュータに実装する際には，条件 1.1 と条件 3.1 に加えて，以下の条件を仮定する必要がある．

条件 4.2.

(i) 問題 (1) の許容領域 F の表現 P_{F} は有限個の２次関数よりなる．

(ii) F を含む既知な有界閉凸集合 C_{0} は有限個の凸２次関数の集合 P̃_{0} を用いて C_{0} ≡ {x ∈ R^{n} : p(x) ≤ 0 (∀p(·) ∈ P̃_{0})} と表現されている．

(iii) ステップ１で選ぶ方向ベクトルの集合 D_{2} は有限集合．

この条件の下では P_{k} ≡ P^{L}(C_{k}, D_{1}) ⋃ P^{2}(C_{k}, D_{1}, D_{2}) は有限になる．また，ステップ２での α(C_{k}, d) (d ∈ D_{1} ∪ D_{2}) の計算は，線形行列不等式系（[2] 等参照）

X ∈ S_{+}^{n}, γ + 2q^{T} x + Q • X ≤ 0 (qf(·; γ, q, Q) ∈ P̃_{0} ∪ P_{F} ∪ P_{k}

の下で， d^{T} x を最大化する半正定値計画問題となり，内点法によって解くことができる（ [4] 等参照）．さらに，１反復で解かなければならない半正定値計画問題の個数も有限個で押さえられている．したがって，条件 1.1， 3.1， 4.2 の下ではアルゴリズム 4.1 は実行可能である．

5 凸包への収束および最大値の近似

この節を通して条件 1.1 と条件 3.1 を仮定する．ステップ１での D_{2} ⊆ \overline{D} の選び方によって異なった収束定理（[9, Corollaries 3.4, 3.6]）と近似定理（[9, Theorems 3.3, 3.5]）を得る．以下では，実現は困難であるが理論的には強力な収束定理，および，実現可能な条件の下で，任意の精度で問題 (1) の最大値を有限回の反復で近似する定理について述べる．定理 5.1. [9, Corollary 3.4] ステップ１で D_{2} = \overline{D} を採用すると，生成される \{C_{k} (k = 0, 1, 2, …)} は上述の (a) および以下の (b)， (c) を満たす．

(b) F = ∅ ⟹ ∃k; F_{k} = ∅.

(c) F ≠ ∅ ⟹ c.hull(F) = ⋂_{k=0}^{∞} C_{k}.

この定理はアルゴリズム 4.1 によって生成される {C_{k} (k = 0, 1, 2, …)} の c.hull(F) への収束を保証している．しかしながら，\overline{D} は連続無限個の方向ベクトルからなるので，定理の仮定 D_{2} = \overline{D} を満たすようにアルゴリズムを実行することはできない．実際， D_{2} = \overline{D} にとった場合には，ステップ２で連続無限個の不等式条件のついた連続無限個の半正定値計画問題を解かなければならず，まさに，実現不能である．アルゴリズム 4.1 を計算機に実装するためには，ステップ２で解くべき問題の個数を有限個にし，かつ，１つの問題の不等式条件を有限個に押さえる必要がある．このためには，条件 4.2 を仮定し，かつ，方向ベクトルの集合\overline{D} を離散化（discretization）する必要がある．

δ > 0, D ⊆ \overline{D} とする． D′ ⊆ \overline{D} は条件

∀d ∈ D, ∃d′ ∈ D′; ∥d′ − d∥ ≤ δ

を満たすとき， D の δ-ネットと呼ぶ． δ > 0 であれば，有限本の方向ベクトルからなる D の δ-ネット D′ が存在する．したがって，そのような D′ は D の一様な離散化を実現していると言える．

上述の定理 5.1 では目的関数 c^{T} x に関する情報が全く使われていない．通常の場合は，許容領域 F の凸法 c.hull(F) そのものが必要なのではなく，問題 (1) の最大値および最大解が得られればよい． c.hull(F) を計算するのはその手段であるが，後者のためには最大解からは遠く離れた箇所では， C_{k} が c.hull(F) を高い精度で近似する必要はない．このことを考慮に入れると， c.hull(F) の近似を目的関数の係数ベクトル c の方向に局所化（localization）することにより，一反復の計算量を減らすことに結びつく．以下の定理では，離散化と局所化の両方を定理 5.1 に組み込んでいる．

定理 5.2. [9, Theorem 3.5] c/∥c∥ ∈ \overline{D} を含む任意の開集合 N と任意の ϵ > 0 を与える．このとき，ある δ > 0 が存在して，ステップ１で D_{2} を方向ベクトルの集合\overline{D} ∩ N の δ-ネットにとれば，生成される {C_{k} (k = 0, 1, 2, …)} は上述の (a), (b) および以下の (c)’ を満たす．

(c)’ F ≠ ∅ ⟹ ∃k; ζ^{∗} ≡ max{c^{T} x : x ∈ F } ≤ max{c^{T} x : x ∈ C_{k}} ≤ ζ^{∗} + ϵ.

6 逐次半無限線形計画緩和法

アルゴリズム 4.1 のステップ３の C_{k+1} の定義において，半正定値条件

( \begin{matrix}1 & x^{T} \\ x & X \end{matrix} ) ∈ S_{+}^{1+n}

を取り除くと，逐次半無限線形計画緩和法を得る．アルゴリズム 6.1.

ステップ０，１，２：アルゴリズム 4.1 のステップ０，１，２と同じ．

ステップ３：

C_{k+1} = F̂^{L}(C_{0}, P_{F} ∪ P_{k}) （半無限線形計画緩和）
≡ \{x ∈ C_{0} : \begin{matrix}∃X ∈ S^{n}; \\ γ + 2q^{T} x + Q • X ≤ 0 (∀qf(·; γ, q, Q) ∈ P_{F} ∪ P_{k}) \end{matrix} \} .

ステップ４： k = k + 1 として，ステップ 2 へもどる．

定義より，半正定値計画緩和と半無限線形計画緩和の間には

F̂(C_{0}, P_{F} ∪ P_{k}) ⊆ F̂^{L}(C_{0}, P_{F} ∪ P_{k})

が成り立つ．すなわち，半正定値計画緩和と比較すると，半無限線形計画緩和は同等か，あるいは，より悪い緩和を与える．しかしながら，定理 5.1 および定理 5.2 はそのまま成立する．特に，条件 4.2 を仮定し，アルゴリズム 6.1 を定理 5.2 の仮定の下で実行する場合には，ステップ２では有限個の線形計画問題

\begin{matrix}目的： & d^{T} x → 最大化 \\ 条件： & γ + 2q^{T} x + Q • X ≤ 0 (∀qf(·; γ, q, Q) ∈ P̃_{0} ∪ P_{F} ∪ P_{k}) \end{matrix} \}

(d ∈ D_{1} ∪ D_{2}) を解けばよいことになる．アルゴリズム 4.1 が各反復で半正定値計画問題を解かなければならないことと比較すると，一反復あたりの計算量はアルゴリズム 6.1 のほうがずっと少ない．さらに，これらの線形計画問題では，制約条件は共通で目的関数の係数ベクトル d のみが異なっているので，パラメトリック単体法を使えばかなり効率よく計算できることが期待できる．

7 おわりに

本稿で紹介した逐次半正定値計画緩和法および逐次半無限線形計画緩和法は，その基礎理論が論文 [8, 9] で提案されたばかりで，それらが今後どのように発展し，実用化へと結びついて行くか，あるいは実用性がないことが判明するかは未知である．いずれにしろ，非凸最適化分野への新しい試みとして研究する価値は十分にあると信じている．

謝辞：大学院学生の武田朗子さんと福田光浩君は原稿を読んで誤字脱字と不明確な箇所の指摘してくれました．ありがとうございました．

参考文献

[1] E. Balas, S. Ceria and G. Cornuéjols, “A lift-and-project cutting plane algorithm for mixed 0-1 programs,” newblock Mathematical Programming 58 (1993) 295–323.

[2] S. Boyd, L. E. Ghaoui, E. Feron and V. Balakrishnan, Linear Matrix Inequalities in System and Control Theory, (SIAM, Philadelphia, 1994).

[3] T. Fujie and M. Kojima, “Semidefinite relaxation for nonconvex programs,” Journal of Global Optimization 10 (1997) 367–380.

[4] K. Fujisawa, M. Kojima and K. Nakata, “SDPA (Semidefinite Programming algorithm) User’s Manual — Version 4.10,” Dept. of Mathematical and Computing Sciences, Tokyo Institute of Technology, May 1998.

[5] M. X. Goemans, “Semidefinite programming in combinatorial optimization,” Mathematical Programming 79 (1997) 143–161.

[6] M. X. Goemans and D. P. Williamson, “Improved approximation algorithms for maximum cut and satisfiability problems using semidefinite programming,” Journal of Assoc. Comput. Mach. 42 (1995) 1115–1145.

[7] N. Karmarkar, “A new polynomial-time algorithm for linear programming,” Combinatorica 4 (1984) 373–395.

[8] M. Kojima and L. Tunçel, “Cones of matrices and successive convex relaxations of nonconvex sets,” Research Report B-338, Dept. of Mathematical and Computing Sciences, Tokyo Institute of Technology, Meguro, Tokyo, Japan, March 1998.

[9] M. Kojima and L. Tunçel, “Discretization and localization in successive convex relaxation methods for nonconvex quadratic optimization problems,” Research Report B-341, Dept. of Mathematical and Computing Sciences, Tokyo Institute of Technology, Meguro, Tokyo, Japan, July 1998.

[10] L. Lovász and A. Schrijver, “Cones of matrices and set functions and 0-1 optimization,” SIAM J. on Optimization 1 (1991) 166–190.

[11] Yu. E. Nesterov, “Semidefinite relaxation and nonconvex quadratic optimization,” CORE Discussion Paper, #9744, June 1997.

[12] Yu. E. Nesterov and A. S. Nemirovskii, Interior-Point Polynomial Algorithms in Convex Programming, (SIAM, Philadelphia, 1994).

[13] H. D. Sherali and C. H. Tuncbilek, A reformulation-convexification approach for solving nonconvex quadratic programming problems, Journal of Global Optimization 7 (1995) 1–31.

[14] N. Z. Shor, “Dual quadratic estimates in polynomial and boolean programming,” Annals of Operations Research 25 (1990) 163-168.

[15] S. J. Wright, Primal-Dual Interior-Point Methods, (SIAM, Philadelphia, 1997).