Discrete optimization problems are ubiquitous in almost any enterprise and many of these are known to be NP-hard (Lenstra and Rinnooy Kan, 1979). The objective of such problems is to find the minimum of a real-valued function f(v ₀, …, v _N−1) (the cost function) over a set of discrete variables v _k . The search space is restricted by hard constraints, which are commonly presented as equalities such as g(v ₀, …, v _N−1) = 0 or inequalities as h(v ₀, …, v _N−1) > 0. Besides using classical heuristics (Dorigo and Di Caro, 1999; Melnikov, 2005) and machine learning methods (Mazyavkina et al., 2021) to solve these problems, there is a growing interest in applying quantum computation (Au-Yeung et al., 2023). A common approach for realizing this consists of first encoding the cost function f in a Hamiltonian H such that a subset of eigenvectors of H represents elements in the domain of f and the eigenvalues are the respective values of f: H ψ = f v 0 , … , v N − 1 ψ . (1) In such an encoding, the ground state of H is the solution to the optimization problem. Having obtained a Hamiltonian formulation, one can use a variety of quantum algorithms to find the ground state including adiabatic quantum computing (Farhi et al., 2000) and variational approaches, for instance, the quantum/classical hybrid quantum approximate optimization algorithm (QAOA) (Farhi et al., 2014) or generalizations thereof such as the quantum alternating operator ansatz (Hadfield et al., 2019). On the hardware side, these algorithms can run on gate-based quantum computers, quantum annealers, or specialized Ising machines (Mohseni et al., 2022).

In the current literature, almost all Hamiltonians for optimization are formulated as Quadratic Unconstrained Binary Optimization (QUBO) problems (Kochenberger et al., 2014). The success of QUBO reflects the strong hardware limitations of current devices, where multiqubit interactions are not available and must be decomposed into two-qubit interactions using ancilla qubits. Moreover, quantum algorithms with the dynamical implementation of hard constraints (Hen and Sarandy, 2016; Hen and Spedalieri, 2016) require driver terms that can be difficult to design and implement on quantum computers. Hence, hard constraints are usually included as energy penalizations of QUBO Hamiltonians. The prevalence of QUBO has also increased the popularity of one-hot encoding, a particular way of mapping (discrete) variables to eigenvalues of spin operators (Lucas, 2014), since this encoding allows for Hamiltonians with low-order interactions which is especially appropriate for QUBO problems.

However, compelling alternatives to QUBO and one-hot encoding have been proposed in recent years. A growing number of platforms are exploring high-order interactions (Chancellor et al., 2017; Lu et al., 2019; Schöndorf and Wilhelm, 2019; Wilkinson and Hartmann, 2020; Menke et al., 2021; 2022; Dlaska et al., 2022; Pelegrí et al., 2022; Glaser et al., 2023), while the Parity architecture (Lechner, 2020; Fellner et al., 2022; Ender et al., 2023) (a generalization of the LHZ architecture (Lechner et al., 2015)) allows the mapping of arbitrary-order interactions to qubits that require only local connectivity. The dynamical implementation of constraints has also been investigated (Hadfield et al., 2019; 2017; Fuchs et al., 2022; Zhu et al., 2023), including the design of approximate drivers (Wang et al., 2020; Sawaya et al., 2022) and compilation of constrained problems within the Parity architecture (Drieb-Schön et al., 2023). Moreover, simulated and experimental results have shown that alternative encodings outperform the traditional one-hot approach (Chancellor, 2019; Sawaya et al., 2020; Chen et al., 2021; Plewa et al., 2021; Tamura et al., 2021; Glos et al., 2022; Stein et al., 2023). Clearly, alternative formulations for Hamiltonians need to be explored but when the Hamiltonian has been expressed in QUBO using one-hot encoding, it is not trivial to switch to other formulations. Automatic tools to explore different formulations would therefore be highly beneficial.

We present a library of more than 20 problems which is intended to facilitate the Hamiltonian formulation beyond QUBO and one-hot encoding. We build upon the recent work of Sawaya et al. (2022) by making use of the encoding-independent approach to revisit common problems in the literature. With this approach, the problems can be encoded trivially using any spin encoding. We also provide a summary of the most popular encodings. Possible constraints of the problems are identified and presented separately from the cost function so that dynamic implementation of the constraints can also be easily explored. Two additional subgoals that are addressed in this library are:

