qat-opt: representation and encoding of combinatorial problems, QAOA

This section describes how to define and solve combinatorial problems using myQLM tools.

Declaring and manipulating combinatorial problems

The most generic class used to describe combinatorial problems is the qat.opt.CombinatorialProblem class. It provides a simple interface to declare boolean variables and clauses, and ways to turn abstract problems into Job or objects that can be fed to solvers.

In the (quite common case) where the problem’s clauses only involve at most two variables at a time (QUBO), we also provide a couple of optimized classes to describe the problem directly in terms of Ising coupling matrix or QUBO \(Q\) matrix (see below). Back-and-forth translations are avaible between the three problem classes using the to_XXX methods.

Generic combinatorial optimization

class qat.opt.CombinatorialProblem(name=None, maximization=False, **kwargs)

Basic interface to describe a combinatorial optimization problem.

The problem declaration is done via methods new_var() (or new_vars() to declare arrays) and add_clause().

from qat.opt import CombinatorialProblem

problem = CombinatorialProblem("MyProblem")
# Declare two fresh variables
var1, var2 = problem.new_vars(2)
# Add a new clause consisting of the logical AND of the two variables
problem.add_clause(var1 & var2)
# Add a new clause consisting of the XOR of the two variables
problem.add_clause(var1 ^ var2)

print(problem)
MyProblem:
 2 variables, 2 clauses

It is possible to add weights to the clauses:

from qat.opt import CombinatorialProblem

problem = CombinatorialProblem()
var1, var2 = problem.new_vars(2)
problem.add_clause(var1 & var2, weight=0.5)

print(problem)
Problem:
 2 variables, 1 clauses

A diagonal Hamiltonian encoding the cost function of the problem can be extracted using the get_observable() method.

from qat.opt import CombinatorialProblem

problem = CombinatorialProblem()
var1, var2 = problem.new_vars(2)
problem.add_clause(var1 & var2, weight=0.5)

obs = problem.get_observable()

print(obs)
0.125 * I^2 +
-0.125 * (Z|[0]) +
0.125 * (ZZ|[0, 1]) +
-0.125 * (Z|[1])

Finally, this class inherits from the CircuitGenerator class, which provides a method to directly generate variational Ansätze to try and minimize the energy of the cost Hamiltonian.

Parameters
  • name (optional, str) – a name to display when the problem is printed

  • maximization (optional, bool) – Used to specify that the problem is a maximization problem (i.e its cost function is the sum of its clauses). In practice, it will simply flip the sign of the generated cost Hamiltonian. Default to false.

add_clause(clause, weight=None)

Adds a new clause to the problem.

Parameters
  • clause (Clause) – a clause object

  • weight (optional, float) – optionally a weight (default to 1)

Returns

the problem itself

get_observable()

Generates a cost Hamiltonian for the problem.

The cost Hamiltonian is diagonal and associate to each bitstring \(|s\rangle\) an energy \(\sum_\alpha w_\alpha C_\alpha(s)\) where \(C_\alpha\) are the clauses of the problem, seen as \(\{0, 1\}\) valued functions and \(w_\alpha\) their corresponding weights.

This encoding is done recursively and is described in the documentation of the Clause class.

If the problem is specified as a maximization problem, the sign of the cost Hamiltonian if flipped. This means that the “best” solution is always encoded in the ground state of the returned Hamiltonian.

Returns

a cost Hamiltonian

Return type

Observable

new_var()

Returns a fresh variable

Returns

a fresh variable

Return type

Var

new_vars(nbvars)

Returns a list of fresh variables.

Parameters

nbvars (int) – the number of fresh variables to declare

Returns

a list of fresh variables of length nbvars

Return type

list

to_ising()

Translates the problem into an Ising problem. Might raise an exception if the problem is not quadratic.

Returns

an Ising object

Return type

Ising

to_qubo()

Translates the problem into a QUBO problem. Might raise an exception if the problem is not quadratic.

Returns

a QUBO object

Return type

QUBO

Clauses (Clause) are declared by combining variables (Var). The cost Hamiltonian extraction is handled by the (Clause) class. The final cost Hamiltonian consists of the weighted sum of the cost Hamiltonian of its clauses.

class qat.opt.boolexpr.Clause(operator, arity, *children)

Class representing a boolean clause (boolean formula).

Clauses are trees whose leaves are Var objects and internal nodes are labelled by boolean operators (&, |, ^, ~).

Clauses are meant to be built either using the .and_clause, .or_clause, .xor_clause, neg_clause static methods, or via boolean operators overloading:

Parameters
  • operator (str) – the operator (in [”&”, “|”, “^”, “~”])

  • arity (int) – the arity of the operator

  • *children (Clause) – the subclauses

get_observable()

