Python + Gurobi による科学技術計算環境の構築

Gurobi のアカデミック・ライセンスは以下の制約の下で利用可能である：

• 学位授与機構を備えた学術機関の教員・学生もしくはスタッフのみが利用できる
• 研究もしくは教育目的にのみ利用できる
• 公認の学術機関ドメインから validate される必要がある．アカデミック・ライセンス の有効期限は12ヶ月だが，適切な教員・学生もしくはスタッフは，この手続きを繰り返 すことでライセンスを更新できる．

1. http://www.gurobi.com/ にアクセスし，Register for Free をクリック

   

2. Account Type は「Academic」とし，Email は tohoku.ac.jp で終わるアドレスを利 用．Company/University には「Tohoku University」とする

3. http://www.gurobi.com/ のメニューバーから Download > Licenses に移動し，FREE ACADEMIC のページに移動．2箇所をチェックして Request License する

   

4. ライセンスが発行される．Gurobi をインストールした後，ここで表示される grbgetkey 以下のコマンドをターミナルで入力するとライセンスファイルが生成さ れる．

   

1. http://www.gurobi.com/ のメニューバーから Download > Gurobi Optimizer のページ に移動し，Current version: 5.6.3 の Mac OS 版と READ ME を選択してDownload する．

   

2. README.txt に従ってインストールする(ダウンロードした pkg ファイルをダブルクリッ クすればインストールできる)．実行ファイルなどは

   /Library/gurobi***/mac64
   

   以下にインストールされる．ここで, *** は Gurobi のバージョンによって変わる． バージョンが 5.6.3 なら，パスは /Library/gurobi563/mac64 になる．

3. ターミナルを立ち上げ，上記のライセンス取得コマンドを実行する：

   grbgetkey ae36ac20-16e6-acd2-f242-4da6e765fa0a # このキーは実際には使用できない
   

   Gurobi のサーバーに接続して認証が完了すると，ライセンスファイルの生成場所を尋 ねてくる：

      Gurobi license key client (version x.x.x) # x.x.x はバージョンによって変わる
   Copyright (c) 2012, Gurobi Optimization, Inc.
   
   -------------------------------
   Contacting Gurobi key server...
   -------------------------------
   
   Key for license ID 45789 was successfully retrieved.
   License expires at the end of the day on 2014-08-26.
   
   ---------------------
   Saving license key...
   ---------------------
   
   In which directory would you like to store the Gurobi license key file?
   [hit Enter to store it in /Users/[username]]]:
   

   デフォルトではホームフォルダにライセンスファイルが生成される．

4. ターミナルで gurobi.sh を実行してエラーが出なければOK:

   $ gurobi.sh                     # gurobi.sh を実行
   Python 2.7.5 (default, Aug 25 2013, 00:04:04) 
   [GCC 4.2.1 Compatible Apple LLVM 5.0 (clang-500.0.68)] on darwin
   Type "help", "copyright", "credits" or "license" for more information.
   
   Gurobi Interactive Shell (mac64), Version 5.6.3
   Copyright (c) 2013, Gurobi Optimization, Inc.
   Type "help()" for help
   

   ここでは Ctrl-D を押して終了させよう．

   ライセンスファイルが無効な場合，以下のようなメッセージが表示される：

   Python 2.7.2 (default, Oct 11 2012, 20:14:37) 
   [GCC 4.2.1 Compatible Apple Clang 4.0 (tags/Apple/clang-418.0.60)] on darwin
   Type "help", "copyright", "credits" or "license" for more information.
   Traceback (most recent call last):
     File "/usr/local/lib/gurobi.py", line 5, in <module>
       from gurobipy import *
     File "/usr/local/lib/python2.7/site-packages/gurobipy/__init__.py", line 1, in <module>
       from gurobipy import *
     File "gurobi.pxi", line 11, in init gurobipy (../../src/python/gurobipy.c:72659)
     File "env.pxi", line 27, in gurobipy.Env.__init__ (../../src/python/gurobipy.c:3099)
   gurobipy.GurobiError: No Gurobi license found (user nagae, host iMac-Late2012.local, hostid b1b0acc5)
   

   この場合も終了させるには Ctrl-D．

ノード集合 とリンク集合 を以下のように定義する：

このとき，最小費用流問題は，以下の3つの要素で構成される：

1. 最小化すべき目的関数(総輸送費用)：

   