• Meta parameters/choices: We present and review the most important choices made in the process of mapping optimization problems in a mathematical formulation to spin Hamiltonians. These mainly include the encodings, which can greatly influence the computational cost and performance of the optimization, but also free meta parameters or the use of auxiliary variables. These degrees of freedom are a consequence of the fact that the optimal solution is typically encoded only in the ground state. Other low-energy eigenstates encode good approximations to the optimal solution and it can be convenient to make approximations so that the solution corresponds to these states (Montanez-Barrera et al., 2022).

In practice, many optimization problems are not purely discrete but involve real-valued parameters and variables. Thus, the encoding of real-valued problems to discrete optimization problems (discretization) as an intermediate step is discussed in Section 7.4.

The focus of this review is on optimization problems which can be formulated as diagonal Hamiltonians written as sums and products of Pauli-z-matrices. This subset of Hamiltonians is usually not applicable to quantum systems or quantum simulations. For an introduction to these more general Hamiltonians, we refer to the reviews of Georgescu et al. (2014) for physics problems and McArdle et al. (2020); Cao et al. (2019) for quantum chemistry simulations.

After introducing the notation used throughout the text in Section 2, we present a list of encodings in Section 3. Section 4 reviews the Parity architecture and Section 5 discusses encoding constraints. Section 6 functions as a manual on how to bring optimization problems into a form that can be solved with a quantum computer. Section 7 contains a library of optimization problems which are classified into several categories and Section 8 lists the building blocks used in the formulation of these (and many other) problems. Section 9 offers conclusions.

A discrete set of real numbers is a countable subset U ⊂ R without an accumulation point. Discrete sets are denoted by uppercase Latin letters, except the letter G, which we reserve for graphs. A discrete variable is a variable ranging over a discrete set R and will be represented by lowercase Latin letters, mostly v or w. Elements of R are denoted by lowercase Greek letters.

If a discrete variable has range 0,1 we call it binary or boolean. Binary variables will be denoted by the letter x. Similarly, a variable with range − 1,1 will be called a spin variable and the letter s will be reserved for these. There is an invertible mapping from a binary variable x to a spin variable s: x ↦ s : = 2 x − 1 . (2) For a variable v with range R we follow Chancellor (2019) and Sawaya et al. (2022) and define the value indicator function to be δ v α = 1 i f v = α 0 i f v ≠ α , (3) where α ∈ R.

We also consider optimization problems for continuous variables. A variable v will be called continuous, if its range is given by R d for some d ∈ Z > 0 .

An optimization problem O is a triple (V, f, C), where.

1. V : = v i i = 0 , … , N − 1 is a finite set of variables.

The goal for an optimization problem O = (V, f, C) is to find an extreme value y _ex of f, such that all of the constraints are satisfied at y _ex .

Discrete optimization problems can often be stated in terms of graphs or hypergraphs (Berge, 1987). A graph is a pair (V, E), where V is a finite set of vertices or nodes and E ⊂ v i , v j | v i ≠ v j ∈ V is the set of edges. An element v i , v j ∈ E is called an edge between vertex v _i and vertex v _j . Note that a graph defined like this can neither have loops, i.e., edges beginning and ending at the same vertex, nor can there be multiple edges between the same pair of vertices. Given a graph G = (V, E), its adjacency matrix is the symmetric binary |V|×|V| matrix A with entries A i j = 1 , if v i , v j ∈ E , 0 , else . (6)

A hypergraph is a generalization of a graph in which we allow edges to be adjacent to more than two vertices. That is, a hypergraph is a pair H = (V, E), where V is a finite set of vertices and E ⊆ ⋃ i = 1 | V | { v j 1 , … , v j i } | v j k ∈ V a n d v j k ≠ v j ℓ ∀ k , ℓ = 1 , ⋯ i (7) is the set of hyperedges.

We reserve the word qubit for physical qubits. To go from an encoding-independent Hamiltonian to a quantum program, binary or spin variables become Pauli-z-matrices which act on the corresponding qubits.

For many problems, the cost function and the problem constraints can be represented in terms of two fundamental building blocks: the value of the integer variable v and the value indicator δ v α defined in Eq. 3. When expressed in terms of these building blocks, Hamiltonians are more compact, and recurring terms can be identified across many different problems. Moreover, quantum operators are not present at this stage: an encoding-independent Hamiltonian is just a cost function in terms of discrete variables, which eases access to quantum optimization to a wider audience. Encoding variables and the choice of quantum algorithms can come later.

