Sampling Solutions

caution

The SampleSet described below will be deprecated at the end of April 2024. Please see What is SampleSet for a description of the new SampleSet.

I need a solution. What is a sampler?

Most solvers, including JijZept and OpenJij, use probabilistic heuristic algorithms to find a solution to the problem. They return candidates for the solution instead of the solution itself, and are therefore sometimes called "samplers". Their output can differ on each execution and may fail to obtain the optimal solution.

SampleSet object

The output format of each sampler is different, and the solver does not have information about the user-entered mathematical model description. Therefore, jijmodeling provides the SampleSet class to store the samples returned from samplers and to integrate with the mathematical model.

To demonstrate the features of SampleSet, let us consider the following mathematical model:

import jijmodeling as jm

N = jm.Placeholder('N')
a = jm.Placeholder('a', ndim=2)
x = jm.BinaryVar('x', shape=(N,N))
i = jm.Element('i', belong_to=(0, N))
j = jm.Element('j', belong_to=(0, N))

problem = jm.Problem('Simple Problem')
problem += jm.sum([i, j], a[i, j]*x[i, j])
problem += jm.Constraint('onehot-row', x[:, j].sum() == 1, forall=j)
problem += jm.Constraint('onehot-col', x[i, :].sum() == 1, forall=i)
problem

$$\begin{array}{cccc} \text{Problem:} & \text{Simple Problem} & & \\ & & \min \quad \displaystyle \sum_{i = 0}^{N - 1} \sum_{j = 0}^{N - 1} a_{i, j} \cdot x_{i, j} & \\ \text{{s.t.}} & & & \\ & \text{onehot-col} & \displaystyle \sum_{\ast_{1} = 0}^{N - 1} x_{i, \ast_{1}} = 1 & \forall i \in \left\{0,\ldots,N - 1\right\} \\ & \text{onehot-row} & \displaystyle \sum_{\ast_{0} = 0}^{N - 1} x_{\ast_{0}, j} = 1 & \forall j \in \left\{0,\ldots,N - 1\right\} \\ \text{{where}} & & & \\ & x & 2\text{-dim binary variable}\\ \end{array}$$

Next, we will solve this problem using OpenJij's simulated-annealing sampler. First, we will compile the problem into a QUBO (Quadratic Unconstrained Binary Optimization) format and then execute the sampling.

import openjij as oj
import jijmodeling as jm
import jijmodeling_transpiler as jmt

instance_data = {
 "N": 3,
 "a": [[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]]
}

# Translate into QUBO
compiled_model = jmt.core.compile_model(problem, instance_data)
pubo_builder = jmt.core.pubo.transpile_to_pubo(compiled_model)
qubo, constant = pubo_builder.get_qubo_dict()

# Setup simulated-annealing sampler
sampler = oj.SASampler()

# Execute Sampling, take 5 samples
response = sampler.sample_qubo(qubo, num_reads=5)

# Decode OpenJij format into jijmodeling format
sample_set = jmt.core.pubo.decode_from_openjij(
    response,       # result of sampler
    pubo_builder,   # How to restore the original problem from QUBO
    compiled_model  # Original problem
)

Here, the sample_set object stores 5 samples since the sampling is executed 5 times with num_reads=5.

Samples

Although the SampleSet stores the samples in a sparse manner for efficient memory usage, a dense representation is better for explanation. One can easily convert it into a dense representation using the to_dense() method:

for solution in sample_set.to_dense().record.solution["x"]:
    print(solution, "\n")

This code yields 5 samples of binary variables $x$ as follows:

[[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. 1.]]

[[0. 1. 0.]
 [0. 0. 1.]
 [1. 0. 0.]]

[[0. 0. 1.]
 [0. 1. 0.]
 [1. 0. 0.]]

[[0. 0. 1.]
 [0. 1. 0.]
 [1. 0. 0.]]

[[0. 0. 1.]
 [1. 0. 0.]
 [0. 1. 0.]]

tip

This may differ from the output of your environment.

There are several solutions (there are $3!$ cases) that satisfy all constraints $\sum_i x_{ij} = 1$, $\sum_j x_{ij} = 1$, and minimize the objective function $\sum_{ij} a_{ij} x_{ij}$ to the same value of $15$. The value of the objective function for each sample is stored in:

sample_set.evaluation.objective
# array([15., 15., 15., 15., 15.])

The violation of constraints, i.e. the absolute value of $\|\sum_i x_{ij} - 1\|$ and $\|\sum_j x_{ij} - 1\|$ for each sample, is stored as a dictionary with their names specified in the Problem.

sample_set.evaluation.constraint_violations
# {'onehot-col': array([0., 0., 0., 0., 0.]),
#  'onehot-row': array([0., 0., 0., 0., 0.])}

Sparse Representation

note

If you are familiar with sparse matrix representations, this corresponds to COO format. Please see scipy.sparse.coo_matrix for more details.

As noted above, SampleSet stores the samples in a sparse manner for efficient memory usage. Since there are only 3 non-zero elements in each sample, we only need to store where the non-zero elements are and what values they have. Sparse representations consist of an array of indices, values, and the shape of the array:

for indices, values, shape in sample_set.record.solution["x"]:
    print(f"{indices=}, {values=}, {shape=}")

indices=([0, 1, 2], [0, 1, 2]), values=[1.0, 1.0, 1.0], shape=(3, 3)
indices=([0, 1, 2], [1, 2, 0]), values=[1.0, 1.0, 1.0], shape=(3, 3)
indices=([0, 1, 2], [2, 1, 0]), values=[1.0, 1.0, 1.0], shape=(3, 3)
indices=([0, 1, 2], [2, 1, 0]), values=[1.0, 1.0, 1.0], shape=(3, 3)
indices=([0, 1, 2], [2, 0, 1]), values=[1.0, 1.0, 1.0], shape=(3, 3)

The shape is the shape of the matrix, and the values are the values of the non-zero elements. The indices are the indices of the non-zero elements. For example, the first sample has non-zero elements at $(0, 0)$, $(1, 1)$, and $(2, 2)$, so the indices are ([0, 1, 2], [0, 1, 2]).

Filtering Samples

jijmodeling.SampleSet has some convenient methods to extract the solutions with specific features.

• SampleSet.infeasible() method extracts the infeasible solutions, which break at least one constraint condition.

  sample_set.infeasible()

• SampleSet.feasible() method extracts the feasible solutions, which satisfy all the constraint conditions.

  sample_set.feasible()

• SampleSet.lowest() method extracts the solutions taking the minimum value of the objective functions among the feasible solutions. If there are multiple minimum values, lowest() method returns all these solutions.

  sample_set.lowest()