日々の流転

λ. 強双対性を使って線形計画を制約充足問題に還元する

Linear and Discrete Optmization の双対性に関する部分(参考: 第五週目メモ)を復習していて気になったのだけど、強双対性の条件 c^T x^* = b^T y^* って線形式で書けてしまっているのよね。 ということは主問題と双対問題の制約条件と強双対性の連言をとっても線形(QF_LRA)の範囲に収まっていて、これの解(実行可能解)を求めれば、目的関数に関する最適化をしなくても最適解が求まることになるはず。 ということで試してみた。

主問題

簡単な問題として、Linear and Discrete Optmization で最初の例として使われていた、以下の例を考える。

MAXIMIZE 
  100 x1 + 125 x2
SUBJECT TO
  3 x1 + 6 x2 <= 30
  8 x1 + 4 x2 <= 44
BOUNDS
  0 <= x1 <= 5
  0 <= x2 <= 4
END

最適解は (x1,x2)=(4,3) で、目的関数値は775。

双対問題

この問題の双対問題は以下のようになる。

MINIMIZE
  30 y1 + 44 y2 + 5 y3 + 4 y4
SUBJECT TO
  3 y1 + 8 y2 + y3 >= 100
  6 y1 + 4 y2 + y4 >= 125
BOUNDS
  0 <= y1
  0 <= y2
  0 <= y3
  0 <= y4
END

最適解は (y1,y2,y3,y4) = (50/3, 25/4, 0, 0) で、目的関数値は775。

両方を合わせたもの

MAXIMIZE 
  0 x1
SUBJECT TO
  3 x1 + 6 x2 <= 30
  8 x1 + 4 x2 <= 44
  3 y1 + 8 y2 + y3 >= 100
  6 y1 + 4 y2 + y4 >= 125
  100 x1 + 125 x2 - 30 y1 - 44 y2 - 5 y3 - 4 y4 = 0
BOUNDS
  0 <= x1 <= 5
  0 <= x2 <= 4
  0 <= y1
  0 <= y2
  0 <= y3
  0 <= y4
END

(x1,x2,y1,y2,y3,y4) = (4, 3, 50/3, 25/4, 0, 0) はこの問題の(実行可能)解で、このときの目的関数値はどちらも775。

Z3の場合

最適化機能を持たない充足可能性ソルバの例として、SMTソルバのZ3で試してみる（オンラインでも Z3 online at RiSE4Fun で試せる）。

(declare-fun x1 () Real)
(declare-fun x2 () Real)
(define-fun obj-primal () Real
  (+ (* 100.0 x1) (* 125.0 x2)))
(assert (<= (+ (* 3.0 x1) (* 6.0 x2)) 30.0))
(assert (<= (+ (* 8.0 x1) (* 4.0 x2)) 44.0))
(assert (and (<= 0.0 x1) (<= x1 5.0)))
(assert (and (<= 0.0 x2) (<= x2 4.0)))

(declare-fun y1 () Real)
(declare-fun y2 () Real)
(declare-fun y3 () Real)
(declare-fun y4 () Real)
(define-fun obj-dual () Real
  (+ (* 30.0 y1) (* 44.0 y2) (* 5.0 y3) (* 4.0 y4)))
(assert (>= (+ (* 3.0 y1) (* 8.0 y2) y3) 100.0))
(assert (>= (+ (* 6.0 y1) (* 4.0 y2) y4) 125.0))
(assert (<= 0.0 y1))
(assert (<= 0.0 y2))
(assert (<= 0.0 y3))
(assert (<= 0.0 y4))

(assert (= obj-primal obj-dual))

(check-sat) ; => sat

(get-value (obj-primal x1 x2))
; => ((obj-primal 775.0) (x1 4.0) (x2 3.0))

(get-value (obj-dual y1 y2 y3 y4))
; => ((obj-dual 775.0) (y1 (/ 50.0 3.0))
;     (y2 (/ 25.0 4.0)) (y3 0.0) (y4 0.0))

感想とか

整数上の場合だと「解が一つ見つかった後に、目的関数値がそれよりも真に良いということを表す制約を追加して、再度解き直す」ということを繰り返すことで、有界な場合には最適化を行えるのに対して、実数上だと有界であってもZeno状態みたいなことに陥って、最適解に到達できない場合があるはず。 なので、実数の場合には充足可能性ソルバでは最適化問題は解けないと思っていたけれど、こういうやり方をすれば解けるんだなぁ。 ということで面白かった。

まあ、浮動小数点数で計算するような場合には、誤差蓄積で解が見つからなくなるような場合もあるような気はするけれど。