Representation of the building blocks in terms of Ising operators depends on the chosen encoding. An encoding is a function that associates eigenvectors of the σ _z operator with specific values of a discrete variable v: | s 0 , … , s N − 1 〉 → v , s i = ± 1 , (8) where the spin variables s _i are the eigenvalues of σ z ( i ) operators. The encodings are also usually defined in terms of binary variables x _i which are related to Ising variables according to Eq. 2.

A summary of the encodings is presented in Figure 1. Some encodings are dense, in the sense that every quantum state |s ₀, …, s _N−1⟩ encodes some value of the variable v. Other encodings are sparse because only a subset of the possible quantum states are valid states. The valid subset is generated by adding a core term ¹ , i.e., a penalty term for constraints that need to be enforced in order to decode the variable uniquely in the Hamiltonian for every sparsely encoded variable. In general, dense encodings require fewer qubits, but sparse encodings have simpler expressions for the value indicator δ v α , and are therefore favorable for avoiding higher-order interactions. This is because δ v α needs to check a smaller number of qubit states to know whether the variable v has value α or not, whereas dense encodings require the state of every qubit in the register (Sawaya et al., 2020; 2022).

The strong hardware limitations of noisy intermediate-scale quantum (NISQ) (Preskill, 2018) devices have made sparse encodings (especially one-hot) the standard approach to problem formulation. This is mainly because the basic building blocks (value and value indicator) are of linear or quadratic order in the spin variables in these encodings. The low connectivity of qubit platforms requires Hamiltonians in the QUBO formulation and high-order interactions are expensive when translated to QUBO Kochenberger et al. (2014). However, different choices of encodings can significantly improve the performance of quantum algorithms (Chancellor, 2019; Sawaya et al., 2020; Chen et al., 2021; Di Matteo et al., 2021; Tamura et al., 2021), by reducing the search space or generating a smoother energy landscape.

One way this difference between encodings manifests itself is in the number of spin flips of physical qubits needed to change a variable into another valid value (Berwald et al., 2023). If this number is larger than one, there are local minima separated by invalid states penalized with a high cost which can impede the performance of the optimization. On the other hand, such an energy-landscape might offer some protection against errors (Pastawski and Preskill, 2016; Fellner et al., 2022). Furthermore, other fundamental aspects of the algorithms, such as circuit depth and energy scales can be greatly improved outside QUBO (Ender et al., 2022; Drieb-Schön et al., 2023; Fellner et al., 2023; Messinger et al., 2023), prompting us to look for alternative formulations. The Parity Architecture is a paradigm for solving quantum optimization problems (Lechner et al., 2015; Ender et al., 2023) that does not rely on the QUBO formulation, allowing a wide number of options for formulating Hamiltonians. The architecture is based on the Parity transformation, which remaps Hamiltonians onto a 2D grid requiring only local connectivity of the qubits. The absence of long-range interactions enables high parallelizability of quantum algorithms and eliminates the need for costly and time-consuming SWAP gates, which helps to overcome two of the main obstacles of quantum computing: limited coherence time and the poor connectivity of qubits within a quantum register.

The Parity transformation creates a single Parity qubit for each interaction term in the (original) logical Hamiltonian: J i , j , … σ z i σ z j ⋯ → J i , j , … σ z i , j , … , (34) where the interaction strength J _i,j,… is now the local field of the Parity qubit σ z ( i , j , … ) . This facilitates addressing high-order interactions and frees the problem formulation from the QUBO approach. The equivalence between the original logical problem and the Parity-transformed problem is ensured by adding three- and four-body constraints and placing them on a 2D grid such that only neighboring qubits are involved in the constraints. The mapping of a logical problem into the regular grid of a Parity chip can be realized by the Parity compiler (Ender et al., 2023). Although the Parity compilation of the problem may require a larger number of physical qubits, the locality of interactions on the grid allows for higher parallelizability of quantum algorithms. This allows constant depth algorithms (Lechner, 2020; Unger et al., 2022) to be implemented with a smaller number of gates (Fellner et al., 2023).

The following toy example, summarized in Figure 2, shows how a Parity-transformed Hamiltonian can be solved using a smaller number of qubits when the original Hamiltonian has high-order interactions. Given the logical Hamiltonian H = σ z 1 σ z 2 + σ z 2 σ z 4 + σ z 1 σ z 5 + σ z 1 σ z 2 σ z 3 + σ z 3 σ z 4 σ z 5 , (35) the corresponding QUBO formulation requires seven qubits, including two ancillas for decomposing the three-body interactions, and the total number of two-body interactions is 14. The embedding of the QUBO problem on quantum hardware may require additional qubits and interactions depending on the chosen architecture. Instead, the Parity-transformed Hamiltonian only consists of six Parity qubits with local fields and two four-body interactions between close neighbors.

