Closed11

簡単な問題でPuLPを触ってみる

nabeyangnabeyang

この問題の1問目を解いてみます(2問目も数値が違うだけで同じ問題)。
https://www.fit.ac.jp/~song/lecture/06/LP0606.pdf
変数名がプログラム上で、不便なので製品Ⅰはp1で製品Ⅱはp2(pはproductのpという気持ち)、原料E, F, Gはそれぞれm1, m2, m3とします(mはmaterialのm)。

このスクラップでやりたいことは、徐々にプログラムを改善していって、製品や原料の種類を変えた場合も変更無しに動くプログラムにしていくことです。

nabeyangnabeyang

全てベタ書きすると、こんな感じになります。

import pulp

problem = pulp.LpProblem(sense=pulp.LpMaximize)

# * 変数
# x['p{i}']: 製品p{i}の製造数
x = {}
x['p1'] = pulp.LpVariable('x_p1', cat='Continuous')
x['p2'] = pulp.LpVariable('x_p2', cat='Continuous')

# * 制約(1)
# 製品p{i}の製造数は0以上
problem += x['p1'] >= 0
problem += x['p2'] >= 0

# * 制約(2)

# 原料m1の使用可能料は8以下
# 製品p1を1kg作るのに原料m1は1kg必用
# 製品p2を1kg作るのに原料m1は2kg必用
problem += x['p1'] * 1 + x['p2'] * 2 <= 8

# 原料m2の使用可能料は12以下
# 製品p1を1kg作るのに原料m2は3kg必用
# 製品p2を1kg作るのに原料m2は0kg必用
problem += x['p1'] * 3 + x['p2'] * 0 <= 12

# 原料m3の使用可能料は12以下
# 製品p1を1kg作るのに原料m3は0kg必用
# 製品p2を1kg作るのに原料m3は4kg必用
problem += x['p1'] * 0 + x['p2'] * 4 <= 12

# * 目的関数
# 製品p1は1kgあたり4万円の利益
# 製品p2は1kgあたり6万円の利益
# 利益の最大値を求める
problem += x['p1'] * 4 + x['p2'] * 6

status = problem.solve()


# 結果
print(f"status = {pulp.LpStatus[status]}")
print(f"x['p1'] = {x['p1'].value()}") #=> 4
print(f"x['p2'] = {x['p2'].value()}") #=> 2
print(f"利益の最大値 = {problem.objective.value()}万円")#=> 28
nabeyangnabeyang

ここからプログラムらしく、製品や原料の数によらずに動くように書き換えていきます。
まず、パラメータを導入します。

--- main1.py    2021-10-10 18:28:02.000000000 +0900
+++ main2.py    2021-10-10 18:27:24.000000000 +0900
@@ -2,6 +2,24 @@
 
 problem = pulp.LpProblem(sense=pulp.LpMaximize)
 
+# パラメータ
+# P: 製品一覧
+P = ['p1', 'p2']
+# M: 原料一覧
+M = ['m1', 'm2', 'm3']
+# stock: 在庫
+stock = {'m1': 8, 'm2': 12, 'm3': 12}
+# require['p{i}', 'm{j}']: 製品p{i}を1kg作るのに必用な材料m{j}の重量
+require = {
+    ('p1', 'm1'): 1,
+    ('p1', 'm2'): 3,
+    ('p1', 'm3'): 0,
+    ('p2', 'm1'): 2,
+    ('p2', 'm2'): 0,
+    ('p2', 'm3'): 4,
+}
+# gain: 製品p{i}を1kg作ると得られる利益
+gain = {'p1': 4, 'p2': 6}
 # * 変数
 # x['p{i}']: 製品p{i}の製造数
 x = {}
nabeyangnabeyang

変数使って、プログラムを書き直します。ただし、ループは使わずまだベタ書きでやります。

--- main2.py	2021-10-10 18:27:24.000000000 +0900
+++ main3.py	2021-10-10 18:44:18.000000000 +0900
@@ -33,26 +33,26 @@
 
 # * 制約(2)
 