Returns a boolean valued diagonal cost observable matching the evaluation of the clause.

The cost observable is built by induction:

\(H(A \wedge B) = H(A) H(B)\)

\(H(A \vee B) = H(A) + H(B) - H(A) H(B)\)

\(H(A \oplus B) = H(A) + H(B) - 2 H(A) H(B)\)

\(H(\neg A) = 1 - H(A)\)

\(H(v) = (1 - \sigma_z^v)/2\)

Consequently: the resulting Hamiltonian is {0, 1} valued and has as 1-eigenstates the states that satisfy the clause.

The variable class is quite simple and overloads logical operators to closely interact with the Clause class:

class qat.opt.boolexpr.Var(index)

Simple class for boolean variables

index

the variable index (unique identifier)

Type

int

Parameters

index (int) – the variable index

get_observable()

Returns a boolean valued diagonal cost observable matching the evaluation of the variable: \((1 - \sigma_z^i) / 2\) where \(i\) is the index of the variable.

The observable will act upon self.index + 1 qubits.

get_variables()

Returns a singleton with the index of the variable

Quadratic problems: QUBO and Ising

class qat.opt.QUBO(Q, offset_q=0, **kwargs)

Class representing Quadratic Unconstrained Binary Optimization (QUBO) problems.

The class allows for the representation of a problem as QUBO - by providing with a symmetric \(Q\) matrix and a QUBO offset energy \(E_Q\), both of which coming from the respective Hamiltonian encoding,

\[H = - x^T Q x - E_Q\]

where \(x\) is the vector of binary values \(\{0,1\}\) we look for, such that \(H\) is minimum.

The class can also translate from a QUBO problem to an Ising problem via to_ising() by returning an Ising object.

QUBO problems can be translated to a CombinatorialProblem object via to_combinatorial_problem().

This class also inherits from the CircuitGenerator class which allows to construct QAOA-Ansätze.

Parameters
  • Q (2D numpy array) – a symmetric array representing the \(Q\) matrix from the Hamiltonian of the problem

  • offset_q (optional, double) – the value of the QUBO offset energy in the Hamiltonian of the problem

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA) of the respective QUBO problem. SQA is available in QLM, but the temperature parameters temp_max and temp_min could also be used for simulated annealing (SA), available in myQLM.

The method should be called from one of the child problem classes (e.g. GraphColouring, VertexCover, etc.) and not directly from the parent QUBO class as this will raise an exception.

get_observable()

Returns an Observable for the problem containing an Ising tuple.

Returns

an Ising Observable representing the problem

Return type

Observable

get_q_and_offset()

This method returns the \(Q\) matrix and QUBO energy offset, which define the QUBO object.

Returns

2-element tuple containing

  • Q (2D numpy array) - a symmetric array representing the \(Q\) matrix from the Hamiltonian of the problem

  • offset_q (double) - the value of the QUBO offset energy in the Hamiltonian of the problem

to_combinatorial_problem()

Translates the QUBO problem into a combinatorial problem.

Returns

a combinatorial problem instance

Return type

CombinatorialProblem

to_ising()

Translates the QUBO problem into an Ising problem over spins.

Returns

an Ising instance

Return type

Ising

to_job(gamma_t=None, tmax=1.0, **kwargs)

Returns a Job for the problem - ready to run on a QPU for simulated annealing (SA) or simulated quantum anealing (SQA). The second one comes in the QLM.

Parameters:
  • tmax (float, optional) - time duration of the annealing. Default is 1.

  • gamma_t (ArithExpression, SQA only) - a function specifying the time dependence of Gamma. It should be produced using the variable ‘t’ created by the class Variable.

Returns

a ready to run Job for the problem

Return type

Job

Note

For the supported NP problems, some well performing max and min Gamma for a linearly decreasing gamma_t have been found and could be accessed via the method get_best_parameters of the respective problem class.

class qat.opt.Ising(J, h, offset_i=0, **kwargs)

Class representing Ising problems.

The class allows for the representation of a problem in the Ising framework - by providing with a symmetric coupling matrix \(J\) with zero diagonal elements, magnetic field \(h\) and an Ising offset energy \(E_I\), all of which coming from the respective Hamiltonian encoding,

\[H = - s^T J s - h^T s - E_I\]

where \(s\) is the spin vector we look for with values \(\{-1,1\}\), such that \(H\) is minimum.

The class can also translate from an Ising problem to a QUBO problem through to_qubo() by returning a QUBO object.

Ising problems can be translated to a CombinatorialProblem object via to_combinatorial_problem().

Similarly to QUBO and CombinatorialProblem, this class inherits from the CircuitGenerator class, thus is able to generate QAOA-Ansätze.

