線形計画法
はじめに
線形計画問題とは,目的関数と制約条件がすべて線形である最適化問題を言う.
Mathematica には実数変数を持つ線形最適化問題を解くためのアルゴリズム群があり,LinearProgramming,FindMinimum,FindMaximum,NMinimize,NMaximize,Minimize,Maximizeでアクセスできる.LinearProgrammingを使うと,線形計画法アルゴリズムに直接アクセスできる.これは使用するメソッドを最も柔軟に指定することができ,大規模な問題に最も効率的である.FindMinimum,FindMaximum,NMinimize,NMaximize,Minimize,Maximizeは方程式および不等式の形式での線形計画問題を解くのに便利である.
| Out[1]= |  |
| Out[2]= |  |
| Out[3]= |  |
LinearProgramming関数
LinearProgramming
は線形計画法の主要関数であり,使用するメソッドの指定が最も柔軟にでき,大規模な問題に最も効率的である.
以下のオプションを取ることができる.
LinearProgrammingのオプション
Methodオプションは,線形計画問題を解くアルゴリズムを指定する.これで使用できる値はAutomatic,
,
,
である.デフォルトはAutomaticであり,問題の大きさと精度に基づいて,他のメソッドの中から自動的に選ばれる.
Toleranceオプションは,収束の許容誤差を指定する.
例題
内点法とシンプレックス法および/または改訂シンプレックス法との違い
シンプレックス法および改定シンプレックス法は,極小値に達するまで目的関数の継続的に小さい方の値で頂点から頂点へと,制約条件により定義されたポリトープの辺に沿って移動することにより線形計画問題を得.これらの関数の Mathematica の実装では,密な線形代数が使われる.この実装のユニークな点は,厳密・拡張精度問題が解けるという点である.したがって,これらの方法は,非機械数の結果が必要ない,または頂点の解が望まれるような小さいサイズの問題に適している.
線形計画法の内点法は,大まかに言うと,制約条件により定義されたポリトープの内部から反復する.これは非常に速く買いに近づくが,シンプレックス・改訂シンプレックス法とは異なり,厳密な解は見付けない.Mathematica における内点法の実装では,機械精度の疎な線形代数を使う.よって,大規模の機械精度線形計画問題では,内点法の方が効率がよいのでこれを使うべきである.
以下では,重解

と

の間の線分上にある任意の点が解となる)を持つ線形計画問題を解く.内点法では,解集合の中央にある解が返される.
| Out[6]= |  |

または

を使うと,解集合の境界の解が返される.
| Out[7]= |  |
以下は,200個の変数のランダムで疎な線形計画問題では,内点法の方がずっと速く,同様な最適値を与えることを示している.
| Out[44]= |  |
| Out[45]= |  |
| Out[46]= |  |
二重変数の探索
一般的な線形計画法の問題
がある場合,その双対問題は次の通りである.
両方の解が分かっていると便利な場合がある.
主問題と双対問題の解の関係は次の表の通りである.
| |
| 実行可能 | 実行可能 |
| 非有界 | 実行不可能または非有界 |
| 実行不可能 | 非有界または実行不可能 |
両方の問題が実行可能である場合,(P)と(D)の最適値は同じであり,主問題の解
および双対問題の解
に対して次の相補条件が成り立つ.
は,リスト
を返す.
次で主問題

の他に,下の双対問題も解く.

| Out[14]= |  |
以下により,主問題と双対問題が同じ最適値を返すことが証明される.
| Out[15]= |  |
| Out[16]= |  |
この制約条件の双対は

である.これは,制約条件の右辺が1単位増加すると,最適値は2単位増加することを意味している.このことは,制約条件の右辺を0.001で摂動させることで確認できる.
| Out[17]= |  |
| Out[18]= |  |
内点法における実行不可能性と非有界性の扱い
主双対内点法は実行可能な問題に対して設計されており,実行不可能/非有界の問題では発散してしまう. その上,現在の実装では発散が実行不可能性によるものか非有界性によるものかを見分けるのは難しい.
ヒューリスティックが実行不可能/非有界の問題を見付け,適切なメッセージを出力する.
| Out[19]= |  |
ヒューリスティックを使っても,問題が実行不可能/非有界のどちらであるかを確実に見分けることができない場合がある.
| Out[20]= |  |
メッセージで提案されているように,

を使うと問題が非有界であることが分かる.
| Out[21]= |  |
"InteriorPoint"のメソッドオプション
は
のメソッドオプションであり, 密な列を別々に扱うかどうかを決定する.密な列とは,多くの非零要素を含む制約行列の列のことである.このメソッドのデフォルト値はAutomaticであり,密な列は別々に扱われる.
密な行列を含む大きな問題では特にこのメソッドが便利である.
| Out[98]= |  |
| Out[99]= |  |
大きいデータをインポートし,大規模問題を解く
線形計画問題を記述するのによく使われる形式はMPS(数理計画問題)形式である.これは多数のセクションで構成されるテキスト形式である.
方程式のMPSファイルをインポートする
Mathematica はMPS形式のファイルをインポートすることができる.MPSデータのImportにより,線形計画問題はデフォルトで方程式形式で返される.この形式はMinimizeあるいはNMinimizeを使って解くことができる.
次で
MPSファイル

で指定された線形計画問題を解く.
| Out[25]= |  |
| Out[26]= |  |
大規模な問題:行列形式およびベクトル形式でのインポート
大規模な問題の場合,MPSデータファイルをインポートし,行列かベクトルで行列問題を表してからLinearProgrammingを使って解いた方が効率的である.
次で,
MPS形式のデータで3つの要素がインポートできることを示す.
| Out[101]= |  |
| Out[104]= |  |
疎な制約行列を求めるときだけ

指定を使うことができる.
| Out[105]= |  |
自由書式のMPSファイル
標準的なMPS形式のファイルでは固定書式が使われており,各欄の文字の位置が厳密に決まっている.しかし,モデリングシステムの中には自由書式のMPSファイルを出力するものもあり,欄の位置が自由である.後者のファイルでは,Importでオプション"FreeFormat"->Trueを指定することができる.
次の文字列は,
MPSファイルを自由形式で記述している.
一時的ファイル名を割り当て,それに文字列をエキスポートする.
オプション
"FreeFormat"->Trueを使ってファイルをインポートする.
| Out[126]= |  |
線形計画テスト問題
NetLib線形計画テスト問題はすべてExampleData関数からアクセスすることができる.
| Out[12]= |  |
問題

をインポートし,それを解く.
| Out[8]= |  |
| Out[9]= |  |

問題に対してインポートすることのできる他の特性を表示する.
| Out[10]= |  |

を方程式形式でインポートする.
| Out[11]= |  |
線形計画法の応用例題
L1ノルム最小化
最小化問題
を,線形計画問題に変換して解くことができる.

最小化問題を解くための関数を定義する.
| Out[38]= |  |

