[数理最適化 / pyomo] 無料ツールで非線形問題を解きたい時 on Google Colaboratory

環境

Google Colaboratory に glpk と pyomo をインストールして最適化コードを実行していた。

!apt-get install -y glpk-utils
!pip install -q pyomo

発生した問題

ところが、非線形な問題を解こうとしたら、以下のようなエラーが出た。(つ'ヮ'c)ｳｯﾋｮｫｫｫｵｵｱｱｱｱｱ！

pyomo ValueError: Model constraint (C5) contains nonlinear terms that cannot be written to LP format

glpk では、非線形な問題は解けないらしい>_<。

対処法

ipopt というソルバーを使えば非線形問題も解けたので、手順を以下に記載する。

1. ipopt のインストール

Google Colaboratory

!wget -N -q "https://matematica.unipv.it/gualandi/solvers/ipopt-linux64.zip"
!unzip -o -q ipopt-linux64

これで ipopt が使えるようになる。

2. 実行コード

Google Colaboratory

import pyomo.environ as pyo
from pyomo.opt import SolverFactory


model = pyo.ConcreteModel()

# ----------------------------------------
# パラメータ
# ----------------------------------------
model.setInv = pyo.Set(initialize=["productA","productB","productC"])
model.Demand = 30000
M = 1e7

# ----------------------------------------
# 決定変数
# ----------------------------------------
model.Cost = pyo.Var(model.setInv, bounds=(0, None))
model.Product = pyo.Var(model.setInv, bounds=(0, None), domain=pyo.NonNegativeIntegers)
model.Binary = pyo.Var(model.setInv, within=pyo.Binary)

# ----------------------------------------
# 目的関数
# ----------------------------------------
model.obj = pyo.Objective(expr = pyo.summation(model.Cost), sense=pyo.minimize)

# ----------------------------------------
# 制約条件
# ----------------------------------------
model.C1 = pyo.Constraint(expr = pyo.summation(model.Product) == model.Demand)
# 非線形な制約
model.C2 = pyo.Constraint(
    expr = model.Cost["productA"] == (
        0.02 * (model.Product["productA"] ** 3)
        + 0.7 * model.Product["productA"]
        + model.Binary["productA"] * 0.3
    )
)
model.C3 = pyo.Constraint(expr = model.Cost["productB"] == 0.2 * model.Product["productB"] + model.Binary["productB"] * 0.7)
# 非線形な制約
model.C4 = pyo.Constraint(expr = model.Cost["productC"] == 0.13 * (model.Product["productC"] ** 2))
model.C5 = pyo.Constraint(expr = model.Product["productA"] <= model.Binary["productA"] * M)
model.C6 = pyo.Constraint(expr = model.Product["productB"] <= model.Binary["productB"] * M)

# ----------------------------------------
# 求解
# ----------------------------------------
# opt = SolverFactory('glpk')
opt = SolverFactory('ipopt')

model.results = opt.solve(model)

model.pprint()
print("Result objective: ", pyo.value(model.obj))

でけた！( "´༥`" )

参考資料

• 非線形問題を解きたいときに使える solver の解説

• 非線形問題を解くときに使える ipopt を Google Colaboratory で install する方法

• 自分は採用しなかったけど、以下のような方法もあるみたい (conda を使って install)