ナップサック問題#

JijSolverとJijModelingを用いて、ナップサック問題を解く方法を紹介します。 この問題はLucas, 2014, “Ising formulations of many NP problems” 5.2. Knapsack with Integer Weightsなどでも取り上げられています。

ナップサック問題とは#

ナップサック問題は、以下のような場合において最適解を見つける問題です。 この問題はNP困難な整数計画問題として知られています。 ナップサック問題の具体例として、以下のストーリーを考えてみましょう。

ある探検家が洞窟を冒険していると、思いがけず宝物を見つけました。

Treasure A

Treasure B

Treasure C

Treasure D

Treasure E

Treasure F

Price

$5000

$7000

$2000

$1000

$4000

$3000

weight

800g

1000g

600g

400g

500g

300g

しかし、探検家が持っているのは容量2キロの小さなナップサックだけでした。当然、探検家はこのナップサックに出来るだけ価値の高い宝物を詰め込みたいと考えます。探検家は、どの宝物をナップサックに詰め込めば良いでしょうか?

このような状況下で最適解を求めるのがナップサック問題です。ナップサック問題はNP困難な整数計画問題として有名なものの1つです。

定式化#

上の例を一般化して数理モデルを定式化してみましょう。ナップサックに入れるアイテム(宝物)の集合を \(\{ 0, 1, \dots, i, \dots, N-1 \}\) とし、各アイテム \(i\) の価値を \(v_i\)、重量を \(w_i\) と書くことにします。

\[ v = \{v_0, v_1, \dots, v_i, \dots, v_{N-1}\} \]
\[ w = \{w_0, w_1, \dots, w_i, \dots, w_{N-1}\} \]

さらに、 アイテム \(i\) を詰め込むかどうかを表現する決定変数として \(x_i\) を導入します。 \(x_i\) はアイテム \(i\) をナップサックに詰め込む場合に1を取り、そうでない場合には0になるバイナリ変数です。また、ナップサックの容量を \(W\) で表すことにします。
ナップサックに入れるアイテムの価値の合計を最大化することを目指します。 そのために、ナップサックに入れるアイテムの価値の合計を数式で表現しましょう。 そしてナップサックの容量制限を数式で表現する必要があります。 最終的に、この問題の数理モデルは以下のようになります。

\[ \max \quad \sum_{i=0}^{N-1} v_i x_i \tag{1} \]
\[ \mathrm{s.t.} \quad \sum_{i=0}^{N-1} w_i x_i \leq W \tag{2} \]
\[ x_i \in \{0, 1\} \quad (\forall i \in \{0, 1, \dots, N-1\}) \tag{3} \]

JijModelingによる定式化#

次に、JijModelingを使って上記の数理モデルを実装してみましょう。まずは、数理モデルに現れる変数およびパラメーターを定義します。

import jijmodeling as jm

v = jm.Placeholder("v", ndim=1)
N = v.len_at(0, latex="N")
w = jm.Placeholder("w", ndim=1)
W = jm.Placeholder("W")
x = jm.BinaryVar("x", shape=(N,))
i = jm.Element("i", belong_to=(0, N))

ここで v はアイテム \(i\) の価値 \(v_i\) のリスト、 w はアイテム \(i\) の重量 \(w_i\) のリストを定義しており、 N はアイテムの個数、 W はナップサックの容量を定義しています。また、\(x_i\) に相当するバイナリ変数として x を定義しています。 i は後の数理モデル構築のために用意された添字で、アイテム \(i\) を相当します。

目的関数の実装#

目的関数として式(1)を実装すると以下のようになります。

problem = jm.Problem("Knapsack", sense=jm.ProblemSense.MAXIMIZE)
problem += jm.sum(i, v[i]*x[i])

ナップサック問題は目的関数を最大化する問題であるため、 sense=jm.ProblemSense.MAXIMIZE を設定しています。

制約条件の実装#

