原油の製油の最適化

これなに

「OR-Toolsで学ぶ最適化　【線形計画：製油しようぜ】」をpandasを使ったモデルで解き直した記事です。

定式化

変数：原油種類ごと、製品種類ごとの作成量
目的関数：購入費用合計 ー 販売額合計 → 最大化
制約条件
・原油の購入上限
・製品ごと需要の上下限
・製品のオクタン価

Pythonで解く

データ設定

原油と製品の表を作成

from io import StringIO
from ortoolpy import pd, model_max, lpDot, lpSum

dfoil = pd.read_csv(StringIO("""\
原油,原油_オクタン価,所有量,費用_バレル
R0,99,782,55.34
R1,94,894,54.12
R2,84,631,53.68
R3,92,648,57.03
R4,87,956,54.81
R5,97,647,56.25
R6,81,689,57.55
R7,96,609,58.21"""))

dfprd = pd.read_csv(StringIO("""\
製品,製品_オクタン価,需要下限,需要上限,売価
F0,88,415,11707,61.97
F1,94,199,7761,62.04
F2,90,479,12596,61.99"""))

変数表作成

変数表は、原油と製品の直積で作成

df = pd.DataFrame(
    [
        [*row1, *row2]
        for row1 in dfoil.itertuples(False)
        for row2 in dfprd.itertuples(False)
    ],
    columns=dfoil.columns.tolist() + dfprd.columns.tolist(),
)
df[:2]  # 先頭2行

	原油	原油_オクタン価	所有量	費用_バレル	製品	製品_オクタン価	需要下限	需要上限	売価
0	R0	99	782	55.34	F0	88	415	11707	61.97
1	R0	99	782	55.34	F1	94	199	7761	62.04

モデル化＆結果表示

m = model_max(df=df)
for _, dfs in df.groupby("原油"):
    m += lpSum(dfs.Var) <= dfs.所有量.iloc[0]
for _, dfs in df.groupby("製品"):
    m += lpSum(dfs.Var) >= dfs.需要下限.iloc[0]
    m += lpSum(dfs.Var) <= dfs.需要上限.iloc[0]
    m += (lpDot(dfs.原油_オクタン価, dfs.Var) ==
          lpDot(dfs.製品_オクタン価, dfs.Var))
m.solve(objs=["-費用_バレル", "売価"])

dfr = df.pivot_table("Val", "原油", "製品").round(2)
dfr = pd.concat([dfr, dfr.sum().to_frame("計").T])
dfr["計"] = dfr.sum(1)
dfr

	F0	F1	F2	計
R0	0	782	0	782
R1	894	0	0	894
R2	465.29	30.08	135.63	631
R3	330.04	0	317.96	648
R4	0	956	0	956
R5	0	621.59	25.41	647
R6	689	0	0	689
R7	0	609	0	609
計	2378.33	2998.67	479	5856

（何故かOR-Tools の元記事とは違いますね）

モデル化のポイント

• ortoolpyの0.2.36以降で可能です。
• 最大化モデルはmodel_max、最小化モデルはmodel_minで作成します。
  • モデル作成時に、下記の引数にDataFrameを指定すると変数としてVar列を作成します。
    • df：非負連続変数
    • dfb：0-1変数
    • dfi：非負整数変数
• solve時にobjsで列名のリストを指定すると、その列を目的関数の係数として使用します。
  • 列名に「-」がついていると、-1倍します。
• モデルにDataFrameが結びついていると、solve後に結果をVal列として追加します。

以上