It is not yet clear what the best Hamiltonian representation is for an optimization problem. The answer will probably depend strongly on the particular use case and will take into account not only the number of qubits needed but also the smoothness of the energy landscape, which has a direct impact on the performance of quantum algorithms (King et al., 2019).

In this section, we review how to implement the hard constraints associated with the problem, assuming that the encodings of the variables have already been chosen. Hard constraints c(v ₁, …, v _N ) = K often appear in optimization problems, limiting the search space and making problems even more difficult to solve. We consider polynomial constraints of the form c v 1 , … , v N = ∑ g i v i + ∑ g i , j v i v j + ∑ g i , j , k v i v j v k + … , (36) which remain polynomial after replacing the discrete variables v _i with any encoding. The coefficients g _i , g _i,j , … depend on the problem and its constraints. Even if the original problem is unconstrained, the use of sparse encodings such as one-hot or domain wall imposes hard constraints on the quantum variables.

In general, constraints can be implemented dynamically (Hen and Sarandy, 2016; Hen and Spedalieri, 2016) (exploring only quantum states that satisfy the constraints) or as extra terms H _c in the Hamiltonian, such that eigenvectors of H are also eigenvectors of H _c and the ground states of H _c correspond to elements in the domain of f that satisfy the constraint. These can be incorporated as a penalty term H _c into the Hamiltonian that penalizes any state outside the desired subspace: H c = A c v 1 , … , v N − K 2 , (37) or H c = A c v 1 , … , v N − K , (38) in the special case that c(x ₁, …, x _N ) ≥ K is satisfied. The constant A must be large enough to ensure that the ground state of the total Hamiltonian satisfies the constraint, but the implementation of large energy scales lowers the efficiency of quantum algorithms (Lanthaler and Lechner, 2021) and additionally imposes a technical challenge. Moreover, extra terms in the Hamiltonian imply additional overhead of computational resources, especially for squared terms such as in Eq. 37. The determination of the optimal energy scale is an important open problem. For some of the problems in the library, we provide an estimation of the energy scales (cf. also Section 8.4.2).

Quantum algorithms for finding the ground state of Hamiltonians, such as QAOA or quantum annealing, require driver terms U _drive = exp(−itH _drive) that spread the initial quantum state to the entire Hilbert space. Dynamical implementation of constraints employs a driver term that only explores the subspace of the Hilbert space that satisfies the constraints. Given an encoded constraint in terms of Ising operators σ z ( i ) : c σ z 1 , … , σ z N = ∑ g i σ z i + ∑ g i , j σ z i σ z ( j ) + ⋯ = K , (39) a driver Hamiltonian H _drive that commutes with c ( σ z ( 1 ) , … , σ z ( N ) ) will be restricted to the valid subspace provided the initial quantum state satisfies the constraint: c | ψ 0 〉 = K | ψ 0 〉 U drive c | ψ 0 〉 = U drive K | ψ 0 〉 c U drive | ψ 0 〉 = K U drive | ψ 0 〉 . (40)

In general, the construction of constraint-preserving drivers depends on the problem (Hen and Sarandy, 2016; Hen and Spedalieri, 2016; Hadfield et al., 2017; Chancellor, 2019; Bärtschi and Eidenbenz, 2020; Wang et al., 2020; Fuchs et al., 2021; Bakó et al., 2022). Approximate driver terms have been proposed that admit some degree of leakage and may be easier to construct (Sawaya et al., 2022). Within the Parity Architecture, each term of a polynomial constraint is a single Parity qubit (Lechner et al., 2015; Ender et al., 2023). This implies that for the Parity Architecture the polynomial constraints are simply the conservation of magnetization between the qubits involved: ∑ i g i σ z i + ∑ i , j g i , j σ z i σ z j + ⋯ = K → ∑ u ∈ c g u σ z u = K , (41) where the Parity qubit σ z ( u ) represents the logical qubits product σ z ( u 1 ) σ z ( u 2 ) … σ z ( u n ) and we can sum over all these products that appear in the constraint c. A driver Hamiltonian based on exchange (or flip-flop) terms σ + ( u ) σ − ( w ) + h . c ., summed over all pairs of products u, w ∈ c, preserves the total magnetization and explores the complete subspace where the constraint is satisfied (Drieb-Schön et al., 2023). The decision tree for encoding constraints is presented in Figure 3.