2. 制約条件その1(各ノードでのフロー保存則)：

   

   ただし， は， が起点(a)なら , が終点なら, それ以外なら0 となる定数．

3. 制約条件その2(非負制約)：

   

• まずは Python モードを起動:

  $ python
  

• リンクとノード集合を定義する：

  >>> links = [('a','b'), ('a','c'), ('b', 'c'), ('b','d'),('c','d')]
  >>> nodes = ['a', 'b', 'c', 'd']
  >>> orig = 'a'                      # 起点
  >>> dest = 'd'                      # 終点
  

• ノード集合は各リンクの起点もしくは終点の集まりなので，以下のように生成してもよ い．

  >>> nodes = set([i for (i, j) in links] + [j for (i, j) in links])
  

• 各リンクの費用と容量を設定する

  >>> cost = dict(zip(links,[2,8,3,8,4]))
  >>> capacity = dict(zip(links, [10, 10, 3, 9, 8]))
  

• 総輸送量を設定する

  >>> total_volume = 17
  

• Python 用 Gurobi パッケージを「grb」という名前で読み込む

  >>> import gurobipy as grb
  

• モデルを格納する grb.Model型変数 m を作成する

  >>> m = grb.Model()
  

• GRB.addVar()という関数を用いてリンクごとの未知変数を定義する．

  >>> x = {}
  >>> for (i,j) in links:
  ...    x[i, j] = m.addVar(vtype=grb.GRB.CONTINUOUS, name="x_{%s,%s}" % (i, j)) # x[i, j]の前にはインデントが必要
  ...                             # 何も入力しないで改行
  

• 目的関数を定義し，最大化/最小化のいずれを行うかを指定する

  >>> m.update()                  # モデルに変数が追加されたことを反映させる
  >>> m.setObjective(grb.quicksum(x[i, j]*cost[i, j] for (i, j) in links)) # 目的関数
  >>> m.setAttr("ModelSense", grb.GRB.MINIMIZE)
  

• ノードごとのフロー保存則を追加する

  >>> b = {}
  >>> for n in nodes:
  ...    # フロー保存則の右辺の定数を決定する
  ...    if n == orig:            # if の前にはインデントが必要
  ...        b[n] = - total_volume # b[n] の前には更にインデントが必要
  ...    elif n == dest:
  ...        b[n] = total_volume
  ...    else:
  ...        b[n] = 0
  ...    # ノードごとのフロー保存則
  ...    # 行末に「\」をつけると複数行にまたがって記述できる
  ...    m.addConstr(grb.quicksum(x[i, j] for (i, j) in links if j == n)\ 
  ...                - grb.quicksum(x[i, j] for (i, j) in links if i == n) \
  ...                == b[n], name="flow reservation at %s" % n)
  ...    # 空白行を入力することで for ブロックを抜けられる
  

• 非負制約と容量制約を追加する

  >>> for (i, j) in links:
  ...    x[i, j].lb = 0.0
  ...    x[i, j].ub = capacity[i, j]
  ...    # 空白行を入力することで for ブロックを抜けられる
  

• モデルを更新し，mincostflow.lp に線形計画問題の形で出力する

  >>> m.update()
  >>> m.write("mincostflow.lp")
  

• 最適化を行う

  >>> m.optimize()
  Optimize a model with 4 rows, 5 columns and 10 nonzeros
  Presolve removed 2 rows and 3 columns
  Presolve time: 0.00s
  Presolved: 2 rows, 2 columns, 4 nonzeros
  
  Iteration    Objective       Primal Inf.    Dual Inf.      Time
         0    1.8000000e+02   2.000000e+00   0.000000e+00      0s
         1    1.8300000e+02   0.000000e+00   0.000000e+00      0s
  
  Solved in 1 iterations and 0.00 seconds
  Optimal objective  1.830000000e+02
  

• 最適値と最適解(最小費用流)を表示させる

  >>> print "optimal value:\t%8.4f" % m.ObjVal
  optimal value:  183.0000
  >>> for (i, j) in links:
  ...    print "%s:\t%8.4f" % (x[i, j].VarName, x[i, j].X)
  ...
  x_{a,b}:         10.0000
  x_{a,c}:          7.0000
  x_{b,c}:          1.0000
  x_{b,d}:          9.0000
  x_{c,d}:          8.0000