PuLPを使ってみる¶

特に実装すべきアルゴリズムはなかったが、せっかくなのでpythonで線形計画問題を解くためのパッケージ(?)であるPuLPを使って実際に線形計画問題を解く方法を調べてみた。

インストール方法は、

sudo pip install pulp

In [2]:

from pulp import *

$\max 100x + 100y$

$s.t. \\ x + 2y \leq 16$

$x + y \leq 18$

In [6]:

m = LpProblem(sense=LpMaximize) # 数理モデル
x = LpVariable('x', lowBound=0) # 変数
y = LpVariable('y', lowBound=0)
m += 100 * x + 100 * y # 目的関数
m += x + 2 * y <= 16 #制約1
m += 3 * x + y <= 18 # 制約2
m.solve() # ソルバー実行
print(value(x), value(y))

4.0 6.0

In [7]:

print(m.objective)

100*x + 100*y

In [8]:

print(value(m.objective))

1000.0

In [9]:

status = m.solve()

In [11]:

print (LpStatus[status])

Optimal

PuLP with pandas¶

輸送最適化問題を解く（参考サイトに問題の詳細あり）

In [12]:

import numpy as np, pandas as pd
from itertools import product
from pulp import *

np.random.seed(1)
nw, nf = 3, 4
pr = list(product(range(nw),range(nf)))
supply = np.random.randint(30, 50, nw)
demand = np.random.randint(20, 40, nf)
cost = np.random.randint(10, 20, (nw,nf))

In [13]:

# pandasを使わずにモデリングすると

m1 = LpProblem() # minimizing
v1 = {(i,j):LpVariable('v%d_%d'%(i,j), lowBound=0) for i,j in pr}
m1 += lpSum(cost[i][j] * v1[i,j] for i,j in pr) # objective function
for i in range(nw):
    m1 += lpSum(v1[i,j] for j in range(nf)) <= supply[i]
for j in range(nf):
    m1 += lpSum(v1[i,j] for i in range(nw)) >= demand[j]
m1.solve()
{k:value(x) for k,x in v1.items() if value(x) > 0}

Out[13]:

{(0, 0): 28.0,
 (0, 1): 7.0,
 (1, 2): 31.0,
 (1, 3): 5.0,
 (2, 1): 22.0,
 (2, 3): 20.0}

In [14]:

v1

Out[14]:

{(0, 0): v0_0,
 (0, 1): v0_1,
 (0, 2): v0_2,
 (0, 3): v0_3,
 (1, 0): v1_0,
 (1, 1): v1_1,
 (1, 2): v1_2,
 (1, 3): v1_3,
 (2, 0): v2_0,
 (2, 1): v2_1,
 (2, 2): v2_2,
 (2, 3): v2_3}

In [16]:

v1[1,1]

Out[16]:

v1_1

In [17]:

v1.items()

Out[17]:

dict_items([((0, 1), v0_1), ((1, 2), v1_2), ((0, 0), v0_0), ((1, 1), v1_1), ((2, 1), v2_1), ((0, 2), v0_2), ((2, 0), v2_0), ((1, 3), v1_3), ((2, 3), v2_3), ((2, 2), v2_2), ((1, 0), v1_0), ((0, 3), v0_3)])

In [18]:

value(v1[1,1])

Out[18]:

0.0

In [20]:

# pandasを使ってモデリング
a = pd.DataFrame([(i,j) for i, j in pr], columns=['warehouse', 'factory'])
a['cost'] = cost.flatten()
a[:3]

Out[20]:

	factory	cost
0	0	10
1	1	10
2	2	11

In [31]:

m2 = LpProblem()
a['Var'] = [LpVariable('v%d'%i, lowBound=0) for i in a.index]
m2 += lpDot(a.cost, a.Var)
for k, v in a.groupby('warehouse'):
    m2 += lpSum(v.Var) <= supply[k]
for k, v in a.groupby('factory'):
    m2 += lpSum(v.Var) >= demand[k]
m2.solve()
a['Val'] = a.Var.apply(value)
a[a.Val > 0]

Out[31]:

	warehouse	factory	cost	Var	Val
0	0	0	10	v0	28.0
1	0	1	10	v1	7.0
6	1	2	12	v6	31.0
7	1	3	14	v7	5.0
9	2	1	12	v9	22.0
11	2	3	12	v11	20.0

In [23]:

for i in a.index:
    print (i)

In [32]:

for k, v in a.groupby('warehouse'):
    print (k)
    print (v.Var)

0
0    v0
1    v1
2    v2
3    v3
Name: Var, dtype: object
1
4    v4
5    v5
6    v6
7    v7
Name: Var, dtype: object
2
8      v8
9      v9
10    v10
11    v11
Name: Var, dtype: object

In [ ]: