matsulibの日記

自分用メモ

GurobiとPythonで数理最適化

数理最適化ソルバーGurobiにはPython含む複数の言語のインターフェースが用意されていて、
gurobipyモジュールを使うと使い慣れたPythonがソルバーのモデリング言語になる。

例題

http://www.fujilab.dnj.ynu.ac.jp/lecture/system2.pdf

あるレストランで,手持ちの材料からハンバーグとオムレツを作って利益を最大にしたいと考
えている.手持ちの材料は,
• ひき肉 3800 [g]
• タマネギ 2100 [g]
ケチャップ 1200 [g]
であり,それぞれの品を作るのに必要な材料の量は,
• ハンバーグ 1個あたり,ひき肉 60 [g],タマネギ 20 [g],ケチャップ 20 [g]
• オムレツ1個あたり, ひき肉 40 [g],タマネギ 30 [g],ケチャップ 10 [g]
であるとする.(他に必要な材料は十分な量があるものとする)
販売価格は,
• ハンバーグ 400 [円/個]
• オムレツ 300 [円/個]
とする.総売上を最大にするには,それぞれハンバーグとオムレツを幾つずつ作れば良いか?

定式化

x1をハンバーグの個数、x2をオムレツの個数として、次のように定式化する。


\begin{eqnarray} 
\mbox{maximize}&~& 400x_1 + 300x_2
\end{eqnarray}

\begin{eqnarray}
\mbox{subject to}
&~& 60x_1 + 40x_2 \leq 3800 \\
&~& 20x_1 + 30x_2 \leq 2100 \\
&~& 20x_1 + 10x_2 \leq 1200 \\
&~& x_i \geq 0, \\
&~& x_i \in {\bf\mbox{Z}},\:\:\: i=1,2
\end{eqnarray}

Pythonによるモデルの記述

gurobipyモジュールのメソッドを使ってモデルに変数や目的関数、制約式などを追加していく。

from gurobipy import *

try:
    # Create a new model
    m = Model("restaurant")
    
    # Create variables
    x = [m.addVar(vtype=GRB.INTEGER, name='x'+str(i+1)) for i in range(2)]
    
    # Integrate new variables
    m.update()
    
    # Set objective
    m.setObjective(400*x[0] + 300*x[1], GRB.MAXIMIZE)
    
    # Add constraints
    m.addConstr(60*x[0] + 40*x[1] <= 3800, "ground meat")
    m.addConstr(20*x[0] + 30*x[1] <= 2100, "onion")
    m.addConstr(20*x[0] + 10*x[1] <= 1200, "ketchup")

    m.optimize()

    # OUTPUT
    for v in m.getVars():
        print v.varName, v.x
    print 'Obj:', m.objVal

except GurobiError:
    print 'Error reported'

コード内に非負制約を明示していないが、それは決定変数を宣言するaddVarメソッドのデフォルト値として設定されている。引数の数やデフォルト値は対話シェルでも確認することができる。

>>> help(m.addVar)
Help on method addVar in module gurobipy:

addVar(...) method of gurobipy.Model instance
    ROUTINE:
      addVar(lb, ub, obj, vtype, name, column)

    PURPOSE:
      Add a variable to the model.

    ARGUMENTS:
      lb (float): Lower bound (default is zero)
      ub (float): Upper bound (default is infinite)
      obj (float): Objective coefficient (default is zero)
      vtype (string): Variable type (default is GRB.CONTINUOUS)
      name (string): Variable name (default is no name)
      column (Column): Initial coefficients for column (default is None)

    RETURN VALUE:
      The created Var object.

    EXAMPLE:
      v = model.addVar(ub=2.0, name="NewVar")

つまり、v=m.addVar() は v=m.addVar(0,inf,0,GRB.CONTINUOUS,'',None) と同義になる。
また、キーワード引数を上手く使うと楽に書ける。

決定変数の添字が定式化では1オリジンなのにソースコードでは(リストを使っているので)0オリジンで気持ち悪いという場合は辞書を使うと良い。

>>> x = {}
>>> for i in range(1, 3):
...     x[i] = m.addVar(vtype =GRB.INTEGER, name='x'+str(i))

>>> m.update()
>>> x
{1: <gurobi.Var x1>, 2: <gurobi.Var x2>}

もちろん、辞書のキーには整数だけでなく文字列やタプルなどを使うこともできる。

実行結果
Optimize a model with 3 rows, 2 columns and 6 nonzeros
Found heuristic solution: objective 24000
Presolve time: 0.00s
Presolved: 3 rows, 2 columns, 6 nonzeros
Variable types: 0 continuous, 2 integer (0 binary)

Root relaxation: objective 2.700000e+04, 3 iterations, 0.00 seconds

    Nodes    |    Current Node    |     Objective Bounds      |     Work
 Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

/*    0     0               0    27000.000000 27000.0000  0.00%     -    0s

Explored 0 nodes (3 simplex iterations) in 0.00 seconds
Thread count was 2 (of 2 available processors)

Optimal solution found (tolerance 1.00e-04)
Best objective 2.700000000000e+04, best bound 2.700000000000e+04, gap 0.0%
x1 30.0
x2 50.0
Obj: 27000.0

よって、ハンバーグ30個、オムレツ50個のとき、総売上27000円で最大となる。

まとめ

なによりソルバーをPythonソースコードの中で使えることが一番うれしい。

  • 一般的なモデリング言語(AMPLなど)を覚えなくて済む
  • ソルバーの実行結果を他の処理に使いたいとき、結果を取り出してソースコードに手打ちしたり別途バッチ処理を書いたりする手間が省ける
  • 決定変数やその添字に名前を付けることで、単に可読性のためだけでなく、文字列処理でマッチするオブジェクトを抜き出して計算するようなことも割と楽に書ける

※ 有償ソフトをアカデミックライセンスで使用中。フリートライアル版ではすごく小規模な問題しか解けない。