制約条件として式(2)を実装すると以下のようになります。

problem += jm.Constraint("weight", jm.sum(i, w[i]*x[i]) <= W)

これで数理モデルの実装は完了です。正しく数理モデルが実装されているかをLaTeX表示を通して確認してみましょう。

problem
\[\begin{split}\begin{array}{cccc}\text{Problem:} & \text{Knapsack} & & \\& & \max \quad \displaystyle \sum_{i = 0}^{N - 1} v_{i} \cdot x_{i} & \\\text{{s.t.}} & & & \\ & \text{weight} & \displaystyle \sum_{i = 0}^{N - 1} w_{i} \cdot x_{i} \leq W & \\\text{{where}} & & & \\& x & 1\text{-dim binary variable}\\\end{array}\end{split}\]

インスタンスデータの準備#

次にインスタンスデータを準備します。ここでは100個のアイテムをランダムに生成し、ナップサックの容量を100とするものとします。

import numpy as np

inst_v = np.random.randint(5,30,100)
inst_w = inst_v + np.random.randint(-2,20,100)
inst_W = 100
instance_data = {"v": inst_v, "w": inst_w, "W": inst_W}

inst_v , inst_w, inst_Wはそれぞれアイテムの価値のリスト、アイテムの重量のリスト、ナップサックの容量を表します。

ナップサック問題を解く#

jijsolverを用いて、この問題を解きましょう。

import jijsolver

