Combinatorial classes¶

Four Quadratic Unconstrained Binary Optimisation (QUBO) example problems are listed here: the Travelling Salesman Problem, the Quadratic Assignment Problem, the Minimal Vertex Cover and the Maximum Independent Set.

Travelling Salesman Problem¶

The Travelling Salesman Problem (sometimes referred to as the Travelling Salesperson Problem), commonly abbreviated as TSP, is a NP-hard problem in combinatorial optimisation.

Briefly, the problem revolves around finding the shortest possible route for a salesman to visit some cities before returning to the origin. TSP is usually formulated as a graph problem with nodes specifying the cities and edges denoting the distances between each city.

The idea behind TSP can be mapped to similar-type problems. For instance, what is the optimal route for the salesman to take in order to minimise something.

In this module, the TSP class follows Hadfield’s 2017 paper.

class qiboopt.combinatorial.combinatorial.TSP(distance_matrix, two_to_one=None, backend=None)[source]¶

Class representing the Travelling Salesman Problem (TSP). The implementation is based on the work by Hadfield.

Parameters:

distance_matrix – a numpy matrix encoding the distance matrix.
two_to_one – Mapping from two coordinates to one coordinate
backend – Backend to use for calculations. If not given the global backend will be used.

Example

from qiboopt.combinatorial.combinatorial import TSP
import numpy as np
from collections import defaultdict
from qibo import Circuit, gates
from qibo.models import QAOA
from qibo.result import CircuitResult


def convert_to_standard_Cauchy(config):
    m = int(np.sqrt(len(config)))
    cauchy = [-1] * m  # Cauchy's notation for permutation, e.g. (1,2,0) or (2,0,1)
    for i in range(m):
        for j in range(m):
            if config[m * i + j] == '1':
                cauchy[j] = i  # citi i is in slot j
    for i in range(m):
        if cauchy[i] == 0:
            cauchy = cauchy[i:] + cauchy[:i]
            return tuple(cauchy)
            # now, the cauchy notation for permutation begins with 0


def evaluate_dist(cauchy):
    '''
    Given a permutation of 0 to n-1, we compute the distance of the tour

    '''
    m = len(cauchy)
    return sum(distance_matrix[cauchy[i]][cauchy[(i+1)%m]] for i in range(m))


def qaoa_function_of_layer(layer, distance_matrix):
    '''
    This is a function to study the impact of the number of layers on QAOA,
    it takes in the number of layers and compute the distance of the mode
    of the histogram obtained from QAOA

    '''
    small_tsp = TSP(distance_matrix)
    obj_hamil, mixer = small_tsp.hamiltonians()
    qaoa = QAOA(obj_hamil, mixer=mixer)
    best_energy, final_parameters, extra = qaoa.minimize(initial_p=[0.1] * layer,
                                         initial_state=initial_state, method='BFGS')
    qaoa.set_parameters(final_parameters)
    quantum_state = qaoa.execute(initial_state)
    circuit = Circuit(9)
    circuit.add(gates.M(*range(9)))
    result = CircuitResult(quantum_state, circuit.measurements,
            small_tsp.backend, nshots=1000)
    freq_counter = result.frequencies()
    # let's combine freq_counter here, first convert each key and sum up the frequency
    cauchy_dict = defaultdict(int)
    for freq_key in freq_counter:
        standard_cauchy_key = convert_to_standard_Cauchy(freq_key)
        cauchy_dict[standard_cauchy_key] += freq_counter[freq_key]
    max_key = max(cauchy_dict, key=cauchy_dict.get)
    return evaluate_dist(max_key)

np.random.seed(42)
num_cities = 3
distance_matrix = np.array([[0, 0.9, 0.8], [0.4, 0, 0.1],[0, 0.7, 0]])
distance_matrix = distance_matrix.round(1)
small_tsp = TSP(distance_matrix)
initial_parameters = np.random.uniform(0, 1, 2)
initial_state = small_tsp.prepare_initial_state([i for i in range(num_cities)])
qaoa_function_of_layer(2, distance_matrix)

Reference:: 1. S. Hadfield, Z. Wang, B. O’Gorman, E. G. Rieffel, D. Venturelli, R. Biswas, From the Quantum Approximate Optimization Algorithm to a Quantum Alternating Operator Ansatz. (arxiv:1709.03489)

hamiltonians()[source]¶

Constructs the phaser and mixer Hamiltonians for the TSP.

Returns:: A tuple containing the phaser and mixer Hamiltonians.
Return type:: tuple

prepare_initial_state(ordering)[source]¶

Prepares a valid initial state for TSP QAOA based on the given city ordering.

Parameters:: ordering (list) – A list representing the permutation of cities.
Returns:: The quantum state representing the initial solution.
Return type:: array

penalty_method(penalty)[source]¶

Constructs the TSP QUBO object using a penalty method for feasibility.

The TSP is formulated as:

\[\min \sum_{u,v,j} d_{u,v} \, x_{u,j} \, x_{v,j+1}\]

Subject to constraints:

\[\sum_j x_{v,j} = 1 \quad \forall v \ \sum_v x_{v,j} = 1 \quad \forall j\]

The penalty method converts this to an unconstrained QUBO:

\[f(x) = \text{objective}(x) + \lambda \left[ \sum_v \left(\sum_j x_{v,j} - 1\right)^2 + \sum_j \left(\sum_v x_{v,j} - 1\right)^2 \right]\]

Parameters:: penalty (float) – The penalty parameter for constraint violations. It should be large enough to enforce constraints but not so large as to cause numerical issues.
Returns:: A QUBO object for the TSP with penalties applied.
Return type:: QUBO
Raises:: ValueError – If penalty is negative.

Quadratic Assignment Problem¶

The Quadratic Assignment Problem, commonly abbreviated as QAP, is a NP-hard problem in combinatorial optimisation, first introduced by Koopmans and Beckmann.

Briefly, the problem concerns assigning a set of facilities to a set of locations in a way that minimises the total assignment cost. The assignment cost is typically determined by the flow between facilities and the distance between their assigned locations. QAP is usually formulated using two matrices: one representing flows between facilities and another representing distances between locations.

To map QAO to a QUBO, we define binary variables \(x \in \{0, 1\}\), where 1 denotes if the \(i\)-th facility is assigned to the \(j\)-th location and 0 otherwise.

The objective function is

\[\min \sum_{i,k=1}^{n} \sum_{j,l=1}^{n} f_{ik} \, d_{jl} \, x_{ij} \, x_{kl}\]

subject to the assignment constraints

\[\sum_{k=1}^{n} x_{ik} = 1 \quad \forall i \in \{1,\dots,n\}\]

and

\[\sum_{i=1}^{n} x_{ik} = 1 \quad \forall k \in \{1,\dots,n\}\]

where \(f_{ik}\) are the elements of the flow matrix \(F \in \mathbb{R}^{n \times n}\), representing the interaction between facilities \(i\) and \(k\), and \(d_{jl}\) are the elements of the distance matrix \(D \in \mathbb{R}^{n \times n}\), representing the distance between locations \(j\) and \(l\).

class qiboopt.combinatorial.combinatorial.QAP(flow_matrix, distance_matrix, two_to_one=None)[source]¶

Class for representing the Quadratic Assignment Problem (QAP).

The QAP assigns a set of facilities to a set of locations in a way that minimises the total assignment cost.

Parameters:

flow_matrix – An \((n,n)\) numpy array describing the flows between the facilities.
distance_matrix – A \((n,n)\) numpy array describing the distance between the locations
two_to_one (optional) – Mapping from 2-index to a single index.

Example

from qiboopt.combinatorial.combinatorial import QAP