Parameters
  • J (2D numpy array) – a symmetric array with zero diagonal elements for the coupling between each two spins - it represents the \(J\) matrix from the Hamiltonian of the problem

  • h (1D numpy array) – an array with the magnetic field acting on each of the spins, coming from the Hamiltonian of the problem

  • offset_i (optional, double) – the value of the Ising offset energy in the respective Hamiltonian

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA) of the respective Ising problem. SQA is available in QLM, but the temperature parameters temp_max and temp_min could also be used for simulated annealing (SA), available in myQLM.

The method should be called from one of the child problem classes (e.g. MaxCut, NumberPartitioning, etc.) and not directly from the parent Ising class as this will raise an exception.

get_j_h_and_offset()

This method returns the \(J\) coupling matrix, the magnetic field \(h\) and Ising energy offset, which define the Ising object.

Returns

3-element tuple containing

  • J (2D numpy array) - a symmetric array with zero diagonal elements for the coupling between each two spins - it represents the \(J\) matrix from the Hamiltonian of the problem

  • h (1D numpy array) - an array with the magnetic field acting on each of the spins, coming from the Hamiltonian of the problem

  • offset_i (double) - the value of the Ising offset energy in the respective Hamiltonian

get_observable()

Returns an Observable for the problem containing an Ising tuple.

Returns

an Ising Observable representing the problem

Return type

Observable

to_combinatorial_problem()

Translates the Ising problem into a CombinatorialProblem.

Returns

a combinatorial problem instance

Return type

CombinatorialProblem

to_job(gamma_t=None, tmax=1.0, **kwargs)

Returns a Job for the problem - ready to run on a QPU for simulated annealing (SA) or simulated quantum anealing (SQA). The second one comes in the QLM.

Parameters:
  • tmax (float, optional) - time duration of the annealing. Default is 1.

  • gamma_t (ArithExpression, SQA only) - a function specifying the time dependence of Gamma. It should be produced using the variable ‘t’ created by the class Variable.

Returns

a ready to run Job for the problem

Return type

Job

Note

For the supported NP problems, some well performing max and min Gamma for a linearly decreasing gamma_t have been found and could be accessed via the method get_best_parameters of the respective problem class.

to_qubo()

Translates the Ising problem into a QUBO problem.

Returns

a QUBO object

Return type

QUBO

Generating the QAOA Ansatz

The Quantum Approximate Optimization Algorithms is a heuristics to design variational Ansätze for combinatorial optimization. It is inspired from the digitalization of an analog evolution using a linear ramp, starting from a simple initial Hamiltonian \(H_0 = - \sum_i \sigma_x^i\) to a diagonal Hamiltonian whose ground state encodes the solution to our problem. This digitalization leads to a layered parametrized quantum circuit consisting of entangling layers sperated by collective \(R_X\) rotations.

QAOA Ansätze are usually parametrized by a depth parameter specifying the number of alternating layers.

It is possible to directly generate ready to run QAOA jobs (containing an Ansatz and the target Hamiltonian) from an instance of qat.opt.CombinatorialProblem/qat.opt.QUBO/qat.opt.Ising instance via the following interface:

class qat.opt.circuit_generator.CircuitGenerator

Class for circuit generation interface from a diagonal observable.

qaoa_ansatz(depth, cnots=True, strategy='coloring', to_circ_args=None, **kwargs)

Generates a QAOA Ansatz using the cost observable returned by the abstract method get_observable.

Warning

When setting the cnots option to False, the circuit might make use of generalized many-qubits Z rotations. In that case, you might want to instantiate your variational plugins using a gate set that contains definition of these gates. If not, some matrices in the circuit structure will be missing and some QPUs may not be able to handle the circuit.

The following piece of code should allow you to link the correct gate set to a variational plugin:

from qat.plugins import ScipyMinimizePlugin
from qat.vsolve.ansatz import get_qaoa_gate_set

# This plugin will no be able to bind variables inside a
# job generated with cnot set to False!
my_plugin = ScipyMinimizePlugin()

# This plugin can now be used with job generated with the
# cnots option sets to False!
my_plugin = ScipyMinimizePlugin(gate_set=get_qaoa_gate_set())
Parameters
  • depth (int) – the depth of the Ansatz

  • strategy (str) – the strategy to adopt to generate the circuit. Possible strategies are “default” or “coloring”. The “coloring” strategy uses a greedy coloring heuristics to try to optimize the overall depth of the Ansatz. Default is “default” which synthesize the circuit without optimizing the term ordering.

  • cnots (optional, bool) – If set to True the Ansatz will only use CNOT gates. If set to False, some abstract gates will be used to generate collective pauli rotations, resulting in a lower gate count. Defaults to True.

  • **kwargs – optional arguments that will be transfered to the job’s constructor (e.g nbshots, etc).

Returns

a qlm job, ready to run

Return type

Job

