PuLPで解く線形計画法と連立方程式とナップサック問題

線形計画法 (Linear Programming, LP)を解くためのPythonライブラリとしてPuLPがあるので、これの利用例として、簡単な線形計画法と、連立方程式とナップサック問題を解いた例を以下記します。

インストール

!pip install pulp

Collecting pulp
[?25l  Downloading https://files.pythonhosted.org/packages/14/c4/0eec14a0123209c261de6ff154ef3be5cad3fd557c084f468356662e0585/PuLP-2.4-py3-none-any.whl (40.6MB)
[K     |████████████████████████████████| 40.6MB 104kB/s 
[?25hCollecting amply>=0.1.2
  Downloading https://files.pythonhosted.org/packages/f3/c5/dfa09dd2595a2ab2ab4e6fa7bebef9565812722e1980d04b0edce5032066/amply-0.1.4-py3-none-any.whl
Requirement already satisfied: docutils>=0.3 in /usr/local/lib/python3.7/dist-packages (from amply>=0.1.2->pulp) (0.17.1)
Requirement already satisfied: pyparsing in /usr/local/lib/python3.7/dist-packages (from amply>=0.1.2->pulp) (2.4.7)
Installing collected packages: amply, pulp
Successfully installed amply-0.1.4 pulp-2.4

基本的な使い方

目的の関数を最小 (あるいは最大) にする変数の値を求める「数理計画法 (Mathematical Programing)」のうち、制約条件と目的関数が線形方程式で表される問題を「線形計画法 (Linear Programming, LP)」と呼びます。

PuLPではまず、次のようにして線形計画法のためのインスタンスを作成します。

import pulp

prob = pulp.LpProblem(name='A_simple_Linear_Programming', sense=pulp.LpMaximize)

この問題で取り扱う変数を次のように定義します。ここでは、最低値がゼロという制約を加えています。

x = pulp.LpVariable('x', lowBound = 0)
y = pulp.LpVariable('y', lowBound = 0)

最大化または最小化したい目的関数を設定します。最大なのか最小なのかは、インスタンス作成時に sense=pulp.LpMaximize のようにして定義済みです。

prob += x + y

制約条件を設定します。

prob += x + 2 * y <= 6
prob += 2 * x + y <= 6

これまで設定した制約条件をチェックしましょう。

prob

A_simple_Linear_Programming:
MAXIMIZE
1*x + 1*y + 0
SUBJECT TO
_C1: x + 2 y <= 6

_C2: 2 x + y <= 6

VARIABLES
x Continuous
y Continuous

次の文を実行すればソルバーが頑張って解を探してくれます。

status = prob.solve()

無事、解が見つかったかどうかは次のようにして分かります。

pulp.LpStatus[status]

'Optimal'

最適解の各変数の値は次のようにして確認できます。

print("Optimal x=", x.value())
print("Optimal y=", y.value())
print("Optimal z=", prob.objective.value())

Optimal x= 2.0
Optimal y= 2.0
Optimal z= 4.0

連立方程式

PuLP は連立方程式にも用いることができます。例としてつぎのような連立方程式を解いてみましょう。

3 x_1 + 2 x_2 + 7 x_3 + x_4 = 8 \\
x_1 + 5 x_2 + x_3 - x_4 = 5 \\
4 x_1 + x_2 + 3 x_3 - 2 x_4 = 7 \\
x_1 + 6 x_2 + 4 x_3 + 3 x_4 = 13

変数の定義

import pulp

x_1 = pulp.LpVariable('x_1')
x_2 = pulp.LpVariable('x_2')
x_3 = pulp.LpVariable('x_3')
x_4 = pulp.LpVariable('x_4')

制約条件の設定

prob = pulp.LpProblem(name='Simultaneous_linear_equations')
prob += 3 * x_1 + 2 * x_2 + 7 * x_3 + x_4 == 8
prob += x_1 + 5 * x_2 + x_3 - x_4 == 5
prob += 4 * x_1 + x_2 + 3 * x_3 - 2 * x_4 == 7
prob += x_1 + 6 * x_2 + 4 * x_3 + 3 * x_4 == 13

制約条件の確認

prob

Simultaneous_linear_equations:
MINIMIZE
None
SUBJECT TO
_C1: 3 x_1 + 2 x_2 + 7 x_3 + x_4 = 8

_C2: x_1 + 5 x_2 + x_3 - x_4 = 5

_C3: 4 x_1 + x_2 + 3 x_3 - 2 x_4 = 7

_C4: x_1 + 6 x_2 + 4 x_3 + 3 x_4 = 13

VARIABLES
x_1 free Continuous
x_2 free Continuous
x_3 free Continuous
x_4 free Continuous

ソルバー実行

status = prob.solve()
pulp.LpStatus[status]

'Optimal'

解の確認。

x_1.value(), x_2.value(), x_3.value(), x_4.value()

(3.5, 1.0, -1.0, 2.5)

ナップサック問題

アイテムの集合 $I$ = {$1, 2, ..., N$} 中の各アイテム $i \in I$ の価値を $v_i$、重みを $w_i$とし、それを入れるナップサックの最大容量を $W_{capacity}$ としたとき、下記の最適化問題をナップサック問題と言います。

max \sum_{i=1}^{N} v_i x_i \\
s.t. 
\sum_{i=1}^{N} w_i x_i \le W_{capacity} \\
\qquad \qquad x_i \in \mathbb{N} \qquad (\forall i \in I)

• $\in$ は英語で「in」、 $\forall$ は「for all」の意味
• $s.t.$ は「such that」の略で「下記を満たす条件の元で上記を求める」くらいの意味
• $\mathbb{N}$ は自然数全体の集合
• $x_i \in \mathbb{N}$ を $x_i \in$ {$0$, $1$} に置き換えたものを特に「0-1ナップサック問題」と言います。

ナップサック問題は「動的計画法」で解ける問題として典型的なものですが、ここでは「線形計画法」としてPuLPで解きます。

データの用意

n =  5
W_capacity =  175
W =  [95, 56, 65, 54, 32]
V =  [79, 9, 57, 55, 27]

制約条件の設定

import pulp

prob = pulp.LpProblem('kp01', sense = pulp.LpMaximize)
xs = [pulp.LpVariable('{}'.format(x), cat='Binary', lowBound = 0) for x in range(n)]
prob += pulp.lpDot(V, xs)
prob += pulp.lpDot(W, xs) <= W_capacity
prob

kp01:
MAXIMIZE
79*0 + 9*1 + 57*2 + 55*3 + 27*4 + 0
SUBJECT TO
_C1: 95 0 + 56 1 + 65 2 + 54 3 + 32 4 <= 175

VARIABLES
0 <= 0 <= 1 Integer
0 <= 1 <= 1 Integer
0 <= 2 <= 1 Integer
0 <= 3 <= 1 Integer
0 <= 4 <= 1 Integer

ソルバーの実行

status = prob.solve()
pulp.LpStatus[status]

'Optimal'

答え

print([x.value() for x in xs])
print(prob.objective.value())

[0.0, 0.0, 1.0, 1.0, 1.0]
139.0

個数制限なしナップサック問題

$x_i \in \mathbb{N}$ をゼロ以上の自然数に拡張したものを、個数制限なしのナップサック問題と言います。動的計画法では、0-1ナップサック問題と個数制限なしのナップサック問題とでアルゴリズムを変更する必要が出てきますが、PuLPでは、次のように少しだけ変更すればOKです。

データの用意

n = 3
W_capacity = 7
W = [3, 4, 2]
V = [4, 5, 3]

制約条件の設定

import pulp

prob = pulp.LpProblem('kp_unlimited', sense = pulp.LpMaximize)
xs = [pulp.LpVariable('{}'.format(x), cat='Integer', lowBound = 0) for x in range(n)]
prob += pulp.lpDot(V, xs)
prob += pulp.lpDot(W, xs) <= W_capacity
prob

kp_unlimited:
MAXIMIZE
4*0 + 5*1 + 3*2 + 0
SUBJECT TO
_C1: 3 0 + 4 1 + 2 2 <= 7

VARIABLES
0 <= 0 Integer
0 <= 1 Integer
0 <= 2 Integer

ソルバーの実行

status = prob.solve()
pulp.LpStatus[status]

'Optimal'

答え

print([x.value() for x in xs])
print(prob.objective.value())

[1.0, 0.0, 2.0]
10.0