[Package] Pulp (Python) と JuMP (Julia) の記法比較アンチョコ (2022年7月版)

はじめに

タイトルの通りですが、Pythonのpulpでモデルを書くときと、JuliaのJuMPでモデルを書くときの両方があります。毎回pulpとJuMPのドキュメントを見に行くのが大変になってきたので、自分のために比較を書いておくことにしました。

以下はこの記事を書くために使用したバージョンです。

	バージョン	備考
Python	3.7.13	Google colaboratory上
pulp	2.6.0	Google colaboratory上
Julia	1.6.7	macOS (Apple M1/Memory 8GBマシン上)
JuMP	1.1.1	-
MathOptInterface	1.6.1	終了判定用

例題1.ナップサック問題

サイゼリヤ問題などでよく使われているあれです。例題はこちらのサイトの例題からパクりました。結果は価値17.0、解が$\{0,3,4\}$です。

Python+PuLP

単体の変数 P.LpVariable("x", 0, 1, "Integer") を使ってもいいですけど、かなりの場面で P.LpVariable.dicts で使う気がしています。

変数宣言に関するドキュメントです: classpulp.LpVariable(name, lowBound=None, upBound=None, cat='Continuous', e=None)

以前、lpDot や lpSum を使うほうが高速と聞いたことがありますが、現代ではどうか知らないです。

import pulp as P

def solve(W, V, Capacity):
    N = len(W)
    model = P.LpProblem(sense=P.LpMaximize)
    x = P.LpVariable.dicts("x", range(N), 0, 1, "Integer")
    solver = P.PULP_CBC_CMD(msg=0)

    ## If we use `lpDot`
    # model += P.lpDot(V, x)
    # model += P.lpDot(W, x) <= Capacity

    # lpSum
    model += P.lpSum(V[i] * x[i] for i in range(N))
    model += P.lpSum(W[i] * x[i] for i in range(N)) <= Capacity
    
    model.solve(solver)
    
    obj = P.value(model.objective)
    resx = [P.value(x[i]) for i in range(N)]
    print(P.LpStatus[model.status])
    print(P.LpStatus)
    if P.LpStatus[model.status] == "Optimal":
        print(obj)
        print(resx)
    
# 例題: https://dev2prod.site/python/napsack/
W = [4, 2, 2, 1, 10]
V = [12, 2, 1, 1, 4]
Capacity = 15
solve(W, V, Capacity)

Julia+JuMP

Juliaの場合、for文を回しても高速と信じられているので、気にせずに数式をfor文で書いても良いです（たぶん）。Julia＋JuMPで最適解が求まったかどうかを判定したい場合には、MathOptInterfaceを使ってください。例えば ts=1e-5 などを入れると、TIME_LIMIT になります。

ソルバーへのパラメータは set_optimizer_attribute という名前に変わりました。また内部のソルバー（ここではCbc）のバージョンによって、渡す名前が変わったりするので、ソルバーのオプションを眺めて渡す必要があります。

using JuMP
using Cbc
using MathOptInterface

# 関数
function solve(W, V, Capacity; ts=1)
    N = length(W)
    model = Model(Cbc.Optimizer)
    set_optimizer_attribute(model, "log", 0)
    set_optimizer_attribute(model, "seconds", ts)
    
    @variable(model, x[1:N], Bin)
    @objective(model, Max, sum(V[i] * x[i] for i in 1:N))
    @constraint(model, sum(W[i] * x[i] for i in 1:N) <= Capacity)
    optimize!(model)
    
    println(MOI.get(model, MOI.TerminationStatus())) # term. status
    if MOI.get(model, MOI.TerminationStatus()) == MOI.OPTIMAL # check
        obj = objective_value(model)
        resx = value.(x)
        println(obj)
        println(resx)
    end
end

# 例題: https://dev2prod.site/python/napsack/
W = [4, 2, 2, 1, 10]
V = [12, 2, 1, 1, 4]
Capacity = 15
solve(W, V, Capacity)

例題2. X

他に自分が見て役に立ちそうな例が思いついたら、更新していきます（例えばPythonのnetworkxを使う場合 vs JuliaのGraphs.jlを使う場合などです）。

まとめ

自動でお互い変換してくれ！（終わりです）。