最小化の解が得られる.
| Out[39]= |  |
| Out[40]= |  |
| Out[41]= |  |
| Out[42]= |  |
最適なアンカーの設計
この例題は[2]から採用したものである.この問題の目的は,負荷をサポートするのにできるだけ少ない材料を使ってアンカーを設計するというものである.
この問題は,まずノードとリンクを使って離散化し,シミュレートすることによりモデル化することができる.モデリングの過程は以下の図で示す.ここで7×10のノードを持つ格子が生成される.それぞれのノードはマンハッタン距離が3以下のすべてのノードにリンクにより接続される. 赤い3つのノードは壁に固定されるもので,他のノードはすべて圧縮力および張力のバランスが取れていなければならない.
各リンクは厚みのある硬い棒を表しており,重さはそれにかかる力に比例している.目標は,使用される材料を最小化することであり,下を意味する.
従って,これは数学的に線形の制約付き最小化問題であり,目的関数が線形関数の絶対値の和となっている.
目的関数の絶対値
は,
を圧縮力と張力(ともに非負)に分けて置き換えることができる.ここで
をリンクの集合,
をノードの集合,
をノード
と
の間のリンクの長さ,
と
をそれぞれリンク上の圧縮力と張力とすると,上記のモデルは下の線形計画問題に変換することができる.
次でモデルを設定し,それを解き,結果をプロットする.これはAMPLモデルに基づく(
[2]を参照).
30個のノードを水平・垂直方向に並べて問題を解く.
| Out[6]= |  |
しかし,このアンカーを壁に固定せず,空間のいくつかの点に固定する場合,結果が葉の形になる.これは多分,進化の過程で葉の構造が最適化されたのであろう.
| Out[11]= |  |
線形計画のアルゴリズム
シンプレックス法と改訂シンプレックス法
シンプレックス法と改訂シンプレックス法は,制約条件により定義された多面体の頂点で実行可能解を構築し,その後目的関数のより小さな値を連続的に使い,最小値に達するまでその多面体の辺を頂点に向かって移動させていくことにより線形計画問題を解く.
シンプレックス法と改訂シンプレックス法の実装が少ないと実際に非常に効率的であり,大域的最適値が見付かることが保障されているが,最悪の場合には動作に問題が生じる場合もある.構築された線形計画問題に対して,シンプレックス法と改訂シンプレックス法がその問題の大きさに指数関数的な数のステップを取らなければならなくなる可能性があるのである.
Mathematica は密な線形代数を使ってシンプレックス法と改訂シンプレックス法を実装する.この実装のユニークな点は,厳密・拡張精度問題が解けるという点である.従って,これらのメソッドは機械数でない結果が必要とされる小さい問題により適している.
以下で20の制約条件と200個の変数を持つランダムな線形計画問題を設定する.
この問題を解く.制約条件の数よりも変数の数の方がかなり多い線形計画問題では,通常改訂シンプレックス法の方が早い.これに反して,変数よりも制約条件の数のほうが多い場合は,シンプレックス法の方が速い.
機械数の結果だけが必要な場合は,問題を機械数の問題に変換し,内点法を使わなければならない.
内点法
シンプレックス法と改訂シンプレックス法は通常非常に効率的であるが,最悪の動作をする場合もある.シンプレックス法と改訂シンプレックス法が問題の大きさに指数関数的な数のステップを取らなければならないような線形計画問題を構築する可能性があるのである.
さらに Mathematica ではシンプレックス法と改訂シンプレックス法の実装には密な線形代数が使われるが,内点法の実装には機械数の疎な線形代数が使われる.従って,機械数の大規模な線形計画問題には,内点法の方が効率的なのでこれを使うべきである.
内点法の公式
標準化された線形計画問題
を考えてみる.ここで
,
,
である.この問題は,正の制約条件を扱う障壁関数式を使って解くことができる.
上の問題の一階必要条件により,以下が与えられる.
がベクトル
からなる対角行列を表し,
であるとすると,
が成り立つ.これは制約条件を持つ
の線形/非線形方程式の集合である.これはニュートン法
であり
で解くことができる.この線形系を解く方法のひとつに,ガウス消去法を使い行列をブロック三角行列形式に簡約するというものがある.
このブロック三角行列を解くためには,この行列が通常の系に含まれる,いわゆる通常の系を解く.
この行列は正定値行列になるが,解に近付くにつれ非常に悪条件になる.よって解法を安定させるために数値的手法が使われる.通常内点法は,「収束の許容誤差」で説明する許容範囲では,ほぼ
付近までの問題を解くよう想定されている.Mathematica ではMehrotraの予測子・修正子法[1]が使われている.
収束の許容誤差
一般の線形計画問題はまず,標準形の
に変換される.ここで対応する双対問題は次の通りである.
内点法の収束条件は
であり,許容範囲はデフォルトで
に設定されている.
参考文献
[1] Vanderbei, R. Linear Programming: Foundations and Extensions. Springer-Verlag, 2001.
[2] Mehrotra, S. "On the Implementation of a Primal-Dual Interior Point Method." SIAM Journal on Optimization 2 (1992): 575-601.