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}A=\{1, 2, 3, 4\}. It is easy to divide this set into equal sums: {1,4},{2,3}\{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 AA be the set to be partitioned and ai(i={0,1,,N1})a_i (i = \{0,1,\dots,N-1\}) be its elements. Here NN is the number of elements in this set. We consider dividing AA into two sets A0A_0 and A1A_1. Let xix_i be a variable whose iith element of AA is 0 when it is contained in the set A0A_0 and 1 when it is contained in A1A_1. Using this variable xix_i, the total value of the numbers in A0A_0 is written as iai(1xi)\sum_i a_i (1-x_i) and iaixi\sum_i a_i x_i for A1A_1. As we find a solution that satisfies the constraint that the sum of the numbers contained in A0A_0 and A1A_1 be equal, this can be expressed as:

iai(1xi)=iaixi\sum_i a_i (1-x_i)=\sum_i a_i x_i

The problem is to find xix_i that satisfies the constraint. By transforming this expression, we can write iai(2xi)=0\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=(i=0N1ai(2xi))2H=\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
Problem: number partitionmin(i=0N1ai(2xi1))2xi0{0,1}\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 NiN_i to NfN_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:

total value=(Ni+Nf)(NfNi+1)4\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 A1A_1 and A0A_0 in AA 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.