GLPK (GNU Linear Programming Kit)
GLPK (GNU Linear Programming Kit)はオープンソースの混合整数最適化法のソルバーである．Rのパッケージ glpk や Rglpk を利用すれば，RからGLPKを使用することができる．モデルとデータを別ファイルに記述することができるので，同一モデルによる計算はデータファイルのみを変更すればよい．また，APIが提供されているので，最適化計算を含むアプリケーションを自作する時，最適化計算の部分でGLPK APIを利用するということもできる．GLPKのインストール
Windowsのバイナリ版の場合
- こちらからwinglpk-xxxx.zip (2020年9月15日時点では winglpk-4.65.zip)をダウンロードする．
- ダウンロードしたzipファイルを展開する．
- Windows が 32 bit版であれば，w32フォルダの glpsol.exe と glpk_4_65.dll を実行パスの通っているフォルダにコピーする．もちろん，DLLファイルはWindows\System32フォルダなどでもよい． Windows が 64 bit版であれば，w64フォルダのglpsol.exe と glpk_4_65.dll を使う．
ソースコードからコンパイルする場合
Linux環境(gcc)を使用することと仮定する．Unix系OSなら同じようにできると思う．- GLPKのページに入り，最新のglpk-xxx.tar.gz(2020/09/15時点では glpk-4.65.tar.gz)をダウンロードする．
- 以下の手順でインストールする．
$ tar zxf glpk-xxx.tar.gz
$ cd glpk-xxx
$ ./configure
$ make
$ make check
$ make install
その他のシステムやインストール方法
Linuxならディストリビューションごとにインストール方法があり，GLPKが対応していれば，そちらの方法でインストールするほうが便利かもしれない．また，Windows で VC++ や Borland C++ によりコンパイルすることもできる．GLPKの基本的な使い方
線形最適化法
次の線形最適化問題を例として用いる．
最大化 z = 2x1 + 3x2 制約条件 1x1 + 3x2 ≤ 24 4x1 + 4x2 ≤ 48 2x1 + 1x2 ≤ 22 x1 ≥ 0, x2 ≥ 0
この問題を解くためのプログラム（モデルファイル）は次のように書く． ダウンロードはこちら
/* LP.mod.txt */数式と似た表記なので，モデルの記述はこれだけで大体分かると思う．
var x1 >= 0 ;
var x2 >= 0 ;
maximize z: 2*x1 + 3*x2 ;
s.t. st1: x1 + 3*x2 <= 24 ;
s.t. st2: 4*x1 + 4*x2 <= 48 ;
s.t. st3: 2*x1 + x2 <= 22 ;
end;
- /* と */ で囲まれた部分はコメントで計算には影響を与えない．
- s.t. は subject to の略で，制約条件であることを表している．
- 決定変数の値の制約として上限値と下限値を両方指定したい時は次のように書く．
var x1, >=0, <=9;モデルファイルの拡張子は mod を用いるのが慣習のようであるが，拡張子を変えても計算はできる．Windowsの場合は txt とするのが使いやすいであろう．ただし，モデルファイルであることは明示しておきたいので， LP.mod.txt のようなファイル名にしておくとよいであろう．
モデルファイルができたら，ターミナル，Windowsならコマンドプロンプトで次のように打つ．
> glpsol -m LP.mod.txt -o LP.out.txtすると，LP.out.txtというファイルに結果が書き込まれる．LP.out.txtの内容は次のようになっている．
Problem: LP
Rows: 4
Columns: 2
Non-zeros: 8
Status: OPTIMAL
Objective: z = 30 (MAXimum)
No. Row name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 z B 30
2 st1 NU 24 24 0.5
3 st2 NU 48 48 0.375
4 st3 B 18 22
No. Column name St Activity Lower bound Upper bound Marginal
------ ------------ -- ------------- ------------- ------------- -------------
1 x1 B 6 0
2 x2 B 6 0
以下略
- 表の1行目は，z の Activity が30とある．これは目的関数の最大値が30であることを意味している．
- 表の2行目は，st1 の Activity が24，Upper bound が24，Marginalが0.5とある． これは，最適解において最初の制約(st1)の値は24で，制約式の右辺の値（上限値）が24で，制約の上限を満たしていることを表している． Marginalの0.5は潜在価格（シャドープライス）と呼ばれるもので，制約式の右辺の値を1大きくすると，zの最大値を0.5大きくできるという意味である．
- 表の4行目は，st3 の Activity が18で Upper bound が22とある． 最適解において，3番目の制約の値は18で，上限の22に対して余裕があるという意味である． Marginalが空欄なのは，制約の上限を満たしていないのに，制約の上限を22からさらに大きくしても目的関数の値は大きくならないからである．
- 下の表はx1, x2 の Activity が各々 6, 6 とあり，これが最適解である．
感度分析を行うには，次のように打つ．
> glpsol -m LP.mod.txt --ranges LP.sens.txtすると，LP.sens.txtというファイルに結果が書き込まれる．LP.sens.txtの内容は次のようになっている．
GLPK 4.65 - SENSITIVITY ANALYSIS REPORT Page 1
Problem: LP
Objective: z = 30 (MAXimum)
No. Row name St Activity Slack Lower bound Activity Obj coef Obj value at Limiting
Marginal Upper bound range range break point variable
------ ------------ -- ------------- ------------- ------------- ------------- ------------- ------------- ------------
1 z BS 30.00000 -30.00000 -Inf 26.00000 -1.00000 . st1
. +Inf 30.00000 +Inf +Inf
2 st1 NU 24.00000 . -Inf 16.00000 -.50000 26.00000 st3
.50000 24.00000 36.00000 +Inf 36.00000 x1
3 st2 NU 48.00000 . -Inf 32.00000 -.37500 24.00000 x1
.37500 48.00000 54.40000 +Inf 32.40000 st3
4 st3 BS 18.00000 4.00000 -Inf 8.00000 -.60000 19.20000 st2
. 22.00000 24.00000 1.00000 48.00000 st1
GLPK 4.65 - SENSITIVITY ANALYSIS REPORT Page 2
Problem: LP
Objective: z = 30 (MAXimum)
No. Column name St Activity Obj coef Lower bound Activity Obj coef Obj value at Limiting
Marginal Upper bound range range break point variable
------ ------------ -- ------------- ------------- ------------- ------------- ------------- ------------- ------------
1 x1 BS 6.00000 2.00000 . -Inf 1.00000 24.00000 st2
. +Inf 10.00000 3.00000 36.00000 st1
2 x2 BS 6.00000 3.00000 . 2.00000 2.00000 24.00000 st1
. +Inf 8.00000 6.00000 48.00000 st2
End of report
Page 1の表は制約式の右辺に関する分析である．
- Aacitivityを見ると，制約条件st 1の左辺の値は24，制約条件st 2の左辺の値は48，制約条件st 3の左辺の値は18であることが分かる． st 1とst 2は制約条件の右辺の値と一致し，st 3は制約式右辺の値である22より4小さい．
- Slack/Marginalを見ると，制約条件の左辺と右辺の値の差(slack)と潜在価格(marginal)が分かる． st 1とst 2については，左辺と右辺の値が等しいのでslackはない．st 3は左辺の値が18で右辺の値が22なので，slackは4である．
- 潜在価格は，左辺と右辺の値が等しい場合にのみ値を持つことがある． st 1の潜在価格は0.5である．これは制約条件st 1の右辺の値を24から1増加して25にすると，目的関数の最大値が0.5増加するという意味である． st 2の潜在価格は0.375である．これは制約条件st 2の右辺の値を48から1増加して49にすると，目的関数の最大値が0.375増加するという意味である． st 3は元々4のslackがあるので，st 3のみを増加しても目的関数をこれ以上増加することはできない．
- Activity rangeの列は潜在価格が有効な制約条件右辺の値の範囲である． st 1の潜在価格を0.5としてよい制約式右辺の値の範囲は16から36である． st 2の潜在価格を0.375としてよい制約式右辺の値の範囲は32から54である． st 3のActivity rangeの解釈は分からない．
- Activityの列は最適解を示す．ここではx1, x2 の activity がともに6で，最適解がx1 = 6, x2 = 6であることを意味する．
- Obj coef/Marginalを見ると，目的関数の係数(obj coef)と限界コスト(marginal)が分かる． ここでは目的関数のx1の係数が2，x2の係数が3である． 限界コストはactivityが0の変数に対して値を持つことがある． ここではx1 = 6, x2 = 6なので，限界コストは示されていない． 仮にx1 = 0で，限界コストが-0.1であったとする． このとき，最適解のx1 = 0でなくx1 = 1とすると，目的関数は最大値から0.1減少する．
- Obj coef rangeは目的関数の係数を変化しても最適解が変化しない範囲を示している． x1の係数を1から3の範囲で変化させても，最適解はx1 = x2 = 6のままである． また，x2の係数を2から6の範囲で変化させても，最適解はx1 = x2 = 6のままである．
整数最適化法
先ほどの線形最適化問題を決定変数が非負整数をとる整数最適化問題に変更する．
最大化 z = 2x1 + 3x2 制約条件 1x1 + 3x2 ≤ 24 4x1 + 4x2 ≤ 48 2x1 + 1x2 ≤ 22 x1, x2 = 0, 1, 2, ...
この整数最適化問題のモデルファイルは，線形最適化問題のモデルファイルを次のように変更すればよい． ダウンロードはこちら
/* IP.mod.txt */0-1変数の場合は，integer の代わりに binary と指定する．その場合は不等号は不要で，単に
var x1 , integer, >= 0 ;
var x2 , integer, >= 0 ;
maximize z: 2*x1 + 3*x2 ;
s.t. st1: x1 + 3*x2 <= 24 ;
s.t. st2: 4*x1 + 4*x2 <= 48 ;
s.t. st3: 2*x1 + x2 <= 22 ;
end;
var x1, binary ;と書けばよい．コマンドラインでのコマンドは同じである．
> glpsol -m LP.mod.txt -o LP.out.txtただし，整数最適化法では感度分析は使えない．