-# 原料m1の使用可能料は8以下
-# 製品p1を1kg作るのに原料m1は1kg必用
-# 製品p2を1kg作るのに原料m1は2kg必用
-problem += x['p1'] * 1 + x['p2'] * 2 <= 8
-
-# 原料m2の使用可能料は12以下
-# 製品p1を1kg作るのに原料m2は3kg必用
-# 製品p2を1kg作るのに原料m2は0kg必用
-problem += x['p1'] * 3 + x['p2'] * 0 <= 12
-
-# 原料m3の使用可能料は12以下
-# 製品p1を1kg作るのに原料m3は0kg必用
-# 製品p2を1kg作るのに原料m3は4kg必用
-problem += x['p1'] * 0 + x['p2'] * 4 <= 12
+# 原料m1の使用可能料はstock['m1']以下
+# 製品p1を1kg作るのに原料m1はrequire['p1', 'm1']kg必用
+# 製品p2を1kg作るのに原料m1はrequire['p2', 'm1']kg必用
+problem += x['p1'] * require['p1', 'm1'] + x['p2'] * require['p2', 'm1'] <= stock['m1']
+
+# 原料m2の使用可能料はstock['m2']以下
+# 製品p1を1kg作るのに原料m2はrequire['p1', 'm2']kg必用
+# 製品p2を1kg作るのに原料m2はrequire['p2', 'm2']kg必用
+problem += x['p1'] * require['p1', 'm2'] + x['p2'] * require['p2', 'm2'] <= stock['m2']
+
+# 原料m3の使用可能料はstock['m3']以下
+# 製品p1を1kg作るのに原料m3はrequire['p1', 'm3']kg必用
+# 製品p2を1kg作るのに原料m3はrequire['p2', 'm3']kg必用
+problem += x['p1'] * require['p1', 'm3'] + x['p2'] * require['p2', 'm3'] <= stock['m3']
 
 # * 目的関数
-# 製品p1は1kgあたり4万円の利益
-# 製品p2は1kgあたり6万円の利益
+# 製品p1は1kgあたりgain['p1']万円の利益
+# 製品p2は1kgあたりgain['p2']万円の利益
 # 利益の最大値を求める
-problem += x['p1'] * 4 + x['p2'] * 6
+problem += x['p1'] * gain['p1'] + x['p2'] * gain['p2']
 
 status = problem.solve()
nabeyangnabeyang

ループや関数を使って、書き直していきます。
pulp.LpVariable.dictsを使うと、xは次のように1行で書けます。

--- main3.py	2021-10-10 18:44:18.000000000 +0900
+++ main4.py	2021-10-10 18:47:49.000000000 +0900
@@ -22,9 +22,7 @@
 gain = {'p1': 4, 'p2': 6}
 # * 変数
 # x['p{i}']: 製品p{i}の製造数
-x = {}
-x['p1'] = pulp.LpVariable('x_p1', cat='Continuous')
-x['p2'] = pulp.LpVariable('x_p2', cat='Continuous')
+x = pulp.LpVariable.dicts('x', P, cat='Continuous')
 
 # * 制約(1)
 # 製品p{i}の製造数は0以上

制約(1)はforループで書き直します。

--- main4.py	2021-10-10 18:47:49.000000000 +0900
+++ main5.py	2021-10-10 18:51:53.000000000 +0900
@@ -26,8 +26,8 @@
 
 # * 制約(1)
 # 製品p{i}の製造数は0以上
-problem += x['p1'] >= 0
-problem += x['p2'] >= 0
+for p in P:
+    problem += x[p] >= 0
 
 # * 制約(2)

制約(2)はまず原料Mのループを使うと次のようになります。

--- main5.py	2021-10-10 18:51:53.000000000 +0900
+++ main6.py	2021-10-10 21:27:45.000000000 +0900
@@ -31,20 +31,11 @@
 
 # * 制約(2)
 
-# 原料m1の使用可能料はstock['m1']以下
-# 製品p1を1kg作るのに原料m1はrequire['p1', 'm1']kg必用
-# 製品p2を1kg作るのに原料m1はrequire['p2', 'm1']kg必用
-problem += x['p1'] * require['p1', 'm1'] + x['p2'] * require['p2', 'm1'] <= stock['m1']
-
-# 原料m2の使用可能料はstock['m2']以下
-# 製品p1を1kg作るのに原料m2はrequire['p1', 'm2']kg必用
-# 製品p2を1kg作るのに原料m2はrequire['p2', 'm2']kg必用
-problem += x['p1'] * require['p1', 'm2'] + x['p2'] * require['p2', 'm2'] <= stock['m2']
-
-# 原料m3の使用可能料はstock['m3']以下
-# 製品p1を1kg作るのに原料m3はrequire['p1', 'm3']kg必用
-# 製品p2を1kg作るのに原料m3はrequire['p2', 'm3']kg必用
-problem += x['p1'] * require['p1', 'm3'] + x['p2'] * require['p2', 'm3'] <= stock['m3']
+for m in M:
+    # 原料mの使用可能料はstock[m]以下
+    # 製品p1を1kg作るのに原料mはrequire['p1', m]kg必用
+    # 製品p2を1kg作るのに原料mはrequire['p2', m]kg必用
+    problem += x['p1'] * require['p1', m] + x['p2'] * require['p2', m] <= stock[m]
 
 # * 目的関数
 # 製品p1は1kgあたりgain['p1']万円の利益