In this section, we present an example of the complete procedure to go from the encoding-independent formulation to the spin Hamiltonian that has to be implemented in the quantum computer, using an instance of the Clustering problem (Section 7.1.1).

Every problem in this library includes a Problem description, indicating the required inputs for defining a problem instance. In the case of the clustering problem, a problem instance is defined from the number of clusters K we want to create, N objects with weights w _i , and distances d _i,j between the objects. Two different types of discrete variables are required, variables v _i = 1, …, K (i = 1, …, N) indicate to which of the K possible clusters the node i is assigned, and variables y _j = 1, …, W _max track the weight in cluster j.

The cost function of the problem only depends on δ v i k , the value indicators of variables v _i . In this case, δ v i k indicates whether the node associated with the discrete variable v _i belongs to the cluster k or not. If two nodes i, j belong to the same cluster k, then the cost function increases by d _i,j , giving the total cost function: f = ∑ k = 1 K ∑ i < j d i , j δ v i k δ v j k . (42) We can choose any encoding for the variables v _i . For example, if we decide to use binary or Gray encoding, a variable v _i requires D = log₂(K) qubits { x u ( i ) } u = 1 D (we assume K = 2 ^D for simplicity, but the resource estimate will not change significantly if this is not the case). From Eq. 12 we see that δ v i k is a polynomial of order D in variables x u ( i ) with D j terms of order j and 2 ^D = K terms in total, so the product δ v i k δ v j k in the cost function will have 2 ^D × 2 ^D = 2^2D = K ² terms, with orders between zero and 2D. This product is summed over all the N ( N − 1 ) 2 pairs i, j, so we can say that O ( N 2 K 2 ) terms are required for the cost function. High-order interactions may be prohibitive for some hardware platforms, but if multiqubit gates are available, binary encoding offers an important reduction in the number of qubits.

Alternatively, we can choose a sparse encoding for the variables v _i . Using one-hot encoding, we need K qubits per node, so NK qubits are required for the cost function. The product δ v i k δ v j k is just a two-qubit interaction x k ( i ) x k ( j ) , so the cost function requires O ( K N 2 ) terms. The K qubits associated with a variable v _i must satisfy the core constraint: c i = ∑ u = 1 K x u i = 1 . (43) if this constraint is implemented as an energy penalization ( c i − 1 ) 2 , then O ( K 2 ) terms are added to the spin Hamiltonian for each variable, so O ( N K 2 ) terms in total are associated with the constraint. The advantage of sparse encodings is that low-order interactions are needed for the value indicator functions (in this case, the maximum order is two).

Besides the core constraints associated with the encodings, the clustering problem includes K additional constraints (one per each cluster). The total weight of nodes in any cluster cannot exceed a problem instance specific maximal value W _max: ∑ i w i δ i k ≤ W max , ∀ k ∈ { 1 , … , K } . (44) For the kth cluster, this constraint can be expressed in terms of auxiliary variables y _k : y k = ∑ i w i δ v i k . (45) Variables y _k are discrete variables in the range 0, …, W _max. Because the value indicators of y _k are not necessary for the constraints, we can use a dense encoding such as binary without dealing with the high-order interactions associated with value indicators of dense encodings. The variables y _k require K log₂(W _max) qubits if we use binary encoding, or KW _max if we use one-hot. These constraints can also be implemented as energy penalizations in the Hamiltonian or can be encoded in the driver term.

The complete procedure for obtaining the spin Hamiltonian is outlined in Figure 4. We emphasize that the optimal encodings and constraint implementations depend on the details of the hardware, such as native gates, connectivity, and the number of qubits. Moreover, the efficiency of quantum algorithms is also related to the smoothness of the energy landscape, and some encodings can provide better results even though they require more qubits (see, for example, Tamura et al., 2021).

The problems included in the library are classified into four categories: subsets, partitions, permutations, and continuous variables. These categories are defined by the role of the discrete variable and are intended to organize the library and make it easier to find problems but also to serve as a basis for the formulation of similar use cases. An additional category in Section 7.5 contains problems that do not fit into the previous categories but may also be important use cases for quantum algorithms. In Section 8 we include a summary of recurrent encoding-independent building blocks that are used throughout the library and could be useful in formulating new problems.