The qat.vsolve.ansatz.AnsatzFactory provides a recipe to produce such a variational circuits from a target Hamiltonian.

class qat.vsolve.ansatz.AnsatzFactory

This class regroups the implementation all the different Ansätze available in the QLM.

static qaoa_circuit(observable, depth, strategy='default', cnots=True, **to_circ_args)

Generates a QAOA Ansatz from an observable and an Ansatz depth

Example

from qat.core import Observable, Term
from qat.vsolve.ansatz import AnsatzFactory

line_obs = Observable(10)
for i in range(9):
    line_obs.add_term(Term(1., "ZZ", [i, i+1]))
ansatz_with_cnots = AnsatzFactory.qaoa_circuit(line_obs, 3)
print("The Ansatz has {} gates".format(len(ansatz_with_cnots.ops)))
ansatz_with_rzz = AnsatzFactory.qaoa_circuit(line_obs, 3, cnots=False)
print("The Ansatz has {} gates".format(len(ansatz_with_rzz.ops)))
The Ansatz has 121 gates
The Ansatz has 67 gates

The synthesis strategy may influence the depth of the circuit:

def depth(circuit):
    ''' Computes the depth of a circuit '''
    slices = [set()]
    for op in circuit:
        qbits = op.qbits
        insert_in = None
        for index, slic in enumerate(reversed(slices)):
            if all(qb not in slic for qb in qbits):
                continue
            insert_in = index
            break
        if insert_in is None:
            for qb in qbits:
                slices[0].add(qb)
        elif insert_in == 0:
            slices.append(set(qbits))
        else:
            for qb in qbits:
                slices[len(slices) - insert_in].add(qb)
    return len(slices)

from qat.core import Observable, Term
from qat.vsolve.ansatz import AnsatzFactory

line_obs = Observable(10)
for i in range(9):
    line_obs.add_term(Term(1., "ZZ", [i, i+1]))
ansatz_default = AnsatzFactory.qaoa_circuit(line_obs, 3, strategy="default")
print("The Ansatz has depth {}".format(depth(ansatz_default)))
ansatz_coloring = AnsatzFactory.qaoa_circuit(line_obs, 3, strategy="coloring")
print("The Ansatz has depth {}".format(depth(ansatz_coloring)))
The Ansatz has depth 43
The Ansatz has depth 22

When considering QAOA instances with large Clauses (i.e clauses with more than 2 variables), the “gray_synth” strategy can often remove lots of CNOTS:

def cnot_count(circuit):
    ''' count cnots in a circuit '''
    return sum(1 if name == "CNOT" else 0
               for name, args, qbits in circuit.iterate_simple())

from qat.core import Observable, Term
from qat.vsolve.ansatz import AnsatzFactory

n = 5
line_obs = Observable(n)
for i in range(n - 2):
    line_obs.add_term(Term(1., "ZZZ", [i, i+1, i+2]))

ansatz_default = AnsatzFactory.qaoa_circuit(line_obs, 3, strategy="default")
ansatz_gray_synth = AnsatzFactory.qaoa_circuit(line_obs, 3, strategy="gray_synth")
print("Cnot count in default:", cnot_count(ansatz_default))
print("Cnot count in gray synth:", cnot_count(ansatz_gray_synth))
Cnot count in default: 36
Cnot count in gray synth: 24

Synthesis strategies:

  • default: uses the default term ordering provided by the input observable

  • coloring: orders terms using a graph coloring technique in order to reduce circuit depth

  • gray_synth: uses Amy et al GraySynth algorithm to synthesize the entangling layer. This might help in reducing the overall CNOT count.

Parameters
  • observable (Observable) – some diagonal observable

  • depth (int) – the depth of the Ansatz

  • strategy (str) – the strategy to adopt to generate the circuit.

  • cnots (optional, bool) – if set to True, the generator will onlt use CNOT gates as entangling gates. Default to True. This argument is ignored for some strategies.

  • **to_circ_args – arguments passed to the to_circ method

Generating QA jobs

It is possible to turn qat.opt.CombinatorialProblem/qat.opt.QUBO/qat.opt.Ising into Quantum Annealing jobs thanks to the ScheduleGenerator interface:

class qat.opt.schedule_generator.ScheduleGenerator

Interface for problems that can be solved by a quantum annealing.

It requires the implementation of method that returns an Hamiltonian formulation of the target cost function.

annealing_job(tmax=None, mixing=None, **kwargs)

Generates a Quantum Annealing job performing a linear interpolation between an initial mixing Hamiltonian and the problem’s Hamiltonian.

Mixing Hamiltonians can be generated through the qat.opt.MixingFactory:

class qat.opt.MixingFactory

A factory to define problem independent mixing Hamiltonians (i.e initial Hamiltonians in quantum annealing schedules).