結果についてもループで書き直します。

--- main6.py	2021-10-10 18:56:03.000000000 +0900
+++ main7.py	2021-10-10 18:58:22.000000000 +0900
@@ -48,6 +48,6 @@
 
 # 結果
 print(f"status = {pulp.LpStatus[status]}")
-print(f"x['p1'] = {x['p1'].value()}")
-print(f"x['p2'] = {x['p2'].value()}")
+for p in P:
+    print(f"x['{p}'] = {x[p].value()}")
 print(f"利益の最大値 = {problem.objective.value()}万円")
nabeyangnabeyang

制約(2)と目的関数はpulp.lpSumPのリスト内記法できれいに書き直せます。

--- main7.py	2021-10-10 21:27:18.000000000 +0900
+++ main8.py	2021-10-10 21:26:45.000000000 +0900
@@ -33,15 +33,13 @@
 
 for m in M:
     # 原料mの使用可能料はstock[m]以下
-    # 製品p1を1kg作るのに原料mはrequire['p1', m]kg必用
-    # 製品p2を1kg作るのに原料mはrequire['p2', m]kg必用
-    problem += x['p1'] * require['p1', m] + x['p2'] * require['p2', m] <= stock[m]
+    # 製品pを1kg作るのに原料m1はrequire[p, m]kg必用
+    problem += pulp.lpSum([x[p] * require[p, m] for p in P]) <= stock[m]
 
 # * 目的関数
-# 製品p1は1kgあたりgain['p1']万円の利益
-# 製品p2は1kgあたりgain['p2']万円の利益
+# 製品pは1kgあたりgain[p]万円の利益
 # 利益の最大値を求める
-problem += x['p1'] * gain['p1'] + x['p2'] * gain['p2']
+problem += pulp.lpSum([x[p] * gain[p] for p in P])
 
 status = problem.solve()
nabeyangnabeyang

ここからstock, require, gainを外部のcsvファイルから読み取るように変更します。
まずstockは次のように変えます。

--- main8.py	2021-10-10 21:06:20.000000000 +0900
+++ main9.py	2021-10-10 21:08:23.000000000 +0900
@@ -1,5 +1,7 @@
+import pandas as pd
 import pulp
 
+stock_df = pd.read_csv('./data/stocks.csv')
 problem = pulp.LpProblem(sense=pulp.LpMaximize)
 
 # パラメータ
@@ -8,7 +10,7 @@
 # M: 原料一覧
 M = ['m1', 'm2', 'm3']
 # stock: 在庫
