Number Partitioning#

Here we show how to solve the number partitioning problem using OpenJij, JijModeling, and JijModeling Transpiler. This problem is also mentioned in 2.1. Number Partitioning in Lucas, 2014, “Ising formulations of many NP problems”.

Overview of the Number Partitioning Problem#

Number partitioning is the problem of dividing a given set of numbers into two subsets such that the sum of the numbers is equal.

Example#

Let us have a set of numbers $A=\{1, 2, 3, 4\}$. It is easy to divide this set into equal sums: $\{1, 4\}, \{2, 3\}$ and the sum of each subset is 5. Thus, when the size of the set is small, the answer is relatively easy to obtain. When the problem is large, however, it is not immediately solvable. Here, we solve this problem using annealing.

from jijmodeling.transpiler.pyqubo import to_pyqubo
import openjij as oj
import jijmodeling as jm
import numpy as np


Mathematical Model#

First, let us model the Hamiltonian of this problem. Let $A$ be the set to be partitioned and $a_i (i = \{0,1,\dots,N-1\})$ be its elements. Here $N$ is the number of elements in this set. We consider dividing $A$ into two sets $A_0$ and $A_1$. Let $x_i$ be a variable whose $i$th element of $A$ is 0 when it is contained in the set $A_0$ and 1 when it is contained in $A_1$. Using this variable $x_i$, the total value of the numbers in $A_0$ is written as $\sum_i a_i (1-x_i)$ and $\sum_i a_i x_i$ for $A_1$. As we find a solution that satisfies the constraint that the sum of the numbers contained in $A_0$ and $A_1$ be equal, this can be expressed as:

$\sum_i a_i (1-x_i)=\sum_i a_i x_i$

The problem is to find $x_i$ that satisfies the constraint. By transforming this expression, we can write $\sum_i a_i (2-x_i)=0$, and by using the Penalty method and squaring this constraint, the Hamiltonian for the number-splitting problem is:

$H=\left( \sum_{i=0}^{N-1} a_i (2-x_i)\right)^2$

Modeling by JijModeling#

Next, we show how to implement the above equation using JijModeling. We first define variables for the mathematical model described above.

problem = jm.Problem("number partition")
a = jm.Placeholder("a",dim = 1)
N = a.shape[0].set_latex("N")
x = jm.Binary("x",shape=(N,))
i = jm.Element("i",(0,N))
s_i = 2*x[i] - 1
problem += (jm.Sum(i,a[i] * s_i))**2
problem

\begin{alignat*}{4}\text{Problem} & \text{: number partition} \\\min & \quad \left( \sum_{ i = 0 }^{ N - 1 } a_{i} \cdot \left( 2 \cdot x_{i} - 1 \right) \right) ^ { 2 } \\& x_{i_{0}} \in \{0, 1\}\end{alignat*}

Instance#

As an example, let us solve an easy problem; let us consider the problem of dividing numbers from 1 to 40. When dividing consecutive numbers from $N_i$ to $N_f$ and keeping the total number of consecutive numbers even, there are several patterns of division. However, the total value of the divided set is:

$\mathrm{total\ value} = \frac{(N_{i} + N_{f})(N_{f} - N_{i} + 1)}{4}$

In this case, the total value is expected to be 410. Let us confirm this.

N = 40
instance_data = {"a":np.arange(1,N+1)}


Conversion to PyQUBO by JijModeling Transpiler#

JijModeling has executed all the implementations so far. By converting this to PyQUBO, it is possible to perform combinatorial optimization calculations using OpenJij and other solvers.

model,cache = to_pyqubo(problem,instance_data,{})
Q,offset = model.compile().to_qubo()


The PyQUBO model is created by to_pyqubo with the problem created by JijModeling and the instance_data we set to a value as arguments. Next, we compile it into a QUBO model that can be computed by OpenJij or other solver.

Optimization by OpenJij#

This time, we will use OpenJij’s simulated annealing to solve the optimization problem. We set the SASampler and input the QUBO into that sampler to get the result of the calculation.

sampler = oj.SASampler(num_reads=100)
response = sampler.sample_qubo(Q=Q)


Decoding and Displaying the Solution#

Let us take a look at the results obtained. We decode the returned calculation results to facilitate analysis.

result = cache.decode(response)


We sum the indices classified as $A_1$ and $A_0$ in $A$ here.

class_1_index = result.record.solution['x'][0][0][0]
class_0_index = [i for i in range(0,N) if i not in class_1_index]

class_1 = instance_data['a'][class_1_index]
class_0 = instance_data['a'][class_0_index]

print(f"class 1 : {class_1} , total value = {np.sum(class_1)}")
print(f"class 0 : {class_0} , total value = {np.sum(class_0)}")

class 1 : [ 2  4  7  8 11 12 13 15 20 22 23 27 29 31 35 36 37 38 40] , total value = 410
class 0 : [ 1  3  5  6  9 10 14 16 17 18 19 21 24 25 26 28 30 32 33 34 39] , total value = 410


As expected, both total values are 410. Above, we dealt with a problem whose solution is known because it is a consecutive number. We recommend you try more complex problems, such as generating numbers randomly.