All methods generate object implementing the InitialStateBuilder interface and returns a pair (such an object, mixing Hamiltonian).

static bit_flip(nbits, restrict_to=None)

Builds a bit flip mixing Hamiltonian:

\[H_0 = - \sum \sigma_x\]

Its ground state is the product state:

\[|\psi_0\rangle = |+\rangle ^ {\otimes n}\]

and is returned in string format.

This is the standard Hamiltonian for most applications.

Parameters
  • nbits (int) – the number of qubits in the system

  • restrict_to (optional, list<int>) – a possible list of qubits to restrict the mixing to. If set, the sum over \(\sigma_x\) is restricted to these qubits. The initial state of other qubits will be set to \(|0\rangle\).

Returns

a pair (structure preparing the initial state, mixing Hamiltonian).

Return type

(SimpleInitialState, Observable)

static bit_move(nbits, hamming_weight, exact=True, tmax_psi_0=50.0)

Builds a bit move mixing Hamiltonian that preserves Hamming weight:

\[H_0 = - \sum_{i, j} S_{i, j}\]

where \(S_{i,j}\) is the following two qubits operator located on qubits \(i, j\):

\[S_{i, j} = \frac{\sigma_y\otimes\sigma_y + \sigma_x\otimes\sigma_x}{2}\]

The main purpose of this mixing is to explore subset of bitstrings of constant Hamming weight. In that regard, it is important to start from an appropriate initial state that corresponds to the ground state of \(H_0\) restricted to the classical states of fixed Hamming weight \(k\).

This ground state is the equi-superposition of classical states of Hamming weight \(k\):

\[|\psi_0\rangle = {n \choose k}^{-\frac{1}{2}} \sum_{x} |x\rangle\]

where the sum ranges over all bit-strings \(x\) of Hamming weight \(k\).

When using a simulator, one can prepare this state using a numpy array (which is fine for small number of qubits).

In real use cases, it is possible to approximate this state by performing a problem independent annealing using:

\[ \begin{align}\begin{aligned}H'_0 = - \sum_i \sigma_x^i\\|\psi'_0\rangle = |+\rangle ^{\otimes n}\end{aligned}\end{align} \]

and