-stock = {'m1': 8, 'm2': 12, 'm3': 12}
+stock = {row.m: row.stock for row in stock_df.itertuples()}
 # require['p{i}', 'm{j}']: 製品p{i}を1kg作るのに必用な材料m{j}の重量
 require = {
     ('p1', 'm1'): 1,

ここでstocks.csvは次の通りです。

m,stock
m1,8
m2,12
m3,12
nabeyangnabeyang

requiresも同様に書き換えます。

--- main9.py	2021-10-10 21:08:23.000000000 +0900
+++ main10.py	2021-10-10 21:14:38.000000000 +0900
@@ -2,6 +2,7 @@
 import pulp
 
 stock_df = pd.read_csv('./data/stocks.csv')
+require_df = pd.read_csv('./data/requires.csv')
 problem = pulp.LpProblem(sense=pulp.LpMaximize)
 
 # パラメータ
@@ -12,14 +13,7 @@
 # stock: 在庫
 stock = {row.m: row.stock for row in stock_df.itertuples()}
 # require['p{i}', 'm{j}']: 製品p{i}を1kg作るのに必用な材料m{j}の重量
-require = {
-    ('p1', 'm1'): 1,
-    ('p1', 'm2'): 3,
-    ('p1', 'm3'): 0,
-    ('p2', 'm1'): 2,
-    ('p2', 'm2'): 0,
-    ('p2', 'm3'): 4,
-}
+require = { (row.p, row.m): row.require for row in require_df.itertuples()}
 # gain: 製品p{i}を1kg作ると得られる利益
 gain = {'p1': 4, 'p2': 6}
 # * 変数

ここで、requires.csvの内容は次の通り。

p,m,require
p1,m1,1
p1,m2,3
p1,m3,0
p2,m1,2
p2,m2,0
p2,m3,4
nabeyangnabeyang

gainも同様です。

--- main10.py	2021-10-10 21:14:38.000000000 +0900
+++ main11.py	2021-10-10 21:18:00.000000000 +0900
@@ -3,6 +3,7 @@
 
 stock_df = pd.read_csv('./data/stocks.csv')
 require_df = pd.read_csv('./data/requires.csv')
+gain_df = pd.read_csv('./data/gains.csv')
 problem = pulp.LpProblem(sense=pulp.LpMaximize)
 
 # パラメータ
@@ -15,7 +16,7 @@
 # require['p{i}', 'm{j}']: 製品p{i}を1kg作るのに必用な材料m{j}の重量
 require = { (row.p, row.m): row.require for row in require_df.itertuples()}
 # gain: 製品p{i}を1kg作ると得られる利益
-gain = {'p1': 4, 'p2': 6}
+gain = {row.p: row.gain for row in gain_df.itertuples()}
 # * 変数
 # x['p{i}']: 製品p{i}の製造数
 x = pulp.LpVariable.dicts('x', P, cat='Continuous')
p,gain
p1,4
p2,6
nabeyangnabeyang

最後にP, Mを次のように書き換えれば、製品や原料の数に依らないプログラムになります。

--- main11.py	2021-10-10 21:18:00.000000000 +0900
+++ main12.py	2021-10-10 21:20:04.000000000 +0900
@@ -8,9 +8,9 @@
 
 # パラメータ
 # P: 製品一覧
-P = ['p1', 'p2']
+P = gain_df['p'].tolist()
 # M: 原料一覧
-M = ['m1', 'm2', 'm3']
+M = stock_df['m'].tolist()
 # stock: 在庫
 stock = {row.m: row.stock for row in stock_df.itertuples()}
 # require['p{i}', 'm{j}']: 製品p{i}を1kg作るのに必用な材料m{j}の重量
nabeyangnabeyang

最終的なプログラムは、次のようになります。

import pandas as pd
import pulp

stock_df = pd.read_csv('./data/stocks.csv')
require_df = pd.read_csv('./data/requires.csv')
gain_df = pd.read_csv('./data/gains.csv')
problem = pulp.LpProblem(sense=pulp.LpMaximize)

# パラメータ
# P: 製品一覧
P = gain_df['p'].tolist()
# M: 原料一覧
M = stock_df['m'].tolist()
# stock: 在庫
stock = {row.m: row.stock for row in stock_df.itertuples()}
# require['p{i}', 'm{j}']: 製品p{i}を1kg作るのに必用な材料m{j}の重量
require = { (row.p, row.m): row.require for row in require_df.itertuples()}
# gain: 製品p{i}を1kg作ると得られる利益
gain = {row.p: row.gain for row in gain_df.itertuples()}
# * 変数
# x['p{i}']: 製品p{i}の製造数
x = pulp.LpVariable.dicts('x', P, cat='Continuous')

# * 制約(1)
# 製品p{i}の製造数は0以上
for p in P:
    problem += x[p] >= 0

# * 制約(2)

for m in M:
    # 原料mの使用可能料はstock[m]以下
    # 製品pを1kg作るのに原料mはrequire[p, m]kg必用
    problem += pulp.lpSum([x[p] * require[p, m] for p in P]) <= stock[m]

# * 目的関数
# 製品pは1kgあたりgain[p]万円の利益
# 利益の最大値を求める
problem += pulp.lpSum([x[p] * gain[p] for p in P])

status = problem.solve()


# 結果
print(f"status = {pulp.LpStatus[status]}")
for p in P:
    print(f"x['{p}'] = {x[p].value()}")
print(f"利益の最大値 = {problem.objective.value()}万円")
このスクラップは2021/10/10にクローズされました