JuliaとJuMPによる数理最適化

Julia はプログラミング言語であり,JuMP はJuliaに埋め込まれた(embedded)数理最適化のためのモデリング言語である. インストールしない場合は,JuliaBox でブラウザから利用できる.ただし,googleでのサインインが必要である.

このページの内容は,JuMP – Julia for Mathematical Programming によっている.

で使用できるようになる.

目的関数のタイプ

使える目的関数のタイプは

  • 線形
  • 凸二次
  • 非線形(凸も非凸でも)

制約式のタイプ

扱える制約式のタイプは

  • 線形
  • 凸二次
  • 二次錐
  • 半正定値
  • 非線形(凸も非凸も)

変数のタイプ

  • 連続
  • 整数
  • 半連続

使い方

JuMPをつかうには,Juliaの命令として,

とする. ここではソルバとしてCbcを使うことにするので,次のようにする:

モデル作成には次のコマンドを実行する:

m = Model(solver=CbcSolver())

ここで,mはモデルの名前で,なんでもよい.

変数の定義

変数は,モデルと関連付けて定義される. upper bound, lower boundは省略してもよい.

@variable(m,0<=x<=2)
@variable(m,0<=y<=30)

リストがあるときにその要素それぞれに対して変数を定義することもできる.

s=["Green","Blue"]
@variable(m,x[-10:10,s],Int)

これにより,2つの添え字を持つ変数x[i,j]が定義される.iは-10,-9,-8,…,-1,0,1,2,…,10をとり,jはGreenかBlueをとる.たとえばx[-4,”Green”]

上限と下限は,一部の添え字の変数のみに設定することができる:

@variable(m,x[i=1:10]>=i)

目的関数の設定

目的関数の設定には次のようにする:

@objective(m,Max,sum{x[i,j],i=1:10,j=["Green","Blue"]})

制約式の設定

たとえば,10次元のベクトル変数zの要素に対して制約式を指定するには, 次のようにする:

ここで注意することは,変数に係数をかけるときは,係数*変数の順番にすることである.逆順,すなわち変数*係数にするとエラーになる.

線形最適化問題の例

using JuMP
using Cbc
m=Model(solver=CbcSolver())
@variable(m,0<=x<=2)
@variable(m,0<=y<=30)
@objective(m,Max,5x+3*y)
@constraint(m,x+5y<=3.0)
print(m)
status=solve(m)
println("Objective value:",getobjectivevalue(m))
println("x=",getvalue(x))
println("y=",getvalue(y))

変数(詳細)

変数のより詳しい使い方を紹介する.

2つの添え字を持つ変数の一部に対する処理をする場合は,次のようなやり方がある

@variable(m,x[i=1:10,j=i:10]>=0)
@variable(m,x[i=1:10,j=1:10;isodd(i+j)]>=0)

初期値を設定することができる.

@variable(m,x[i=1:10],start=(i/2))

これは,次の処理と同じ効果がある

@variable(m,x[i=1:10])
for i in 1:10
        setvalue(x[i],i/2)
end

変数のタイプを指定するには,タイプを指定する.たとえば,0-1変数を定義するには,次のようにする:

半正定値行列変数も簡単に定義することができる.ただし,問題を解くには,適切な(半正定値最適化問題に対応した)ソルバを用いる必要がある:

半正定値制約は,次で定義することができる:

アフィン表現

係数ベクトル,変数ベクトル,定数を指定することで,アフィン表現を定めることができる:

凸2次制約を表すのに,QuadExprを用いることができる:

quad=QuadExpr([x,y],[x,z],[3.0,4.0],aff)
#3x^2+4yz+3x+4z+2

式の表現

制約式で用いる式の表現を定めることができる.たとえば,次のような使い方をする:

@expression(m,shared,sum{i*x[i],i=1:5})
@constraint(m,shared+y>=5)
@constraint(m,shared+z<=10)

制約式の参照

制約を生成した後に操作するには,参照が必要である.これには,@constraint()を実行する際に適当な引数を用いるのが簡単である:

m=Model(solver=CbcSolver())
@variable(m,x[1:5])
@variable(m,y[1:5])
@objective(m,Min,sum{x[i]+y[j],i=1:5,j=1:5})
@constraint(m,constr[i=1:5,j=1:5;i+j>=3],x[i]-y[j]==1)

ベクトル操作

線形代数演算を表す演算子(たとえば掛け算,足し算)は,適宜定義されている.たとえば,行列Aとベクトルxの掛け算は,A*xで求められる:

これを用いると,乱数によって生成した制約式を追加することができる:

制約式の変更

制約式における変数の係数を変更することはできない.上限下限を変更することはできる.

目的関数の変更

@objectiveを再度実行すると,置き換えられる.

非線形モデル

非線形モデルを解くには,非線形最適化問題に対応したソルバが必要である.いまは,Ipoptを用いることにし,次を実行しておく: