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

Combinatorial Problem

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 * (Z|[1]) +
0.125 * (ZZ|[0, 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.

Clause

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:

Var

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

evaluate(assignment)

Returns the boolean value of the variable in an assignment list

Parameters: assignment (list<bool>) – A list of boolean representing the value of each variables
Returns: The result of the evaluation
Return type: bool

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

QUBO

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.

Ising

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:

Circuit Generator

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.

Ansatz Factory

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:

Schedule Generator

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:

Mixing Factory

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:

Initial State Builder

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:

Simple Initial State

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

Independent Annealing

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.

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.
<stdin>:6: FutureWarning: adjacency_matrix will return a scipy.sparse array instead of a matrix in Networkx 3.0.

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

parse_result(result, inverse=False)

Returns the best approximated solution of the Max Cut problem from a list of samples

Parameters: result (BatchResult) – BatchResult containing a list of samples
Returns: The best partition among the samples with the maximum cut size
Return type: GraphPartitioningResult

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.
<stdin>:13: FutureWarning: adjacency_matrix will return a scipy.sparse array instead of a matrix in Networkx 3.0.

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

parse_result(result, inverse=False)

Returns the best approximated solution of the Graph Partitioning problem from a list of samples

Parameters: result (BatchResult) – BatchResult containing a list of samples
Returns: The best balanced partition among the samples with the minimum cut size
Return type: GraphPartitioningResult

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.
<stdin>:10: FutureWarning: adjacency_matrix will return a scipy.sparse array instead of a matrix in Networkx 3.0.

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

parse_result(result, inverse=False)

Returns the best approximated solution of the Graph Colouring problem from a list of samples

Parameters: result (BatchResult) – BatchResult containing a list of samples
Returns: The best partition among the samples thatrepresents the colour of each node
Return type: GraphPartitioningResult

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.
<stdin>:13: FutureWarning: adjacency_matrix will return a scipy.sparse array instead of a matrix in Networkx 3.0.
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

parse_result(result, inverse=False)

Returns the best approximated solution of the K-Clique problem from a list of samples

Parameters

result (BatchResult) – BatchResult containing a list of samples

Returns

The best partition among the samples with a K-clique as the first subset

Return type

GraphPartitioningResult

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.
<stdin>:11: FutureWarning: adjacency_matrix will return a scipy.sparse array instead of a matrix in Networkx 3.0.

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

parse_result(result, inverse=False)

Returns the approximated solution of the Vertex cut problem from a list of samples

Parameters: result (BatchResult) – BatchResult containing a list of samples

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

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

Declaring and manipulating combinatorial problems

Generic combinatorial optimization

Quadratic problems: QUBO and Ising

Generating the QAOA Ansatz

Generating QA jobs

Encoding NP-hard Problems

Unconstrained Graph Problems

Constrained Graph Problems

Other problems