interpreter = jm.Interpreter(instance_data)
instance = interpreter.eval_problem(problem)
solution = jijsolver.solve(instance, time_limit_sec=1.0)
Process ForkPoolWorker-1:
Traceback (most recent call last):
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/process.py", line 315, in _bootstrap
    self.run()
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/process.py", line 108, in run
    self._target(*self._args, **self._kwargs)
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py", line 125, in worker
    result = (True, func(*args, **kwds))
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py", line 48, in mapstar
    return list(map(*args))
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/pathos/helpers/mp_helper.py", line 15, in <lambda>
    func = lambda args: f(*args)
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/jijsolver/solver.py", line 186, in _process
    return _sample_impl(
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/jijsolver/solver.py", line 46, in _sample_impl
    result = weighted_local_search_solver(
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/jijsolver/solver_components.py", line 82, in weighted_local_search_solver
    result = wls.sample(
pyo3_runtime.PanicException: In 4-flip, the variables are assumed to be all different from each other.
thread '<unnamed>' panicked at src/system/weighted_local_search_system/weighted_local_search.rs:488:9:
In 4-flip, the variables are assumed to be all different from each other.
note: run with `RUST_BACKTRACE=1` environment variable to display a backtrace
Process ForkPoolWorker-2:
Process ForkPoolWorker-3:
Traceback (most recent call last):
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/process.py", line 315, in _bootstrap
    self.run()
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/process.py", line 108, in run
    self._target(*self._args, **self._kwargs)
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py", line 114, in worker
    task = get()
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/queues.py", line 368, in get
    res = self._reader.recv_bytes()
Traceback (most recent call last):
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/connection.py", line 219, in recv_bytes
    buf = self._recv_bytes(maxlength)
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/process.py", line 315, in _bootstrap
    self.run()
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/connection.py", line 417, in _recv_bytes
    buf = self._recv(4)
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/process.py", line 108, in run
    self._target(*self._args, **self._kwargs)
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py", line 114, in worker
    task = get()
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/queues.py", line 367, in get
    with self._rlock:
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/connection.py", line 382, in _recv
    chunk = read(handle, remaining)
KeyboardInterrupt
  File "/opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/synchronize.py", line 101, in __enter__
    return self._semlock.__enter__()
KeyboardInterrupt
---------------------------------------------------------------------------
KeyboardInterrupt                         Traceback (most recent call last)
Cell In[6], line 5
      3 interpreter = jm.Interpreter(instance_data)
      4 instance = interpreter.eval_problem(problem)
----> 5 solution = jijsolver.solve(instance, time_limit_sec=1.0)

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/jijsolver/solver.py:388, in solve(ommx_instance, time_limit_sec, initial_ommx_state, seed, options)
    310 def solve(
    311     ommx_instance: Instance,
    312     *,
   (...)
    316     options: Optional[Union[JijSolverOption, SolverOptions]] = None,
    317 ) -> Solution:
    318     """
    319     Solve a problem using the weighted local search algorithm or the local ILP algorithm.
    320 
   (...)
    386     ```
    387     """
--> 388     entries = _sample_entries(
    389         ommx_instance=ommx_instance,
    390         time_limit_sec=time_limit_sec,
    391         initial_ommx_state=initial_ommx_state,
    392         num_samples=1,
    393         seed=seed,
    394         options=options,
    395     )
    397     assert len(entries) == 1
    399     return ommx_instance.evaluate(entries[0].state)

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/jijsolver/solver.py:218, in _sample_entries(ommx_instance, time_limit_sec, initial_ommx_state, num_samples, seed, options)
    216 else:
    217     with ProcessPool(nodes=num_processes) as pool:
--> 218         results = pool.map(
    219             _process,
    220             [
    221                 (
    222                     index,
    223                     converted_ommx_instance,
    224                     converted_ommx_instance_buf,
    225                     id_translator.to_data(),
    226                     initial_ommx_state,
    227                     num_samples,
    228                     seed,
    229                     options,
    230                 )
    231                 for index, options in enumerate(options.processes)
    232             ],
    233         )
    235 states: list[State] = []
    236 for result in results:

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/pathos/multiprocessing.py:154, in ProcessPool.map(self, f, *args, **kwds)
    152 if not OLD312a7:
    153     warnings.filterwarnings('ignore', category=DeprecationWarning)
--> 154 return _pool.map(star(f), zip(*args), **kwds)

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py:364, in Pool.map(self, func, iterable, chunksize)
    359 def map(self, func, iterable, chunksize=None):
    360     '''
    361     Apply `func` to each element in `iterable`, collecting the results
    362     in a list that is returned.
    363     '''
--> 364     return self._map_async(func, iterable, mapstar, chunksize).get()

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py:765, in ApplyResult.get(self, timeout)
    764 def get(self, timeout=None):
--> 765     self.wait(timeout)
    766     if not self.ready():
    767         raise TimeoutError

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/site-packages/multiprocess/pool.py:762, in ApplyResult.wait(self, timeout)
    761 def wait(self, timeout=None):
--> 762     self._event.wait(timeout)

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/threading.py:581, in Event.wait(self, timeout)
    579 signaled = self._flag
    580 if not signaled:
--> 581     signaled = self._cond.wait(timeout)
    582 return signaled

File /opt/hostedtoolcache/Python/3.9.22/x64/lib/python3.9/threading.py:312, in Condition.wait(self, timeout)
    310 try:    # restore state no matter what (e.g., KeyboardInterrupt)
    311     if timeout is None:
--> 312         waiter.acquire()
    313         gotit = True
    314     else:

KeyboardInterrupt: 

解の可視化#

得られた解を用いて、詰め込んだアイテムの価値と重量を表示してみましょう。

df = solution.decision_variables
obj = solution.objective
indices = np.ravel(df[df["value"]==1]["subscripts"].to_list())
sum_w = np.sum(inst_w[indices])

print("Values of chosen items: ", inst_v[indices])
print("Weights of chosen items: ", inst_w[indices])
print("Total value from objective: ", obj)
print("Total weight: ", sum_w)
Values of chosen items:  [27 24 18 21 17]
Weights of chosen items:  [26 22 16 19 17]
Total value from objective:  107.0
Total weight:  100

参考資料#

[1] Lucas, 2014, “Ising formulations of many NP problems”