flow_matrix = np.array([[1, 2],[3, 4]])
distance_matrix = np.array([[4, 3],[2, 1]])
qap = QAP(flow_matrix, distance_matrix)
penalty = 10
qp = qap.penalty_method(penalty)

penalty_method(penalty)[source]¶: construct the QUBO instance when the penalty method is being used.

suggest_penalty()[source]¶: This function suggest a penalty coefficient for a QAP instance. Currently, we use the dimension multiplied by the maximum of flow and maximum of the distance matrix. This is a heuristic.

Minimum Weighted Vertex Cover¶

The Minimum Weighted Vertex Cover Problem, commonly abbreviated as MWVC, is a NP-hard problem in combinatorial optimisation.

Briefly, the problem involves selecting a subset of vertices in a graph such that every edge in the graph is incident to at least one selected vertex. Each vertex is assigned a weight and the objective is to minimise the total weight of the selected vertices (rather than their count). MWVC is usually formulated as a graph problem where nodes represent entities and edges represent relationships between them.

Given a graph \(G = (V, E)\), the goal of MWVC is to select a subset of vertices in such a way that every edge is covered while minimising total vertex weight.

Define binary variables \(x_i \in \{0, 1\}\) for each vertex \(i \in V\), where \(x_i = 1\) indicates that vertex \(i\) is selected and 0 otherwise.

To map MWVC to a QUBO, constraints are enforced via penalty terms:

\[\min \sum_{i \in V} w_i x_i + P \sum_{(i,j) \in E} (1 - x_i - x_j)^2\]

where \(w_i\) is the weight of vertex \(i\) and \(P\) is a sufficiently large penalty coefficient enforcing coverage.

class qiboopt.combinatorial.combinatorial.MWVC(graph)[source]¶

Class for representing the Minimum Weighted Vertex Cover (MWVC) problem.

The MWVC problem selects a subset of vertices in a graph such that every edge in the graph is incident to at least one selected vertex.

Parameters:: graph – a networkx graph where the weight is encoded with keyword “weight”.

Example

import networkx as nx
from qiboopt.combinatorial.combinatorial import MWVC

g = nx.Graph()
g.add_nodes_from([(0, {"weight":2}), (1, {"weight":3}), (2, {"weight":4})])
g.add_edges_from([(0,1), (1,2), (2,1)])
mwvc = MWVC(g)
penalty = 10
qp = mwvc.penalty_method(penalty)

penalty_method(penalty=2)[source]¶: this function implement the penalty method for minimum vertex cover. note that for unweighted version of minimum vertex cover, setting 2 as the penalty coefficient works. :param penalty: the penalty is a real positive number.

suggest_penalty(epsilon)[source]¶: This function suggest a penalty coefficient fo the MWVC instance. We set it to be the largest weight plus a small perturbation. :param epsilon: a positive number that should be a small positive number, serving as a perturbation

Maximum Independent Set¶

The MIS problem involves selecting the largest subset of non-adjacent vertices in a graph.

class qiboopt.combinatorial.combinatorial.MIS(g)[source]¶

Class for representing the Maximal Independent Set (MIS) problem.

The MIS problem involves selecting the largest subset of non-adjacent vertices in a graph.

Parameters:: g (networkx.Graph) – A graph object representing the problem.

Example

import networkx as nx
from qiboopt.combinatorial.combinatorial import MIS

g = nx.Graph()
g.add_edges_from([(0, 1), (1, 2), (2, 0)])
mis = MIS(g)
penalty = 10
qp = mis.penalty_method(penalty)

penalty_method(penalty)[source]¶

Constructs the QUBO Hamiltonian for the MIS problem using a penalty method.

Parameters:: penalty (float) – The penalty parameter for constraint violations.
Returns:: A QUBO object for the Maximal Independent Set (MIS) problem.
Return type:: QUBO (qiboopt.opt_class.opt_class.QUBO)