\[H'_1 = \left(\sum_i \frac{1 - \sigma_z^i}{2} - k\right)^2\]

The ground state of \(H'_1\) is exactly the state \(|\psi_0\rangle\). Moreover, it is argued in [] that the min-gap of this evolution is located close to the end of the evolution and that its width is a \(O(\frac{1}{T})\). In particular, this entails that one can pick a polynomially large \(T\) and be sure to build a large overlap between the state of our system and the perfect \(|\psi_0\rangle\).

Parameters
  • nbits (int) – the number of qubits in the system

  • hamming_weight (int) – the Hamming weight of the bit strings to mix

  • exact (optional, bool) – if set to True, the initial state will be described using a numpy array (which is not a scalable solution, but is faster and simpler for small systems).

  • tmax_psi_0 (optional, float) – the tmax for the first annealing (in the case of an inexact preparation). Default to 50.

Example:

from qat.opt import MixingFactory

# Mixes bit-strings of length 10 and Hamming weight 5
# This generates a numpy array of length 2^10 to describe the initial state (be careful :))
mixing = MixingFactory.bit_move(10, 5)
# Mixes bit-strings of length 10 and Hamming weight 5
# and prepares the initial state via a problem independent annealing
mixing = MixingFactory.bit_move(10, 5, exact=False)
# Mixes bit-strings of length 10 and Hamming weight 5
# and prepares the initial state via a problem independent annealing
# of duration 70
mixing = MixingFactory.bit_move(10, 5, exact=False, tmax_psi_0=70)
<stdin>:5: DeprecationWarning: `np.complex` is a deprecated alias for the builtin `complex`. To silence this warning, use `complex` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.complex128` here.
Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations
Returns

a pair (structure preparing the initial state, mixing Hamiltonian). The return type of the first entry depends on the exact parameter.

Return type

(SimpleInitialState/IndependentAnnealing, Observable)

You can also define you own mixing by creating a class that implements the following (very straightforward) interface:

class qat.opt.mixing_factory.InitialStateBuilder

A concept for an object that can modify an analog job in order to inject some particular initial state in the computation.

This interface requires implementation of a single method inject_initial_state that can modify a Job object in order to setup the correct initial state for a quantum annealing.

This method might simply modify the psi_0 field of the job, or prepend instructions to the job’s schedule.

abstract inject_initial_state(annealing_job)

Modifies in place an annealing job in order to prepare a particular initial state.

Parameters

annealing_job (qat.core.Job) – the job to modify

This interface is used by the following specializations:

class qat.opt.mixing_factory.SimpleInitialState(psi_0)

A very simple initial state builder that simply forces the initial state to some value.

Parameters

psi_0 (str, numpy.ndarray) – the initial state to prepare

inject_initial_state(annealing_job)

Modifies in place an annealing job in order to prepare a particular initial state.

Parameters

annealing_job (qat.core.Job) – the job to modify

class qat.opt.mixing_factory.IndependentAnnealing(start_ham, end_ham, tmax, true_psi_0)

An implementation of an InitialStateBuilder that prepend a problem independent annealing to the current annealing in order to prepare the correct initial state.

Parameters
  • start_ham (Observable) – the mixing Hamiltonian of the annealing

  • end_ham (Observable) – the target Hamiltonian of the annealing

  • tmax (float/expression) – the annealing time (or a variable/expression)

  • true_psi_0 (str/numpy.ndarray) – the initial state of the initial annealing

inject_initial_state(annealing_job)

Modifies in place an annealing job in order to prepare a particular initial state.

Parameters

annealing_job (qat.core.Job) – the job to modify

Encoding NP-hard Problems

We present here classes for encoding some of the famous NP problems. An instrinsic feature of these problems is that they can be formulated as minimization or maximization problems, i.e. with a cost function. At the same time finding the lowest energy of a physical system, represented by a cost Hamiltonian, is also a minimization problem. Therefore, we can represent the cost function of an NP problem by a cost Hamiltonian. Such a Hamiltonian, given in an Ising form can then be annealed using Simulated Annealing (SA).

Furthermore, problems formulated as Quadratic Unconstrained Binary Optimisation (QUBO), can also be annealed, since we can translate them to Ising via our to_ising() method.

Encoding NP-hard Problems

We present here classes for encoding some of the famous NP problems. An instrinsic feature of these problems is that they can be formulated as minimization or maximization problems, i.e. with a cost function. At the same time finding the lowest energy of a physical system, represented by a cost Hamiltonian, is also a minimization problem. Therefore, we can represent the cost function of an NP problem by a cost Hamiltonian. Such a Hamiltonian, given in an Ising form can then be annealed using Simulated Annealing (SA).

Furthermore, problems formulated as Quadratic Unconstrained Binary Optimisation (QUBO), can also be annealed, since we can translate them to Ising via our to_ising() method.

Note

The Simulated Quantum Annealing is not available in myQLM. QUBO and Ising problems can still be used to construct QAOA Ansätze using the .qaoa_ansatz method or solved via Simulated Annealing.

Unconstrained Graph Problems

Max cut

class qat.opt.max_cut.MaxCut(graph, **kwargs)

Specialization of the Ising class for Max Cut.

This class allows for the encoding of a Max Cut problem for a given graph. The method produce_j_h_and_offset() is automatically called. It calculates the coupling matrix \(J\), magnetic field \(h\) and Ising energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class Ising and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance - see the Max Cut notebook.

Reference

Maximum cut, Theoretical physics, Wikipedia.

import networkx as nx
from qat.opt import MaxCut

graph = nx.full_rary_tree(2, 2**8)

maxcut = MaxCut(graph)

print("To anneal the problem, the solver would need "
      + str(len(graph.nodes())) + " spins.")
To anneal the problem, the solver would need 256 spins.
Parameters

graph (networkx.Graph) – a networkx graph

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.max_cut.produce_j_h_and_offset(graph)

Returns the \(J\) coupling matrix of the problem, along with the magnetic field \(h\) and the Ising energy offset.

Parameters

graph (networkx.Graph) – a networkx graph

Graph Partitioning

class qat.opt.graph_partitioning.GraphPartitioning(graph, A, B=1, **kwargs)

Specialization of the Ising class for Graph Partitioning.

This class allows for the encoding of a Graph Partitioning problem for a given graph. The method produce_j_h_and_offset() is automatically called. It computes the coupling matrix \(J\), magnetic field \(h\) and Ising energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class Ising and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance - see the Graph Partitioning notebook.

For right encoding we need \(\frac { A } { B } \geq \frac { min(2D, N) } { 8 }\) with \(D\) - the maximal degree of a node in the graph and \(N\) - the number of nodes.

Reference

“Ising formulations of many NP problems”, A. Lucas, 2014 - Section 2.2.

import numpy as np
import networkx as nx
from qat.opt import GraphPartitioning

graph = nx.Graph()
graph.add_nodes_from(np.arange(10))
graph.add_edges_from([(0,1), (0,4), (0,6), (1,2), (1,4), (1,7), (2,3), (2,5), (2,8),
                      (3,5), (3,9), (4,6), (4,7), (5,8), (5,9), (6,7), (7,8), (8,9)])

B = 2
A = 5

graph_partitioning_problem = GraphPartitioning(graph, A, B=B)

print("To anneal the problem, the solver would need "
      + str(len(graph.nodes())) + " spins.")
To anneal the problem, the solver would need 10 spins.
Parameters
  • graph (networkx.Graph) – a networkx graph

  • A (double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.graph_partitioning.produce_j_h_and_offset(graph, A, B=1)

Returns the \(J\) coupling matrix of the problem, along with the magnetic field \(h\) and the Ising energy offset. For right encoding we need \(\frac{A}{B} \geq \frac{min(2D, N)}{8}\) with \(D\) - the maximal degree of a node in the graph and \(N\) - the number of nodes.

Parameters
  • graph (networkx.Graph) – a networkx graph

  • A (double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

Constrained Graph Problems

Graph Colouring

class qat.opt.graph_colouring.GraphColouring(graph, number_of_colours, **kwargs)

Specialization of the QUBO class for Graph Colouring.

This class allows for the encoding of a Graph Colouring problem for a given graph and a number of colours. The method produce_q_and_offset() is automatically called. It computes the \(Q\) matrix and QUBO energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class QUBO and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance.

Reference

“Ising formulations of many NP problems”, A. Lucas, 2014 - Section 6.1.

import numpy as np
import networkx as nx
from qat.opt import GraphColouring

graph = nx.Graph()
graph.add_nodes_from(np.arange(4))
graph.add_edges_from([(0,1), (0,2), (1,2), (1,3), (2,3)])
number_of_colours = 3

graph_colouring_problem = GraphColouring(graph, number_of_colours)

print("To anneal the problem, the solver would need "
      + str(len(graph.nodes()) * number_of_colours) + " spins.")
To anneal the problem, the solver would need 12 spins.
Parameters
  • graph (networkx.Graph) – a networkx graph

  • number_of_colours (int) – the number of colours

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.graph_colouring.produce_q_and_offset(graph, number_of_colours)

Returns the \(Q\) matrix and the offset energy of the problem.

Parameters
  • graph (networkx.Graph) – a networkx graph

  • number_of_colours (int) – the number of colours

K-Clique

class qat.opt.k_clique.KClique(graph, K, A, B=1, **kwargs)

Specialization of the QUBO class for K-Clique.

This class allows for the encoding of a K-Clique problem for a given graph and positive factors \(K, A, B\). The method produce_q_and_offset() is automatically called. It computes the \(Q\) matrix and QUBO energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class QUBO and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance - see the KClique notebook.

For a right encoding, one should ensure that \(A > B * K\).

Reference

“Ising formulations of many NP problems”, A. Lucas, 2014 - Section 2.3.

import numpy as np
import networkx as nx
from qat.opt import KClique

graph = nx.Graph()
graph.add_nodes_from(np.arange(6))
graph.add_edges_from([(0,1), (0,2), (0,3), (0,4), (0,5), (1,2), (1,5)])

B = 1
K = 3
A = B*K + 1

k_clique_problem = KClique(graph, K, A, B=B)

print("To anneal the problem, the solver would need "
      + str(len(graph.nodes())) + " spins.")
To anneal the problem, the solver would need 6 spins.
Parameters
  • graph (networkx.Graph) – a networkx graph

  • K (int) – the size of the clique

  • A (double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.k_clique.produce_q_and_offset(graph, K, A, B=1)

Returns the \(Q\) matrix and the offset energy of the problem. The constant \(A\) should be bigger than \(K*B\) for a right encoding. They are also all positive.

Parameters
  • graph (networkx.Graph) – a networkx graph

  • K (int) – the size of the clique

  • A (double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

Vertex Cover

class qat.opt.vertex_cover.VertexCover(graph, A=2, B=1, **kwargs)

Specialization of the QUBO class for Vertex Cover.

This class allows for the encoding of a Vertex Cover problem for a given graph and positive constants \(A\) and \(B\). The method produce_q_and_offset() is automatically called. It computes the \(Q\) matrix and QUBO energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class QUBO and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance - see the Vertex Cover notebook.

For a right encoding, one should ensure that \(A > B\).

Reference

“Ising formulations of many NP problems”, A. Lucas, 2014 - Section 4.3.

import numpy as np
import networkx as nx
from qat.opt import VertexCover

graph = nx.Graph()
graph.add_nodes_from(np.arange(6))
graph.add_edges_from([(0,1), (0,2), (0,3), (0,4), (0,5), (1,5)])
A = 2
B = 1

vertex_cover_problem = VertexCover(graph, A=A, B=B)

print("To anneal the problem, the solver would need "
       + str(len(graph.nodes())) + " spins.")
To anneal the problem, the solver would need 6 spins.
Parameters
  • graph (networkx.Graph) – a networkx graph

  • A (optional, double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied, default is 2. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.vertex_cover.produce_q_and_offset(graph, A=2, B=1)

Returns the \(Q\) matrix and the offset energy of the problem. The constant \(A\) should be bigger than \(B\) for a right encoding. They are also both positive.

Parameters
  • graph (networkx.Graph) – a networkx graph

  • A (optional, double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied, default is 2. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

Other problems

Number Partitioning

class qat.opt.number_partitioning.NumberPartitioning(array_of_numbers, **kwargs)

Specialization of the Ising class for Number Partitioning.

This class allows for the encoding of a Number Partitioning problem for a given array of numbers. The method produce_j_h_and_offset() is automatically called. It computes the coupling matrix \(J\), magnetic field \(h\) and Ising energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class Ising and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance - see the Number Partitioning notebook.

Reference

“Ising formulations of many NP problems”, A. Lucas, 2014 - Section 2.1.

import numpy as np
from qat.opt import NumberPartitioning

array_of_numbers_size = np.random.randint(low=1, high=10000, size=1)[0]
array_of_numbers = np.random.randint(low=1, high=10000, size=array_of_numbers_size)

number_partitioning_problem = NumberPartitioning(array_of_numbers)

print("To anneal the problem, the solver would need "
      + str(array_of_numbers_size) + " spins.")
To anneal the problem, the solver would need 2536 spins.
Parameters

numbers_array (1D numpy array) – an array with all the numbers we want to partition

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.number_partitioning.produce_j_h_and_offset(array_of_numbers)

Returns the \(J\) coupling matrix of the problem, along with the magnetic field \(h\) and the Ising energy offset.

Parameters

numbers_array (1D numpy array) – an array with all the numbers we want to partition

Binary Integer Linear Programming

class qat.opt.binary_linear_integer_programming.BILP(c, S, b, A, B=1, **kwargs)

Specialization of the QUBO class for Binary Integer Linear Programming (BILP).

This class allows for the encoding of a BILP problem from a given matrix \(S\), vectors \(b\) and \(c\) and positive constants \(A\) and \(B\). The aim is to maximise \(c * x\) subject to \(x\) obeying \(S * x = b\). The method produce_q_and_offset() is automatically called. It computes the \(Q\) matrix and QUBO energy offset corresponding to the Hamiltonian representation of the problem, as described in the reference. These are stored in the parent class QUBO and would be needed if one wishes to solve the problem through Simulated Annealing (SA), for instance.

For a right encoding, one should ensure that \(A \gg B\) and \(A > 0, B > 0\).

Reference

“Ising formulations of many NP problems”, A. Lucas, 2014 - Section 3.

import numpy as np
from qat.opt import BILP

c = np.array([0, 1, 1, 1])
S = np.array([[1, 0, 1, 1], [0, 1, 0, 1]])
b = np.array([1, 1])

B = 1
A = 10 * B

bilp_problem = BILP(c, S, b, A, B=B)

print("To anneal the problem, the solver would need " + str(len(c)) + " spins.")
To anneal the problem, the solver would need 4 spins.
Parameters
  • c (1D numpy array of size N) – a specified vector \(c\). We want to maximize \(c * x\).

  • S (2D numpy array of size m*N) – the matrix, for which \(S * x = b\). This equation is our constraint.

  • b (1D numpy array of size m) – a specified vector \(b\) obeying the constraint \(S * x = b\)

  • A (double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1

get_best_parameters()

This method returns a dictionary with the best found parameters (after benchmarking) for simulated quantum annealing (SQA), available in the QLM. However, the temperature parameters could also be used for simulated annealing (SA).

Returns

6-key dictionary containing

  • n_monte_carlo_updates (int) - the number of Monte Carlo updates

  • n_trotters (int) - the number of “classical replicas” or “Trotter replicas”

  • gamma_max (double) - the starting magnetic field

  • gamma_min (double) - the final magnetic field

  • temp_max (double) - the starting temperature

  • temp_min (double) - the final temperature

qat.opt.binary_linear_integer_programming.produce_q_and_offset(c, S, b, A, B=1)

Returns the \(Q\) matrix and the offset energy of the problem. For right encoding \(A \gg B\) and \(A > 0, B > 0\).

Parameters
  • c (1D numpy array of size N) – a specified vector \(c\). We want to maximize \(c * x\).

  • S (2D numpy array of size m*N) – the matrix, for which \(S * x = b\). This equation is our constraint.

  • b (1D numpy array of size m) – a specified vector \(b\) obeying the constraint \(S * x = b\)

  • A (double) – a positive constant by which the terms inside \(H_A\) from \(H = H_A + H_B\) are multiplied. This equation comes from the Hamiltonian representation of the problem.

  • B (optional, double) – similar to \(A\), \(B\) is a positive factor for the \(H_B\) terms, default is 1