PRX QUANTUM 5, 020327 (2024)

Quantum-Informed Recursive Optimization Algorithms

Jernej Rudi Finžgar ,1,2,*,‡ Aron Kerschbaumer,1,2,†,‡ Martin J.A. Schuetz,3,4 Christian B. Mendl,2,5
and Helmut G. Katzgraber3
1
BMW AG, Munich, Germany
2
Department of Computer Science, Technical University Munich, School of CIT, Garching, Germany
3
Amazon Quantum Solutions Lab, Seattle, Washington, USA
4
AWS Center for Quantum Computing, Pasadena, California, USA
5
Technical University of Munich, Institute for Advanced Study, Garching, Germany

(Received 21 September 2023; accepted 21 February 2024; published 3 May 2024)

We propose and implement a family of quantum-informed recursive optimization (QIRO) algorithms

for combinatorial optimization problems. Our approach leverages quantum resources to obtain informa-
tion that is used in problem-specific classical reduction steps that recursively simplify the problem. These
reduction steps address the limitations of the quantum component (e.g., locality) and ensure solution fea-
sibility in constrained optimization problems. Additionally, we use backtracking techniques to further
improve the performance of the algorithm without increasing the requirements on the quantum hard-
ware. We showcase the capabilities of our approach by informing QIRO with correlations from classical
simulations of shallow circuits of the quantum approximate optimization algorithm, solving instances
of maximum independent set and maximum satisfiability problems with hundreds of variables. We also
demonstrate how QIRO can be deployed on a neutral atom quantum processor to find large independent
sets of graphs. In summary, our scheme achieves results comparable to classical heuristics even with rela-
tively weak quantum resources. Furthermore, enhancing the quality of these quantum resources improves
the performance of the algorithms. Notably, the modular nature of QIRO offers various avenues for mod-
ifications, positioning our work as a template for a broader class of hybrid quantum-classical algorithms
for combinatorial optimization.

DOI: 10.1103/PRXQuantum.5.020327

I. INTRODUCTION the solution space properties of certain optimization prob-

lems such as the Z2 symmetry, or the overlap gap property
Quantum optimization has been identified as a promis-
[2], locality was shown to severely limit the performance
ing area of research towards practical quantum advantage.
of QAOA [3–6].
On noisy intermediate-scale quantum (NISQ) devices,
Because introducing nonlocal updates in quantum algo-
much effort has been dedicated to studying hybrid
rithms comes with additional hardware requirements,
quantum-classical algorithms such as the quantum approx-
applying nonlocal updates classically has been proposed
imate optimization algorithm (QAOA) [1]. Importantly,
as an alternative in recursive QAOA (RQAOA) [3,7].
QAOA is a local algorithm. This means that at any con-
Here, values of variables are iteratively frozen by round-
stant circuit depth, only qubits that are separated by less
ing the correlations between variables as measured in the
than a certain distance in the interaction graph of the opti-
quantum state prepared by QAOA. Hence, as variables
mization problem are able to communicate. Together with
are removed from the optimization problem, the distances
between nodes in the interaction graph are reduced iter-
*
jernej-rudi.finzgar@[Link] atively. As such, the nonlocal effects introduced via the
†
Present address: Institute of Science and Technology Austria new connections between previously unconnected nodes
(ISTA), Am Campus 1, 3400 Klosterneuburg, Austria. counterbalance the locality inherent to QAOA.
‡
These authors contributed equally. Building upon RQAOA, here we propose a family of
hybrid quantum-classical algorithms, dubbed quantum-
Published by the American Physical Society under the terms of
the Creative Commons Attribution 4.0 International license. Fur- informed recursive optimization (QIRO). In QIRO, infor-
ther distribution of this work must maintain attribution to the mation generated by quantum resources is used to recur-
author(s) and the published article’s title, journal citation, and sively reduce the size of the optimization problem by
DOI. means of problem-specific classical optimization routines;

2691-3399/24/5(2)/020327(18) 020327-1 Published by the American Physical Society

JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

Simplified problem performance of QIRO by means of large-scale simulations

of low-depth QAOA circuits with up to two hundred vari-
Quantum resources Update rules ables. Furthermore, we use QIRO to solve the maximum
independent set problem on the QuEra Aquila neutral atom
quantum device accessed via Amazon Braket. Finally,
we compare the performance of QIRO against commonly
used classical and quantum optimization techniques. This
allows us to assess the role of the quantum resources used
in QIRO, and showcase the value of classical subroutines.
The remainder of the paper is organized as follows.
FIG. 1. Schematic visualization of the core principles of the
quantum-informed recursive optimization (QIRO) algorithm pre- In Sec. II we introduce the two problem classes studied
sented in this work. Quantum resources (e.g., QAOA) are used in this work, in Sec. III we provide details about QIRO,
to obtain information, e.g., in the form of one-point and two- and in Sec. IV we study the performance of QIRO on
point correlations. This information is then used to simplify the the maximum independent set and maximum satisfiability
problem through problem-specific classical update rules. Here problems. We discuss the implications of our findings in
we exemplify the update rules with simple examples for the max- Sec. V and we conclude by suggesting potential extensions
imum independent set and satisfiability problems—details on the to our framework in Sec. VI.
update rules can be found in Sec. III B. The (simplified) problem
obtained by means of the classical update rules is then used to
restart the cycle. The algorithm terminates when the problem has II. PROBLEM CLASSES
been fully simplified.
Because problem-specific update rules are at the core of
QIRO, we first introduce the specific optimization prob-
see Fig. 1 for a schematic illustration. This scheme lems considered in this work. We formalize these prob-
allows us to leverage decades of research in (classical) lems, and provide mappings to quantum Hamiltonians Ĥc ,
combinatorial optimization [8], and tailor the classical sub- which can in turn be implemented on different quantum
routines to the particular optimization problem of interest, devices [10,11].
thereby enhancing the algorithm’s performance.
Moreover, the inclusion of problem-specific update A. Maximum independent set
rules comes with additional benefits. It allows us to
The maximum independent set (MIS) problem is a
broaden the scope of our algorithm beyond gate-based
paradigmatic NP-hard combinatorial optimization prob-
architectures to analog devices, as shown in our experi-
lem. It is of commercial relevance in areas such as network
ments with the QuEra Aquila neutral atom analog quan-
design or traffic optimization [12]. For a graph G = (V, E)
tum processor. This approach is reminiscent of previous
with vertex set V and an edge set E, the MIS problem
approaches leveraging spin-freezing schemes on quantum
amounts to finding the largest subset S ⊆ V such that no
annealers and in various Monte Carlo methods [9]. Fur-
two vertices in S are adjacent [see Fig. 5(a) below for an
thermore, for problems with hard constraints, the update
example]. For a binary vector x ∈ {0, 1}|V| , finding the MIS
rules can offer options to enforce feasibility by design,
is equivalent to finding the ground state of the (classical)
even in the presence of noise.
Hamiltonian
Finally, we propose backtracking as a strategy to further
improve the performance of QIRO, by attempting to iden-  
tify and correct nonideal decisions made at earlier stages H (x) = − xi + λ xi xj . (1)
of the algorithm. Backtracking provides a way to enhance i∈V (i,j )∈E
algorithmic performance without necessitating an increase
in circuit depth, as is commonplace in quantum optimiza- Here, λ > 1 is a penalty term to enforce the independence
tion algorithms. Taken together, the robustness of QIRO, constraint [10]. However, one should be mindful of the
its applicability to analog devices, and the possibility to pitfalls of enforcing the independence constraint as a soft
improve its efficacy without demanding better quantum penalty term, because low-energy states (other than the
resources make QIRO a promising NISQ algorithm. ground states) of H (x) might not necessarily obey the inde-
In this work we provide concrete implementations pendence constraint. Performing the mapping zi = 2xi − 1
for two paradigmatic NP-hard combinatorial optimiza- to Ising variables zi = ±1, and promoting the spin vari-
tion problems, namely, the maximum independent set ables to quantum operators zi → Ẑi translates Eq. (1) to a
and maximum satisfiability. However, our work should quantum Hamiltonian suitable for quantum devices [10].
be seen mainly as a proof-of-concept demonstration of a Moreover, the MIS problem on a family of graphs called
general template for designing hybrid quantum-classical unit disk graphs (UDGs) has been shown to have a natural
algorithms for combinatorial optimization. We study the mapping to neutral atom quantum devices [13,14]. There,

020327-2
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

the Rydberg blockade mechanism—which prevents simul- (classical) Hamiltonian

taneous excitations of nearby atoms—enforces the inde-
pendence constraint. Finding the MIS of a general UDG, 
m  Kj
m 
despite being a specialized family of graphs, remains an H (φ) = Cj = (1 − jk ). (2)
NP-hard problem [13] and has commercial applications j =1 j =1 k=1
(e.g., in network design [15]). For additional details, we
refer the reader to Appendix D. It is important to keep in mind that the jk are nothing but
the Boolean variables xjk (or their negations x̄jk ). Hence,
B. Maximum satisfiability by using the transformation into Ising variables outlined in
the case of the MIS problem in Sec. II A, we can likewise
The satisfiability (SAT) problem is arguably one of the
encode the MAX-2-SAT problem onto quantum devices.
best-studied combinatorial optimization problems and was
the first problem proven to be NP-complete [16]. Given a
propositional logic formula φ of n Boolean variables xi ∈ III. ALGORITHM
{0, 1}, the SAT problem is a decision problem that asks We now provide details about QIRO. Essentially, in
for the existence of an assignment of the variables xi that QIRO quantum resources are utilized to inform update
satisfies φ. Each clause can then be written as rules that recursively simplify the problem at hand (see
Fig. 1). The idea is a generalization of previous work,
Kj
 in particular RQAOA [3,7] (also described in Appendix
C̃j = jk , B). However, in contrast to RQAOA, the update steps
k=1 used in QIRO are problem specific. This comes with sev-
eral advantages and facilitates cross-pollination with ideas
where literals jk are either a variable xjk or its negation x̄jk from classical optimization [8].
for jk ∈ {1, . . . , n}, and Kj is the length of the j th clause. At each iteration of the QIRO algorithm, for a given
For a fixed Kj = K, the corresponding problem of whether problem Hamiltonian, we first prepare a low-energy quan-
or not the Boolean formula with m clauses tum state, generically in the form of a superposition of
low-energy candidate solutions to the optimization prob-

m lem. Next, we use this quantum state to extract information
φ= C̃j (e.g., correlations between variables), which is then passed
j =1 to a classical update step with the goal of simplifying
the problem. We design this update step to complement
has a satisfying assignment is referred to as the K-SAT the quantum part of the algorithm by performing opera-
problem. tions better suited to classical hardware. This may include
Here, we focus on the optimization version of satisfia- nonlocal steps that address the identified limitations of
bility, referred to as maximum satisfiability (MAX SAT). (local) quantum algorithms, or enforcing hard constraints.
In MAX SAT, the goal is to find an assignment of the vari- As such, the classical update step (and the specific design
ables xi that minimizes the number of violated clauses. The thereof) is crucial to the performance of QIRO, and is
MAX-K-SAT problem deals with propositional logic for- typically problem specific. As such, it might need to be
mulae with at most K literals per clause. An important modified by the end user, depending on their needs. How-
parameter of a (MAX-)K-SAT problem is the clause-to- ever, templates for the update steps may be borrowed from
variable ratio α := m/n, which strongly influences the the literature on classical combinatorial optimization [25].
algorithmic hardness of random problem instances [17]. The procedure of quantum state preparation and the sub-
Research indicates that there is a phase transition in com- sequent classical update step is repeated until the size of
putational complexity [18]. While for the K-SAT problem the problem is reduced sufficiently, such that, for exam-
an easy-hard-easy transition is observed [19,20], MAX-K- ple, we are able to solve it exactly by a brute-force search
SAT exhibits a transition from easy underconstrained to or other means if the reduction does not yield problem
hard overconstrained problems at a certain αc [21]. In what sizes amenable to exact solvers. Finally, we propose the
follows we limit ourselves to K = 2 (clauses contain either use of backtracking to identify and rectify problem reduc-
one or two literals), as MAX-2-SAT is already NP-hard tions that led to an unfavorable outcome, thus providing
[22], despite its decision version (2-SAT) being solvable in an additional strategy to obtain enhanced solutions without
linear time [23]. Here, the phase transition occurs at αc = 1 requiring deeper quantum circuits.
[24]. In what follows, we provide details on the individual
We use the fact that, for Boolean (binary) variables x̄i ≡ components of QIRO. We begin by describing the quantum
(1 − xi ), solving the MAX-2-SAT problem for a given state preparation methods utilized in this work, followed
formula φ is equivalent to finding a ground state of the by the proposed update rules for the MIS and MAX-2-SAT

020327-3
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

problems. Next, we introduce the general strategy of back- and defer the details pertaining to our implementation of
tracking and describe the particular implementation used the MIS problem on neutral atom quantum processors to
in this work. Finally, we combine the building blocks to Appendix D.
present the complete QIRO algorithm. At the beginning of an adiabatic protocol, we start in
an easy-to-prepare ground state of a simple Hamiltonian,
A. Quantum state preparation e.g., Ĥmix . We then slowly transform the Hamiltonian from
It is desirable that the quantum state that is prepared at Ĥmix to Ĥc . If the rate of change is slow enough com-
the beginning of each iteration has a low expectation value pared to the smallest instantaneous spectral gap (between
with respect to the cost Hamiltonian Ĥc . Therefore, we the instantaneous ground and first excited states), the adia-
hope that the information guiding the classical reduction batic theorem guarantees that at the end of the protocol,
steps is extracted from (a superposition of) configurations the ground state of Ĥc is obtained [26,27]. Inspired by
corresponding to good candidate solutions to the opti- adiabatic protocols, the quantum annealing computational
mization problem. A plethora of protocols for preparing paradigm has been the subject of widespread interest in
low-energy states may be found in the literature; in this recent decades [29]. In quantum annealing one forgoes the
manuscript we consider QAOA [1] and adiabatic state strict requirement of adiabaticity and, in general, uses ana-
preparation [26,27], as outlined next. log protocols to prepare near-ground states. These can, in
turn, be used to inform the classical update steps.
1. QAOA
Let us first introduce QAOA. For both the MAX-2- B. Update rules
SAT and MIS problems considered in this manuscript, We now discuss example update rules for the considered
the cost Hamiltonian
 can be
 written in the general (Ising) optimization problems. We begin by preparing a low-
form Ĥc = i,j Jij Ẑi Ẑj + i hi Ẑi , where we have omitted energy quantum state |ψ from which information is dis-
any constant terms. Using the cost Hamiltonian Ĥc from tilled. In QIRO this information comes in the form of one-
 II, and defining the mixer Hamiltonian as Ĥmix =
Sec. point ψ| Ẑi |ψ and two-point correlations ψ| Ẑi Ẑj |ψ,
− ni=1 X̂i , the QAOA circuit is given by an alternating which we store in matrix M (see Algorithm 1 below).
application of the time evolution generated by these two Potential generalizations, such as using expectation values
operators. Specifically, taking the initial state as |ψ0  := of other observables, are left for future research.
|+⊗n , the state produced with QAOA at depth p is If the lowest depth (p = 1) of QAOA is used to gen-
erate these correlations, the only nonzero values of M in
|ψ(β, γ ) = e−iβp Ĥmix e−iγp Ĥc · · · e−iβ1 Ĥmix e−iγ1 Ĥc |ψ0  . an ideal numerical simulation as described in Appendix
A (i.e., free of shot and environmental noise) correspond
The 2p parameters (β, γ ) are then determined by a classi- to pairs of variables that are connected in the original
cal optimization routine such that ψ(β, γ )| Ĥc |ψ(β, γ ) interaction graph, i.e., to pairs of variables with nonzero
is minimized. coupling coefficients Jij . Therefore, we design our update
In this work, we predominantly focus on the lowest rules to only consider correlations between such connected
QAOA depth (p = 1), with the goal to establish a lower variables. Inspired by RQAOA [3,7], we next use these
bound for the performance of QIRO. It is natural to expect correlations to successively simplify the problem with each
that the improved quality of the quantum information update step. However, in contrast to RQAOA, the QIRO
extracted from higher-depth (p > 1) QAOA circuits (in update steps are problem specific in that they exploit the
the absence of noise) would enhance the performance of structure of the given optimization problem. As such, at
the algorithm. Our p > 1 experiments on small problem any stage of the algorithm, the problem that remains to be
instances shown in Appendix E corroborate this intuition. solved is a valid instance of the same class of optimization
For the state produced by QAOA at p = 1, the expec- problems. This will be crucial when deploying QIRO to
tation values of the single- and two-body terms compris- solve the MIS problem on neutral atom quantum hardware.
ing Ĥc can be calculated classically in polynomial time, In the two specific cases detailed in the remainder of this
enabling us to perform numerical simulations of systems section, only the largest correlation in terms of its absolute
with several hundred variables [3,7,28]. We note that there value is used to inform the update step, with ties broken
is no method known (to us) how this could be done effi- at random. Intuitively, the largest correlation corresponds
ciently for p > 1. For details about our implementation, to a variable (pair) that appears in a certain configura-
we refer the reader to Appendix A. tion with high probability in the low-energy state |ψ.
Thus, it is reasonable to use this information to inform
2. Adiabatic state preparation the update step. For example, if a large positive correla-
We next turn our attention to adiabatic state prepa- tion Ẑi Ẑj  ≈ 1 is measured, it is likely that most of the
ration. We present a general description of the protocol low-energy configurations have xi and xj take the same

020327-4
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

ALGORITHM 1. Reduce. (a) (b)

Input: Problem P , partial solution S. +1 –1
Output Simplified prob. P  , updated partial solution S  .
1: Prepare low-energy quantum state |ψ.
2: Store correlations in M ∈ Rn×n : (d)
(c)
∀i ∈ [n] : Mii = ψ| Ẑi |ψ, +1 –1
∀(i, j) ∈ [n] × [n] s.t. Jij = 0 : Mij = ψ| Ẑi Ẑj |ψ.
3: P  , S  ← Simplification (P, S, M)
 Problem-specific simplification is informed by the cor-
relations stored in M .
4: return P  , S  FIG. 2. MIS update rules described in the main text carried
out on an example graph. Panels (a) and (b) show reduction
rules for the case of a positive and negative one-point correlation,
respectively. Analogously, panels (c) and (d) show the respective
value. Generalizations of this deterministic scheme to ran- reduction rules for the case of a positive and negative two-point
domized (e.g., evolutionary) schemes will be studied in the correlation.
future.
As we lay out next, these update steps are carried out
in the original formulation of the respective optimization j simultaneously, i.e., we remove every node k for
problems, facilitating the adoption of techniques from the which (k, i) ∈ E and (k, j ) ∈ E. Intuitively, we know
extensive literature on classical combinatorial optimiza- that if the variables are negatively correlated, one of
tion. We call this update step Reduce and describe it in them will be assigned to be in the IS. Thus, we only
Algorithm 1. remove variables that are connected to both of them
The problem-specific updates employed for the opti- simultaneously.
mization problems considered in this work constitute the
If no nodes were removed from the graph using the
Simplification function in Algorithm 1, and are
described simplifications [only possible in case (d)], we
described next. However, we emphasize that the design
repeat the procedure using the next largest correlation in
choices we make are by no means the only option, with
terms of its absolute value. Finally, we identify all con-
many modifications and extensions possible.
nected components (i.e., connected subgraphs that are not
part of any other connected subgraph) of the graph. If the
1. MIS simplification connected component contains less than a certain number
At the beginning of the Simplification routine of vertices nc (for concreteness, we use nc = 15 here), we
for MIS, we find the entry of M with the largest abso- compute its maximum independent set by a brute-force
lute value. Depending on the sign of the correlation with search, and remove it from the graph. This step is mostly
the largest absolute value and on whether it lies on the performed in an effort to save resources when the algorithm
diagonal (corresponding to a one-point correlation) or not is deployed on real quantum hardware.
(corresponding to a two-point correlation), we perform a Crucially, the rules delineated above explicitly ensure
different simplification. The simplifications are designed that the independence constraint is obeyed. Hence,
such that the correlation-informed reductions are consis- the solutions obtained by repeated applications of
tent with the independence constraint. The four cases are Simplification are definitely feasible. This holds
visualized in Fig. 2 and described next. irrespective of the quality of the (quantum) information
used to perform the simplification, thus bestowing the
(a) If Mii ≥ 0 was selected, we set the ith vertex to scheme with the resilience required in the NISQ era. How-
be in the independent set (IS). We then remove all ever, although the update rules improve the performance of
vertices connected to the ith node from the graph, the algorithm by ensuring solution feasibility, the quality
as including them would violate the independence of the obtained solutions ultimately hinges on the quality
constraint. of the supplied correlations.
(b) If Mii < 0 was selected, we remove the ith vertex Interestingly, in contrast to the RQAOA update rules the
from the graph. QIRO-MIS update rules proposed here are local according
(c) If Mij > 0 [for (i, j ) ∈ E] was selected, we remove to the definition of locality in Sec. I, as no new connec-
both nodes from the graph—this is because the only tions between variables are created. However, the RQAOA
positively correlated assignment of connected ver- update rules lead to intermediate Hamiltonians that gener-
tices consistent with the independence constraint is ically do not take the form of an MIS Hamiltonian as in
not to include either node in the IS. Eq. (1). Conversely, at each stage of the QIRO algorithm
(d) If Mij < 0 [for (i, j ) ∈ E] was selected, we remove the reduced problem is a valid MIS problem (see also Fig.
nodes that are connected to both node i and node 2). Thus, if one can encode the original MIS problem onto

020327-5
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

an analog quantum device (e.g., a neutral atom quantum ALGORITHM 2. QIRO.

processor), the same can be done with intermediate, sim-
Input: Problem P .
plified problems. As such, one might be able to harness the
Output Complete solution S.
nonlocal effects generated by the quantum many-body time
1: Initialize empty solution S.
evolution on a neutral atom quantum processor, circum- 2: while size(P ) > 0 do
venting the issues of local quantum algorithms. Previously, 3: P, S ← Reduce (P, S)  See Algorithm 1.
such adiabatic protocols have been shown to generate solu- 4: end while
tions surpassing those obtained by running QAOA, even at 5: return S  Return the solution.
depths higher than p = 1 used in this manuscript [14].

2. MAX-2-SAT simplification the MAX-2-SAT problem), the solution of this simplified

problem is found by a brute-force search of the remaining
In the MAX-2-SAT case the Simplification step
solution space.
likewise begins by identifying the entry of M with the
largest absolute value. Depending on whether the largest
value corresponds to a one-point or a two-point correlation, C. Quantum-informed recursive optimization
we apply a different reduction. algorithm
We now possess all of the ingredients for the QIRO
(a) If Mii was selected, we set the variable xi → algorithm. At each step of QIRO, we begin by prepar-
[sign(Mii ) + 1]/2, i.e., a positive (negative) one- ing a quantum state (as described in Sec. III A). We
point correlation corresponds to assigning xi to proceed by using the correlation information by follow-
TRUE (FALSE). ing the problem-specific prescriptions from Sec. III B.
(b) If Mij (for i < j ) was selected, we replace xi → xj The simple end-to-end QIRO procedure is summarized in
if sign(Mij ) = 1, and xi → x̄j if sign(Mij ) = −1. Algorithm 2.
Intuitively, this update step captures the relationship
between variables xi and xj as inferred from the cor-
D. Backtracking
relation. If there is a positive (negative) correlation,
we assign xi and xj the same (opposite) values. We now turn our attention to backtracking. QIRO—as
specified in Algorithm 2—is a polynomial-time algorithm.
We note that only the two-point correlation reduction Therefore, under common complexity theoretical assump-
step described in (b) constitutes a nonlocal update because tions, there exist instances of NP-hard problems that QIRO
it potentially introduces a new link between previously will not be able to solve to optimality. This means that
unconnected variables. Because we only consider formulae some of the reduction steps performed by QIRO to arrive
in conjunctive normal form (i.e., conjunctions of disjunc- at a candidate solution might be wrong. A similar obser-
tions), setting a certain variable to the value required by vation was made for RQAOA in Ref. [33], prompting
the clause is sufficient to satisfy the clause (see the exam- the authors to consider reinforcement learning to steer
ple update step in Fig. 1). Hence, we can remove the clause RQAOA update steps. Here, we propose backtracking to
from consideration, thus simplifying the problem. Specif- derive improved solutions starting from an initial candi-
ically, all clauses in which a given literal evaluates to date solution obtained by QIRO. Through backtracking,
TRUE can be removed. Conversely, if the literal evaluates we obtain a method to refine solutions by using the same
to FALSE, only the literal itself can be removed from the quantum resources more frequently. This approach stands
clause. If a certain clause becomes empty (i.e., if both liter- in stark contrast to most other quantum optimization algo-
als were removed from it), it corresponds to an unsatisfied rithms, which usually require enhanced quantum resources
clause. (e.g., increasing the depth p in QAOA) to deliver improved
Next, we apply inference rules, which is a standard tech- solutions.
nique in the field of satisfiability solvers [30,31]. Inference In this work, we only employ this strategy for the
rules deduce information from the structure of the cur- MAX-2-SAT problem. However, extending this approach
rent Boolean formula, simplifying it by assigning values to other optimization problems such as the MIS problem is
to certain variables and thus speeding up further computa- straightforward.
tion. The rules guarantee that the optimal solutions of the As visualized in Fig. 3 the simplification process of
current and simplified problems are of equal quality. The QIRO can be represented by means of a search tree. Deci-
chosen inference rules (e.g., the pure literal rule [32]) are sions (i.e., update steps) made by the algorithm correspond
conceptually simple rules that incur a small computational to arrows in the search tree, which lead to simplified prob-
overhead. Further details on the specific inference rules we lems represented by nodes. Following the update rules in
employ are provided in Appendix C. Finally, if the number Sec. III B (i.e., successively applying Reduce) traces out
of remaining variables falls below nc (we use nc = 10 for a path through the search tree (see the red arrows and

020327-6
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

ALGORITHM 3. QIRO + backtracking (QIRO + BT).

Input: Problem P .
Output Complete solution S.
1: Initialize empty solution S.
2: L ← {}  List to store partial problems & solutions.
3: while size(P ) > 0 do
4: Append (P, S) to L.
5: P, S ← Reduce (P, S )  See Algorithm 1.
6: end while
7: for P  , S  in L do   Backtracking.

8: P  , S  ← Reduce P  , S  ; Solution(0···00)
 Make alternative decision.
9: while size(P  ) > 0 do
10: P  , S  ← Reduce (P  , S  )
11: end while
12: Store the obtained solution.
13: end for
14: return Best solution found.

simply perform the opposite reduction to that we per-

formed initially. For example, if at first the decision was
made to assign the value of a variable xi = TRUE, we
FIG. 3. Schematic visualization of the backtracking procedure reverse that decision and set xi = FALSE instead. Thus,
used to obtain improved solutions. Starting from the original we diverge from the original path through the search tree as
problem P (0) and an empty solution S (0) , we perform a series visualized in Fig. 3. We note that Reduce requires infor-
of reductions (red arrows and nodes) that recursively simplify mation about the solution it is deviating from in order to be
the problem, until the first candidate solution (Solution(0···00) ) able to perform the reversed decision.
is reached. This corresponds to applying QIRO as outlined in In principle, exploring the entire search tree essentially
Algorithm 2. Backtracking, as visualized by the dashed teal corresponds to a brute-force search of the complete solu-
arrows, then amounts to revisiting reduction steps made along
tion space. Because we are interested in designing heuristic
the path. At the chosen ancestor node we have backtracked to, we
can then make an alternative reduction, here denoted by Reduce algorithms that run in polynomial time, we limit the per-
(see the black arrows). Concretely, Reduce amounts to revers- mitted exploration. We find that a simple rule for which the
ing the initially made decision, not needing extra calls to the number of QAOA circuit executions scales quadratically
quantum device. In this work, we limit ourselves to making a in the problem size is sufficient to produce satisfactory
binary decision at each node. As such, we can label the reduced results for the considered problem sizes of MAX-2-SAT.
problem and partial solution at each node with a binary string However, more (or less) sophisticated approaches can be
encoding the path connecting the node to the original problem designed and tailored to the particular needs and resources
node. Each 0 (1) in the superscript labeling P and S corre- of the end user. Next, we describe the specific backtrack-
sponds to one application of Reduce (Reduce) starting from
ing strategy used in our simulations and refer the reader to
the original problem.
Algorithm 3 for the pseudocode.
First, we generate the initial candidate solution
Solution(0···00) by successively applying Reduce, as visu-
nodes in Fig. 3), which eventually leads to a candidate alized in Fig. 3. Put differently, we first apply QIRO
solution (see Solution(0···00) in Fig. 3). Here, the superscript without backtracking. Next, we backtrack to each of the
(0 · · · 00) indicates that this is the solution obtained in the ancestor nodes of the initial candidate solution (as visu-
initial passage through the search tree. alized by the dashed teal lines in Fig. 3). Each of the
We would like to consider alternative solutions, obtain- revisited parent nodes corresponds to a different simplified
able by deviating from the initial path through the search problem and a corresponding partial solution. The decision
tree that led to the initial candidate solution Solution(0···00) . made in the initial passage through the search tree is then
To this end, we define the Reduce function, which reversed (by applying Reduce), whereby a new node in
reverses the decision previously made at the revisited node the search tree (e.g., P (01) , S (01) in Fig. 3) is created. Here,
(i.e., partial problem and solution). We emphasize that this each 0 (1) in the superscript labeling P and S corresponds
initial reversal of the previously made decision requires no to one application of Reduce (Reduce) starting from
additional measurements from the quantum device—we the original problem. Subsequently, if the problem has not

020327-7
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

been fully reduced yet, from this new node a new candi- algorithm, which iteratively assigns a node with the min-
date solution is generated by repeatedly applying Reduce. imal degree to the independent set and likewise removes
In contrast to the initial application of Reduce the sub- adjacent nodes to ensure the validity of the final solution
sequent applications of Reduce require additional calls [42]. We note that the minimal degree greedy algorithm
to the quantum device—see the definition of Reduce in was designed chiefly for sparse graphs (i.e., of low average
Algorithm 1. In other words, each revisited node leads to a degree). This is indicated by the performance guarantee of
new candidate solution that is identical to the initial candi- the greedy algorithm, which is inversely proportional to the
date solution in the values of the variables fixed before the graph’s average degree [42]. However, we note that, even
revisited node. From the revisited node onward, the new for the sparse graphs considered throughout this section,
candidate solution differs from the initial solution. Finally, stronger algorithms exist in the literature—see Ref. [43]
the output of the algorithm is the best candidate solution for benchmarks.
found in this way. We note that the order in which the Importantly, because all of the algorithms considered
parent nodes are revisited is not important. here operate by iteratively simplifying the problem, their
performance is indicative of the fitness of the information
used to perform these simplifications.
IV. RESULTS For each system size and average degree, we generated a
Next, we present the results obtained by applying set of 50 random Erdős-Rényi graphs and solved these with
QIRO to random instances of the MIS and MAX-2-SAT the ReduMis heuristic solver from the KaMIS project
problems. In numerical simulations we compare the per- [44,45]. We use the solution obtained by ReduMis as (a
formance of quantum algorithms (QIRO, RQAOA) and proxy for) the optimal solution, and report the percent-
two prominent classical heuristics (simulated annealing, age of graphs for which each of the algorithms considered
greedy algorithms). Of course, more complex tailored performs on par with the ReduMis reference.
classical approaches can be used, such as the parallel The results for graph sizes with n = 40–200 nodes and
tempering [34–36], borealis [37], and extremal optimiza- average degrees of 3, 5, and 12 are shown in Fig. 4. We
tion [38] heuristics, as well as a plethora of exact solvers used λ = 1.1 in the MIS Hamiltonian [see Eq. (1)] for
benchmarked at the annual (MAX-)SAT competition [39]. RQAOA and QIRO. The performance of the two quantum
Efficient classical simulations of p = 1 QAOA allow us algorithms (QIRO and RQAOA) is found to lie between
to study the system size scaling on problem instances with the random and minimal degree greedy algorithms. This
various connectivities of the underlying interaction graph. highlights the limitations of the quantum correlations
Results using the lowest depth of QAOA should mainly be obtained from QAOA at p = 1, which only “sees” the
interpreted as a lower bound on the performance of QIRO. immediate neighbors of every node. As such, it is unsur-
Furthermore, we investigate the robustness of QIRO and prising that the minimal degree greedy algorithm—which
RQAOA against perturbations of the optimal parameters of makes the optimal greedy move based on information
the QAOA circuits generating the correlations. This allows from nearest neighbors—outperforms QIRO and RQAOA
us to emulate failures of the classical optimizer to find the informed by p = 1 QAOA. However, our results in
optimal parameters. Notably, we also deploy QIRO on the Appendix E indicate that (at least on small graphs), QIRO
QuEra Aquila neutral atom device (as available on Ama- with QAOA at p > 1 outperforms the greedy benchmark.
zon Braket) to (approximately) solve the MIS problem on Furthermore, RQAOA (slightly, but consistently) out-
unit disk graphs with more than a hundred nodes. performs QIRO. This is likely due to the fact that the QIRO
simplification rules as defined in Sec. III B 1 are manifestly
local. The fact that the gap between QIRO and RQAOA
A. MIS seems to diminish for denser graphs (i.e., for increas-
We begin by presenting the results of QIRO and ing average degree) corroborates this intuitive explanation,
RQAOA as applied to the MIS problem. First, we consider as we expect locality to be a weaker limitation if the
the algorithmic performance as a function of the system underlying interaction graph has a higher connectivity.
size on the well-studied family of Erdős-Rényi graphs, A major advantage of QIRO (as compared to, e.g.,
with fixed average degrees [40]. As a classical comparison, RQAOA) is that the problem-specific classical subroutines
we use two variants of a greedy algorithm. This choice is can be used to enforce feasibility of the solutions generated
motivated by recent findings that information from quan- (i.e., all obtained candidate solutions are guaranteed to be
tum devices can be used to enhance the performance of independent sets). To showcase this feature, we performed
greedy solvers [41]. First, we consider a simple random numerical simulations of QIRO and RQAOA with p = 1
greedy algorithm, which iteratively assigns random nodes QAOA with suboptimal parameters. The goal of these
to the independent set, and removes nodes adjacent to the experiments is to consider the robustness of the algorithms
selected node such that no violations are possible. As the against lower-quality quantum information. The parame-
second comparison, we use the minimal degree greedy ters are generated by sampling a uniform 30 × 30 grid in

020327-8
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

100 of valid solutions for RQAOA. Because finding optimal

Solved instances (%)
degree = 3
parameters in an actual experimental setting is challeng-
ing [46,47], the robustness of QIRO is a valuable feature,
50 in particular for realistic experimental scenarios. Addition-
ally, results in Fig. 5(b) suggest that increasing the quality
of the quantum correlations has a positive impact on the
0 performance of both QIRO as well as RQAOA. This is in
agreement with the outcome of our experiments at p > 1
100
Solved instances (%)

degree = 5 in Appendix E, where a clear improved QIRO performance

with increasing p is observed.
Surprisingly, results shown in Fig. 5(b) indicate that the
50
number of valid solutions produced by RQAOA falls off
at a comparatively better parameter quality (lower q) for
the larger value of λ, defying the expectations that a larger
0
penalty term λ should favor valid solutions. In general, we
100 observe complex behaviour in the (λ, q) parameter plane
Solved instances (%)

degree = 12
(see Fig. 10 and the associated discussion in Appendix F).
Min. greedy RQAOA We note that determining the magnitude of penalty terms
Rnd. greedy QIRO
50 is a well-known issue when incorporating hard constraints
into an unconstrained optimization setup [10,11].
Importantly, at each step of the QIRO algorithm, what
0 remains to be solved is a valid MIS problem. Moreover,
40 60 80 100 120 140 160 180 200 if the starting graph is a UDG [as in Fig. 5(a)] then so
Number of vertices are all successive (simplified) graphs—see also Fig. 2.
Because UDGs can be naturally embedded on neutral atom
FIG. 4. System size scaling for the success rate of the ran- quantum processors [13,14], we can use such devices to
dom greedy (“Rnd. greedy”), minimal degree greedy (“Min. produce the required quantum correlations. In what fol-
greedy”), RQAOA, and QIRO for Erdős-Rényi graphs with dif- lows, we present results from the QuEra Aquila neutral
ferent expectation values of the average degree. At each system
atom quantum processor, as accessed through Amazon
size and average degree 50 graphs were generated, and ten runs
of each algorithm were performed. We report the median per-
Braket. The details on the (quasi)adiabatic protocols used
centage of graphs for which we find independent sets of the size to generate the correlations on the QuEra Aquila device
found by the classical ReduMis benchmark algorithm. The error can be found in Appendix D. Because of the limited acces-
bars denote the minimum and maximum performance observed sibility of the device, we report the results on a per-instance
across the ten runs. Note that the results for QIRO and RQAOA basis. We sort the instances in increasing order of their
are obtained using simulations of QAOA at depth p = 1. classical hardness [see Eq. (D1)], as proposed by Ebadi
et al. [14]. We provide details about the hardness of the
considered graph instances in Appendix D.
the (β, γ ) parameter plane, and choosing the parameters In Fig. 6 we display the approximation ratios |IS|/|MIS|
corresponding to the qth quantile in the energy values, obtained by QIRO with correlations from the QuEra device
among the 30 × 30 energies obtained. Given their rele- (QIRO QuEra) and classical simulations of p = 1 QAOA
vance for neutral atom quantum devices [13,14], we use (QIRO QAOA). On the QuEra device, we generate the cor-
UDGs (see Appendix D) for these experiments. relations using the protocols outlined in Appendix D. We
The results for an ensemble of 50 UDGs with 137 compare these to the sizes of the independent sets found by
nodes and two different values of the penalty term λ are the minimal degree greedy algorithm. The reference maxi-
shown in Fig. 5(b). For every graph and parameter qual- mum independent set is computed using the tensor network
ity q, ten runs of each algorithm were performed. We algorithm from Ref. [48]. For clarity, results of RQAOA
show the mean approximation ratio of valid solutions pro- are omitted from the plot, because the approximation ratios
duced by RQAOA and QIRO as plotted against the quality found for UDGs were generally comparable to those found
quantile q of the parameters used to generate the correla- by QIRO QAOA [as seen in Fig. 5(b)].
tions by the QAOA. The parameters at q = 0 correspond Our results indicate that the quantum correlations from
to the fully optimized parameters. As expected, QIRO the QuEra machine are better at guiding the QIRO proce-
only produces valid solutions. Remarkably, if perfect (i.e., dure than numerical simulations of QAOA at p = 1. This
optimized) parameters are supplied to RQAOA, it like- is in line with previous results for solving the MIS problem
wise produces exclusively valid solutions. However, for using analog protocols on such hardware [14,49]. A larger
both values of λ, there is a sharp drop-off in the number sample size would be required to solidify such claims.

020327-9
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

(a) (b) λ = 1.1 λ = 1.5

1.0

Validity
0.5 RQAOA
QIRO
0.0

0.99

|IS|/|MIS|
0.9

0
0.00 0.02 0.04 0.06 0.08 0.10 0.00 0.02 0.04 0.06 0.08 0.10
QAOA parameter quality q QAOA parameter quality q

FIG. 5. (a) An example unit disk graph with 137 nodes. Teal nodes indicate a maximum independent set of size |MIS| = 45. (b) Data
for 500 random unit disk graphs with 137 nodes, comparing the robustness of RQAOA and QIRO, for two different magnitudes of the
penalty term λ [see Eq. (1)]. The upper row shows the fraction of valid solutions obtained by each algorithm. The lower row displays
the ratio between the size of the largest independent set (IS) found and the maximum independent set (MIS) size. Both the validity
and |IS|/MIS| are plotted with respect to the quality of the p = 1 QAOA parameters, with parameters at q = 0.0 corresponding to the
optimal parameters, and higher values of q corresponding to worse parameters (see the main text). We note that we only report the
approximation ratio achieved by RQAOA if the solution produced by it is valid (i.e., corresponds to an independent set).

However, initial results appear to suggest that QIRO instances with clauses of length two. We compare our
QuEra performs comparably to the minimal degree greedy results to those obtained with RQAOA and simulated
algorithm, reinforcing our belief that improved quantum annealing (SA) [50], an established and relatively sim-
correlations boost the performance of QIRO. Strikingly, ple classical heuristic. More detailed simulations with
the minimal degree greedy algorithm fails dramatically other, more powerful heuristics are left for future studies.
on some graphs (e.g., instance 1), where it achieves an Details about our implementation of SA can be found in
approximation ratio of around 0.6, while QIRO QuEra is Appendix G.
able to solve the same problems (nearly) optimally. This We first analyze the system size scaling of the algo-
preliminary evidence showcases the potential utility of rithms’ performance for random MAX-2-SAT instances
quantum correlations for guiding optimization algorithms. with up to 160 variables, at three different values of the
clause-to-variable ratio α. We use the MAX-SAT solver
RC2 [51] to compute the optimal solutions, and use paral-
B. MAX-2-SAT lel tempering (PT) [34] to compute (approximately) opti-
Next, we analyze the performance of QIRO and QIRO mal solutions for instances where the RC2 solver timed
with backtracking (QIRO + BT) on random MAX-2-SAT out. Details of our implementation of PT can be found

QIRO QuEra QIRO QAOA Min. greedy

1.000

0.975
|IS|/|MIS|

0.950

0.925

0.900
1 2 3 4 5 6 7 8 9 10
Graph instance

FIG. 6. Approximation ratio comparing the size of the independent set |IS| found by QIRO with correlations from the QuEra neutral
atom quantum processor (QIRO QuEra), and correlations from numerical simulations of p = 1 QAOA (QIRO QAOA), with the size
of the maximum independent set |MIS|. For comparison, the ratios obtained by the minimal degree greedy algorithm (“Min. greedy”)
are shown. In instance 1 the obtained approximation ratio of the greedy algorithm is outside of the range of the displayed values with
|IS|/|MIS| ≈ 0.63. The median value over five runs is plotted, with error bars showing the best and worst solutions found.

020327-10
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

RQAOA QIRO QIRO+BT SA than half of the instances for almost all system sizes and all
Solved instances (%) Solved instances (%) Solved instances (%) values of α considered here. For α = 2, it even slightly out-
100
performs our SA implementation. It should be noted that
while the number of calls to the quantum device increases
50 significantly (i.e., from linear to quadratic in the problem
size) if backtracking is included, the requirements on the
α=2 quantum hardware in terms of the gate and qubit count
0 are kept constant. This fundamentally differs from QAOA,
100
where an improved performance can only be obtained
by increasing the circuit depth p, and thus the quantum
resources required.
50

V. DISCUSSION
α=4
0 In this work we have extended previous efforts on hybrid
100 quantum-classical algorithms for combinatorial optimiza-
tion, by introducing a family of QIRO algorithms. We pro-
pose the use of quantum-informed classical update steps
50 to recursively simplify the original optimization problem,
allowing for the incorporation of problem-specific clas-
α=6 sical subroutines (e.g., satisfiability inference rules) from
0
the rich literature on classical combinatorial optimization
40 60 80 100 120 140 160
[25]. We have shown via numerical simulations for ran-
Number of variables
dom instances of two NP-hard combinatorial optimization
FIG. 7. System size scaling of the performance of the differ- problems that such classical subroutines can significantly
ent algorithms, at three different values of α. At each system size strengthen the algorithm’s performance. Additionally, we
and α, 50 random instances were generated with two variables compare the performance of QIRO to RQAOA, providing,
per clause, each variable and its polarity chosen uniformly at to the best of our knowledge, the first systematic numerical
random. Parallel tempering (PT) was used to find nearly opti- study of RQAOA on problems that are not Z2 symmetric.
mal solutions; we then report the median percentage of instances Importantly, we think that QIRO offers a promising
on which an algorithm performs at least as well as PT across ten approach for solving optimization problems with hard con-
independent runs. The error bars indicate the minimal and maxi- straints, which have previously proven challenging for
mal performance across these ten runs. Simulation details about
quantum algorithms [10,11,52]. Using the MIS problem as
SA and PT can be found in Appendix G.
a guiding example, we have shown how classical update
steps can be used to enforce hard constraints and only
in Appendix G. We then report the number of times the produce feasible solutions. This is true regardless of the
algorithms perform at least as well as PT or the RC2 solver. quality of the quantum resources, endowing QIRO with the
Figure 7 showcases the results of our experiments. robustness essential for algorithms within the NISQ era.
QIRO performs better than RQAOA for all clause-to- Interestingly, we observe that keeping the quantum
variable ratios α and over the entire range of the considered resources equal, RQAOA outperforms QIRO on the MIS
problem sizes. This shows that performance gains can be problem. We believe that this is a consequence of the
made by including very simple problem-specific update local nature of the QIRO update rules for MIS. This is
rules that require only modest classical resources. We note, corroborated by our finding that improved quantum cor-
however, that the improvement is especially prominent relations from the neutral atom device enhance QIRO’s
for lower values of α, and diminishes for problems with performance compared to when QAOA at p = 1 is used
higher α (connectivity). Moreover, due to the inclusion of to generate the correlations. The intuitive explanation is
inference rules, the Boolean formula is simplified in fewer that the quantum many-body evolution on a neutral atom
update steps, thus requiring fewer calls to the quantum quantum processor provides sufficient nonlocal effects to
device (and fewer optimizations of the QAOA circuit) than counteract the locality of the update steps. In fact, adiabatic
in RQAOA. protocols analogous to those used to generate correlations
As expected, we note a decrease in performance with here have previously been shown to produce candidate
increasing problem sizes, with both QIRO and RQAOA solutions surpassing those produced by QAOA at depths
often struggling to find optimal solutions at the largest higher than p = 1 [14]. Taken together, our results may
problem sizes considered. However, with the addition of serve as yet another indication that locality is a limiting
backtracking, QIRO + BT is able to optimally solve more factor for (quantum) optimization algorithms [53].

020327-11
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

Furthermore, we have confirmed the intuitive expec- just the maximal correlation. This could be accomplished
tation that higher-quality quantum correlations should by a randomized strategy, with the selection probabilities
enhance the performance of the algorithm, as evidenced determined as some function of the correlation matrix.
by means of extensive numerical simulations, as well as Furthermore, devising update rules incorporating corre-
experiments on an actual neutral atom quantum device. lations beyond those between nearest neighbors in the
This result suggests that advances in quantum hardware interaction graph might be useful, especially in combi-
should lead to further performance gains. nation with deeper QAOA circuits or on analog devices
We note that we performed no optimization of the pro- with long-range interactions (e.g., Rydberg atom arrays).
tocols used to generate the correlations in our experiments Additionally, one could exploit the locality of QAOA, and
on the QuEra Aquila neutral atom device. Previously, simultaneously round multiple correlations at each update
dramatic improvements in the quality of candidate solu- step, if the corresponding variables are outside of each
tions have been obtained by first optimizing the schedules other’s sphere of influence. Thus, one could save resources
[14,49]. Moreover, the (unoptimized) protocols and device by requiring fewer calls to the quantum routines, and thus
used in this study (and specified in Appendix D) have mitigate one of the limitations of the presented scheme. On
been shown to produce candidate solutions typically corre- the other hand, if deeper quantum circuits (e.g., QAOA at
sponding to either very small independent sets, or possess- p > 1) are available, extracting correlations beyond those
ing multiple violations of the independence constraint [49]. between variables adjacent in the interaction graph might
Thus, it is reasonable to expect that such schedule opti- be beneficial. Finally, considering problem-tailored alter-
mization schemes could greatly benefit the quantum state natives to a brute-force search of the remaining solution
preparation step in QIRO. Moreover, the ability of QIRO space at a threshold problem size nc may significantly
to find good (feasible) solutions despite the limitations of improve the performance and resource efficiency of QIRO.
the quantum hardware serves as yet another demonstration It is unreasonable to expect that the simple backtrack-
of its robustness. ing strategy presented here would be equally successful for
In addition, we have introduced the use of backtracking larger problem instances, or instances from different prob-
to obtain improved results via additional calls to existing lem classes. Hence, one might need to devise improved
quantum hardware, rather than requiring deeper circuits, as backtracking strategies, both in terms of the resources
is usually the case in quantum (optimization) algorithms. required, as well as in terms of the performance. One
This scheme has led to promising results for the MAX-2- possibility is the use of reinforcement learning to guide
SAT problem, where QIRO enhanced with backtracking the reduction process, as proposed by Patel et al. [33].
could optimally solve most random problem instances, Resources might be saved by efficient pruning of the search
even at comparatively large instance sizes. Because large tree using existing methods [30], identifying and ignoring
instances are accessible in numerical simulations (at p = 1 unfavorable parts of the solution space. Therefore, fewer
of QAOA), we were able to sensibly compare the perfor- quantum state preparations would be required. Finally, one
mance of QIRO(+BT) with widely used classical heuristic could use backtracking to enforce constraints, by identify-
solvers such as simulated annealing. ing steps where a constraint was violated, backtracking to
Finally, we would like to stress that our work should it, and altering the decision made at that step.
be seen as a first demonstration of a general template for It would likewise be interesting to consider alterna-
a larger class of hybrid quantum-classical algorithms for tive methods of generating the quantum correlations. This
combinatorial optimization. We firmly believe that hybrid could include efficient numerical methods to obtain corre-
techniques will be necessary to unlock the advantages of lations with the QAOA at p > 1, by means of, e.g., tensor
quantum devices. The modular nature of QIRO offers a network simulations [54,55]. Naturally, it would be inter-
plethora of ways in which the presented toolbox can be esting to acquire correlations from higher-depth QAOA
modified; we touch upon some potential extensions in the circuits executed on real digital quantum computers, or
next section. by generating them using quantum annealers based on
superconducting flux qubits [56,57]. One could consider
optimizing analog protocols (using, e.g., Bayesian opti-
VI. OUTLOOK
mization [49]) to obtain better correlations, and enhance
Several modifications to the QIRO algorithm as pre- the effectiveness of the algorithm. Moreover, modify-
sented here can be envisioned. Exploring update rules ing the cost function with respect to which the quantum
for optimization problems beyond those examined in this resources are optimized might be helpful, similarly to the
manuscript would be a compelling avenue for future approach employed by Caha et al. [58]. Finally, one could
research. Likewise, it would be intriguing to devise attempt to generate correlations using classical approxima-
the update rules for weighted versions of the problems tion algorithms, as proposed by Wagner et al. [59]. Long
addressed here. Moreover, one could attempt to extend story short, there is a lot to be done, but this paper is already
the information extracted from the quantum state beyond long enough.

020327-12
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

Conveniently, the modular nature of QIRO allows us to all qubits involved are required. This makes QAOA expo-
analyze the promise of quantum computers for combinato- nentially difficult to simulate classically for any depth p
rial optimization, e.g., by comparing QIRO with different [61]. In contrast, the sampling step is replaced by a classi-
quantum and classical methods of computing the correla- cal update step in QIRO. Hence, using the depth-(p = 1)
tions. Should shortcomings of quantum approaches (e.g., QAOA as state preparation for QIRO allows us to simulate
locality) be identified in this manner, the classical sub- the algorithm efficiently.
routines in QIRO can be designed to rectify them. As One may wonder why this should be a quantum
such, we hope that our work can bring us closer to deter- algorithm of interest, if it is amenable to efficient classical
mining whether quantum devices can provide value for simulations. The answer to this lies in the fact that the situ-
combinatorial optimization problems. ation changes dramatically at higher depths of QAOA with
p > 1 where no method for analytically computing the
ACKNOWLEDGMENTS one-point and two-point correlations classically in poly-
J.R.F. and A.K. thank Libor Caha and Alexander Kli- nomial time is known (to us). The same holds at p = 1
esch for insightful discussions. The authors thank Lilly when calculating the expectation values of k-body interac-
Palackal, Maximilian Passek, Carlos Riofrío, and Gili tion terms, with k ≥ 3. Therefore, these cases could lead
Rosenberg for thorough reviews of the manuscript, and to an algorithm that provides quantum advantage. In gen-
the Amazon Braket, BMW, and QuEra teams for their eral, we expect an increasing performance of QIRO for
support. C.M. thanks the Munich Quantum Valley ini- higher depths of QAOA, as indicated by the MIS experi-
tiative, which is supported by the Bavarian State Gov- ments in Sec. IV A. We note that similar observations were
ernment with funds from the Hightech Agenda Bayern previously made by Bravyi et al. [3].
Plus. H.G.K. would like to thank Am Platzl 1A for pro-
viding the necessary environment for creative thinking. APPENDIX B: RQAOA
An open-source implementation of QIRO is available To overcome the limitations of the QAOA due to
online [60]. its local structure [3,53], Bravyi et al. [3,7] introduced
RQAOA. The algorithm optimizes a QAOA state based
APPENDIX A: LOW-DEPTH QAOA on the current problem Hamiltonian and uses this infor-
SIMULATIONS mation to reduce the problem via the elimination of a
single variable. Recursive application of this step reduces
Here, we provide details of our implementation of the
the problem until a specified threshold of remaining vari-
numerical simulations of QAOA at p = 1. If we limit our-
ables nc is reached. The remaining problem can be solved
selves to problems with general Ising-like (i.e., quadratic)
by a classical solver (e.g., by a brute-force search). The
cost Hamiltonians of the form
elimination step is motivated by the rounding of fractional
  correlations in relaxations of linear programs [3]. While
Ĥc = hi Ẑi + Jij Ẑi Ẑj , (A1) hardly a fair comparison in terms of resources required,
i i<j RQAOA has been shown to significantly outperform
QAOA [7,62].
it was shown that the expectation values of Ẑi and Ẑi Ẑj It was first applied to the Z2 -symmetric Max-Cut prob-
according to the depth-(p = 1) QAOA state |ψ(β1 , γ1 ) = lem, which exclusively comprises two-body terms. Here,
e−iβ1 Ĥmix e−iγ1 Ĥc |+⊗n can be calculated in O(n) time via we (trivially) extend the approach to also include single-
simple analytical expressions [3,28]. Because the cost body terms. We additionally note that a generalization to
Hamiltonian comprises O(n2 ) terms, the complexity of Hamiltonians with an arbitrary degree of interaction can be
simulating the calculation of the correlation matrix M in envisaged. However, for Hamiltonians with K-body terms
the update step of Algorithm 1 is O(n3 ). Therefore, assum- with K ≥ 3, one faces an additional obstacle, as the degree
ing O(n) update steps, this yields an overall complexity of the interactions can increase when an elimination step is
of O(n4 ) for the full QIRO simulation and O(n5 ) for the performed. This comes with additional requirements on the
QIRO + BT simulation, when O(n) backtracking steps are quantum hardware; as we are mostly interested in NISQ-
performed, as is the case following the prescription of era algorithms, we do not consider such extensions here.
Algorithm 3. Assuming a quadratic Hamiltonian of form (A1), the
We comment here on the difference in the measurement pseudocode for the RQAOA, as proposed in Refs. [3,7],
processes between QIRO and vanilla QAOA. If the cost is given in Algorithm 4.
Hamiltonian is k local, the process of finding the optimal
parameters according to the expectation value of the cost
APPENDIX C: MAX-2-SAT SIMPLIFICATION
Hamiltonian in principle requires only measurements with
DETAILS
k qubits involved for both algorithms. However, after the
optimal parameters are found, QAOA samples classical bit Here we provide the details on the update rules
strings from the quantum state for which measurements of employed in the QIRO update step for MAX-2-SAT.

020327-13
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

ALGORITHM 4. RQAOA. ALGORITHM 5. Reduce (MAX-2-SAT).

Input: Problem Hamiltonian Ĥ, threshold nc . Input: Formula φ, partial solution S.

Output Solution S. Output Simpl. formula φ, updated partial solution S.
1: L ← {}  Initialize list to store frozen correlations. 1: Prepare low-energy quantum state |ψ.
2: for k = 1 to n − nc do 2: Store correlations in M ∈ Rn×n :
3: Find optimal variational parameters of QAOA: ∀i ∈ [n] : Mii = ψ| Ẑi |ψ,
β ∗ , γ ∗ = arg minβ,γ ψ(β, γ)| Ĥ |ψ(β, γ) ∀(i, j) ∈ [n] × [n] s.t. Jij = 0 : Mij = ψ| Ẑi Ẑj |ψ.
4: |ψ ← |ψ(β ∗ , γ ∗ ) 3: if i = j then  One-point correlation.
5: Store correlations in M ∈ Rn×n , initialize Mij = 0: 4: Assign variable xi → 12 (sign(Mii ) + 1).
∀i ∈ [n] : Mii = ψ| Ẑi |ψ, 5: else  Two-point correlation.
∀(i, j) ∈ [n] × [n] s.t. Jij = 0 : Mij = ψ| Ẑi Ẑj |ψ. 6: if sign(Mij ) = 1 then
6: (i, j) = arg max(i,j)∈[n]×[n] |Mij | 7: Replace xi → xj in φ and extend S accordingly.
7: Append {(i, j), sign(Mij )} to L. 8: else
8: if i = j then  One-point correlation. 9: Replace xi → x̄j in φ and extend S accordingly.
9: Ĥ ← Repl. Ẑi in Ĥ with sign(Mii ) · 1 10: end if
10: else  Two-point correlation. 11: end if
11: Ĥ ← Repl. Ẑi in Ĥ with sign(Mij ) · Ẑj 12: φ ← {}  Initialize empty formula for comparison.
12: end if 13: while φ = φ and size(φ) > nc do  Inference rules.
13: end for 14: φ ← φ
14: Snc ← Bruteforce(Ĥ)  Solve the remaining problem. 15: φ, S ← PureLiteral(φ, S)
15: S ← Reconstruct(Snc , L)  Extract the solution from L. 16: φ, S ← DominatingUnitClause(φ, S)
return S 17: φ, S ← AlmostCommonClause(φ, S)
18: φ, S ← ComplementaryUnitClause(φ, S)
19: end while
20: if size(φ) ≤ nc then
21: φ, S ← BruteForce(φ, S)
Inference rules are a commonly used building block of 22: end if
many classical satisfiability solvers [31]. They are used to 23: return φ, S
simplify the formula by fixing the assignments of specific
variables and therefore accelerate the remaining computa-
tion. Here, we use only very simple inference rules that (d) The dominating unit clause rule [63]: if the total
have low computational (and conceptual) complexity. The count of clauses, regardless of their length, that con-
selected rules ensure that the optimal solutions of both tain a variable xi (x̄i ) is not greater than the count of
the original and the simplified problems are equally good unit clauses that contain x̄i (xi ), then the variable xi
in quality. Despite their simplicity, the inference rules is set to FALSE (TRUE).
improve the performance of QIRO, both in terms of qual-
ity (see Sec. IV B) and in terms of resource efficiency. This The complete QIRO update step Reduce for MAX-2-
is because fewer calls to the quantum device are required, SAT is shown in Algorithm 5.
due to the simplifications performed by the inference rules. Importantly, the MAX-2-SAT QIRO update rules pre-
The selected rules, listed next, were proposed in Refs. sented here preserve the problem structure in that only
[30,31]. clauses with two or less variables appear in the formula.
Consequently, the Hamiltonian remains in the same form
(a) The pure literal rule [32]: if a literal exclusively as in Eq. (A1).
occurs with positive (negative) polarity then the
value of the corresponding variable is set to TRUE
APPENDIX D: MIS ON A NEUTRAL ATOM
(FALSE).
QUANTUM DEVICE
(b) The almost common clause rule [32]: if a MAX-2-
SAT instance includes the clauses xi ∨ xj and x̄i ∨ Here, we provide details about the adiabatic protocols
xj then both clauses can be replaced with a single- on a neutral atom quantum processor that we used to pre-
literal clause (i.e., unit clause) xj . If neither xi nor x̄i pare the quantum states used to guide the QIRO algorithm
appears in the formula anymore, the variable can be for the MIS problem. The MIS problem on UDGs naturally
arbitrarily set to either TRUE or FALSE. emerges in the context of neutral atom quantum processors
(c) The complementary unit clause rule [63]: if the for- based on Rydberg atom arrays [13,14]. UDGs are a fam-
mula contains single-literal clauses xi and x̄i then ily of graphs, where two nodes are connected only if the
these two clauses are removed. Again, if neither xi (Euclidean) distance between them is smaller than some
nor x̄i appears in the formula thereafter, the assign- threshold radius Rd .
ment of the variable can be arbitrarily set to either Here, we generate such graphs by positioning nodes
TRUE or FALSE. on an underlying square lattice with lattice constant a.

020327-14
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

TABLE I. Values of the hardness parameter HP for the graph instances presented in Fig. 6.

Graph instance 1 2 3 4 5 6 7 8 9 10
HP 3 8 14 39 49 95 197 383 863 1435

√
We set 2a < Rd < 2a, such that we get nearest-neighbor Here, is the Rabi frequency, is the laser detuning,
and diagonal connectivity. We chose a lattice constant and the interaction term Vij = C6 xi − xj −6 2 , with C6 the
a = 5.3 µm, inspired by a previous study [49]. The graphs van der Waals coefficient. The particular value C6 depends
presented here were generated by positioning 137 nodes on the atom species used in the experimental setup. The
on such an underlying lattice with 14 × 14 sites. This cor- number operator n̂i := 12 (12 − Ẑi ) counts the number of
responds to a filling factor of roughly 70%. An example Rydberg excitations on the ith site.
graph is shown in Fig. 5(a). An essential property of the Hamiltonian in Eq. (D2)
Furthermore, we classified the generated instances in is the so-called Rydberg blockade phenomenon, in which
terms of their classical hardness parameter two atoms cannot simultaneously be in the (excited) state
|1 if the distance between them is smaller than the Ryd-
N|MIS|−1
HP := , (D1) berg blockade radius Rb ≡ (C6 / )1/6 [14]. Therefore,
|MIS| N|MIS| if one sets the UDG radius to be equal to the Rydberg
blockade radius (Rd = Rb ), ground states of Ĥcost obey the
where NM denotes the number of independent sets of size
independence constraint on UDGs if one assigns atoms
M . The hardness parameter was shown to be related to
that are in state |1 to the MIS [13]. Additionally, if we
the performance of quantum and classical algorithms alike
set 0 < < Vij , the ground state of Hcost maximizes the
[14]. However, we note that evidence has recently been
number of excitations, while not violating independence,
put forward that HP only influences the performance of
hence corresponding to solutions to the original MIS prob-
Markov chain Monte Carlo–based algorithms (e.g., SA)
lem. We thus start the adiabatic protocol with i < 0, such
[64].
that the ground state is |0⊗|V| and then (slowly) increase
Here, we consider graphs spanning several orders
the detuning until some final positive value, as shown in
of magnitude in HP , which we computed using the
Fig. 8, such that the ground states of the final Hamiltonian
tropical tensor network algorithm provided within the
correspond to the solutions of the MIS problem.
GenericTensorNetworks library [48,65]. Specifi-
We ran such protocols on the QuEra Aquila device,
cally, the graphs considered in Fig. 6 have been sorted in
which is available through the Amazon Braket service. The
order of increasing hardness, with the exact values of HP
detailed parameters can be found in the caption of Fig. 8,
given in Table I.
and are partially motivated by a previous study [49]. At
We now briefly summarize how neutral atom arrays
each stage of the QIRO algorithm we estimate the correla-
can be used to solve the MIS problem on UDGs [13,14].
tions using a hundred measurements of the quantum state
The dynamics of neutral atom quantum devices are gov-
prepared by the adiabatic protocol described in Fig. 8.
erned by the Hamiltonian Ĥ (t) = Ĥdr + Ĥcost , where the
individual terms read
 
Ĥdr = (t) X̂i ,
2 i APPENDIX E: HIGHER-DEPTH QAOA
  (D2)
EXPERIMENTS
Ĥcost = − (t) n̂i + Vij n̂i n̂j .
i i<j To provide additional support to the claim that a higher
quality of correlations yields better performance of QIRO,
Ωmax Δf we performed numerical simulations of QIRO using cor-
relations from QAOA at p > 1. In contrast to p = 1, no
analytical formulae for the required correlators are known
Δ(t)
Ω(t)

0 at p > 1, forcing us to resort to state-vector simulations

Ω(t) of QAOA circuits. Consequently, we had to limit our
Δ(t) experiments to solving the MIS problem on smaller graph
0 Δi
instances with 12 nodes, randomly chosen from the Erdős-
0 τΩ tf − τΩ tf Rényi ensemble as in Sec. IV A. For each fixed expected
degree, we generated a hundred random instances, and
FIG. 8. Time dependency of the Rabi frequency (t) and the found their MIS using a brute-force search of the solution
laser detuning (t). We use tf = 4 µs, max = 15.8 MHz, τ = space. We then ran QIRO informed by correlations from
0.1tf , i = −30 MHz, and f = 60 MHz. the QAOA implementation in the PennyLane library

020327-15
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

100
degree = 3 degree = 6 0.00
Solved instances (%)
0.01 1
95

QAOA parameter quality q

0.02
0.03

RQAOA validity
90
0.04
0.05
85
0.06
y 1 2 3 y 1 2 3
ed = = = ed = = = 0.07
re p p p re p p p
G G
0.08
FIG. 9. The fraction of optimally solved instances by QIRO 0.09 0
using QAOA at depths p ∈ {1, 2, 3}, and by the minimal degree 0.10
greedy algorithm. The height of the bars show the mean fraction
across ten runs of the algorithms, with the error bars indicating

1.01
1.02
1.05
1.1
1.2
1.5
2.0
3.0
5.0
10.0
the best and worst performance across those runs.
Penalty factor λ

[66], using the gradient descent optimizer with 300 iter- FIG. 10. Ratio of valid solutions obtained by RQAOA for
ations and 15 restarts of the optimizer at different ran- different qualities of the QAOA parameters q (at p = 1) and dif-
dom initial values of the variational parameters. For these ferent penalty factors λ. As in the results in Fig. 5(b), at each
experiments, we set nc = 1 because of the comparatively (q, λ) the same 50 UDGs of size 137 were generated and ten runs
smaller problem instance sizes. of RQAOA were performed. We here report the ratio valid (i.e.,
independent) sets found by RQAOA. Black stars indicate values
The results in Fig. 9 show the fraction of optimal inde-
of (q, λ) where all found solutions were feasible.
pendent sets found by QIRO with correlations from QAOA
at depths p ∈ {1, 2, 3}, and by the minimal degree greedy
algorithm. The performance of QIRO clearly improves by an intricate interplay between the deformations of the
with increased QAOA depth, for both considered graph spectrum of the cost Hamiltonian (1) when λ is varied,
densities. Together with the results in Sec. IV A, this pro- with the particularities of QAOA states at p = 1. As λ
vides yet another indication that enhanced correlations is varied, only the energy levels of states that correspond
lead to better performance of QIRO. Moreover, these to configurations with violations of the independence con-
experiments indicate that correlations from deeper circuits straint are changed. Thus, the nature of the excited states
enable QIRO to surpass the performance of the greedy can change dramatically, i.e., states corresponding to con-
benchmark. Although this matches the intuitive expecta- figurations with violations might shift relative to states
tions, further experiments on larger instances and quantum corresponding to independent sets. It is difficult to predict
hardware are required to strengthen these claims. the amplitudes of such excited states in the optimal QAOA
states at p = 1, let alone at perturbed parameters. As such,
APPENDIX F: FURTHER RQAOA VALIDITY a pattern similar to that seen is not inconceivable. How-
EXPERIMENTS ever, more experiments are required to pin down the exact
underlying mechanism.
In Sec. IV A we assessed the performance of QIRO and
RQAOA at different qualities of QAOA parameters. As the
APPENDIX G: PARALLEL TEMPERING AND THE
results presented in Fig. 5(b) exhibited interesting behav-
SIMULATED ANNEALING SETUP
ior, we here report an extended set of experiments, where
RQAOA validities for additional magnitudes of the penalty For a MAX-2-SAT instance with n variables, starting
term λ [cf. Eq. (1)] are shown. from a random initial configuration, we perform simu-
Interestingly, the results partially defy the naïve expecta- lated annealing with 600n attempted variable flips. The
tion that a larger penalty term should favor enforcing con- annealing schedule consists of evenly spaced inverse tem-
straints. This is especially prominent in the region between peratures β = T−1 , with an initial value β = 0 and a final
λ = 1.2 and λ = 5, where the robustness of RQAOA inverse temperature of β = 6.
seems to decrease. However, the trend reverses again for For PT [34–36], we initialize 12 replicas of the prob-
λ = 10. lem in random configurations with different temperatures.
While a thorough analysis of the artefacts observed in We chose the fine-tuned temperatures given in Table II.
Fig. 10 is beyond the scope of this paper, we offer some These were chosen in accordance with the prescription
insights into what might underpin the peculiar behavior at that replica exchanges should be accepted with probabil-
intermediate values of λ. We suspect that this is caused ities between 0.2 and 0.8 [67]. We tested this empirically

020327-16
QUANTUM-INFORMED RECURSIVE OPTIMIZATION ALGORITHMS PRX QUANTUM 5, 020327 (2024)

TABLE II. Temperatures of the replicas in our implementation of PT.

Replica index 1 2 3 4 5 6 7 8 9 10 11 12
Temperature 0.10 0.20 0.29 0.39 0.50 0.62 0.75 0.90 1.09 1.33 1.67 2.20

for a collection of representative MAX-2-SAT instances. [13] H. Pichler, S.-T. Wang, L. Zhou, S. Choi, and M. D. Lukin,
We then used the following PT scheme. In each run of PT Quantum optimization for maximum independent set using
we carried out 15 000 cycles, with one simulation cycle Rydberg atom arrays, ArXiv:1808.10816 (2018).
consisting of one Monte Carlo sweep (n single spin flips) [14] S. Ebadi, et al., Quantum optimization of maximum inde-
pendent set using Rydberg atom arrays, Science 376, 1209
per replica followed by an attempted replica exchange
(2022).
between all replica pairs with neighboring temperatures, [15] J. Wurtz, P. Lopes, N. Gemelke, A. Keesling, and S. Wang,
starting from the lowest temperature. The configuration Industry applications of neutral-atom quantum comput-
with the lowest energy found in any of the replicas dur- ing solving independent set problems, ArXiv:2205.08500
ing the course of the algorithm was finally returned as the (2022).
solution. [16] S. A. Cook, in Proceedings of the Third Annual ACM Sym-
posium on Theory of Computing—STOC ’71 (ACM Press,
Shaker Heights, Ohio, United States, 1971), p. 151.
[1] E. Farhi, J. Goldstone, and S. Gutmann, A quantum approx- [17] P. Cheeseman, B. Kanefsky, and W. M. Taylor, in Pro-
imate optimization algorithm, ArXiv:1411.4028 (2014). ceedings of the 12th International Joint Conference on
[2] D. Gamarnik, The overlap gap property: A topological bar- Artificial Intelligence—Volume 1, IJCAI’91 (Morgan Kauf-
rier to optimizing over random structures, Proc. Natl. Acad. mann Publishers Inc., San Francisco, CA, USA, 1991),
Sci. 118, e2108492118 (2021). p. 331.
[3] S. Bravyi, A. Kliesch, R. Koenig, and E. Tang, Obsta- [18] L. Zeng, J. Zhang, and M. Sarovar, Schedule path optimiza-
cles to variational quantum optimization from symmetry tion for adiabatic quantum computing and optimization,
protection, Phys. Rev. Lett. 125, 260505 (2020). J. Phys. A: Math. Theor. 49, 165305 (2016).
[4] E. Farhi, D. Gamarnik, and S. Gutmann, The quantum [19] M. Mezard, G. Parisi, and R. Zecchina, Analytic and algo-
approximate optimization algorithm needs to see the whole rithmic solution of random satisfiability problems, Science
graph: A typical case, ArXiv:2004.09002 (2020). 297, 812 (2002).
[5] E. Farhi, D. Gamarnik, and S. Gutmann, The quantum [20] J. Sleegers, R. Olij, G. van Horn, and D. van den Berg,
approximate optimization algorithm needs to see the whole Where the really hard problems aren’t, Oper. Res. Perspect.
graph: Worst case examples, ArXiv:2005.08747 (2020). 7, 100160 (2020).
[6] C.-N. Chou, P. J. Love, J. S. Sandhu, and J. Shi, in [21] G. Ochoa, F. Chicano, and M. Tomassini, in Parallel Prob-
49th International Colloquium on Automata, Languages, lem Solving from Nature–PPSN XVI (Springer Interna-
and Programming (ICALP 2022), Leibniz International tional Publishing, Leiden, The Netherlands, 2020), p. 125.
Proceedings in Informatics (LIPIcs), Vol. 229, edited by [22] M. R. Garey and D. S. Johnson, Computers and Intractabil-
M. Bojańczyk, E. Merelli, and D. P. Woodruff (Schloss ity: A Guide to the Theory of NP-Completeness (W. H.
Dagstuhl–Leibniz-Zentrum für Informatik, Dagstuhl, Ger- Freeman & Co., USA, 1990).
many, 2022), p. 41:1. [23] B. Aspvall, M. F. Plass, and R. E. Tarjan, A linear-time
[7] S. Bravyi, A. Kliesch, R. Koenig, and E. Tang, Hybrid algorithm for testing the truth of certain quantified Boolean
quantum-classical algorithms for approximate graph color- formulas, Inf. Process. Lett. 8, 121 (1979).
ing, Quantum 6, 678 (2022). [24] V. Chvatal and B. Reed, in Proceedings 33rd Annual Sym-
[8] J. R. McClean, M. P. Harrigan, M. Mohseni, N. C. Rubin, posium on Foundations of Computer Science (SFCS ’92)
Z. Jiang, S. Boixo, V. N. Smelyanskiy, R. Babbush, and H. (IEEE Computer Society, Pittsburgh, PA, USA, 1992),
Neven, Low-depth mechanisms for quantum optimization, p. 620.
PRX Quantum 2, 030312 (2021). [25] C. H. Papadimitriou and K. Steiglitz, Combinatorial Opti-
[9] H. Karimi, G. Rosenberg, and H. G. Katzgraber, Effective mization, Dover Books on Computer Science (Dover Pub-
optimization using sample persistence: A case study on lications, Mineola, NY, 1998).
quantum annealers and various Monte Carlo optimization [26] E. Farhi, J. Goldstone, S. Gutmann, and M. Sipser, Quan-
methods, Phys. Rev. E 96, 043312 (2017). tum computation by adiabatic evolution, ArXiv:quant-ph/
[10] A. Lucas, Ising formulations of many NP problems, Front. 0001106 (2000).
Phys. 2, 1 (2014). [27] T. Kadowaki and H. Nishimori, Quantum annealing in the
[11] F. Glover, G. Kochenberger, and Y. Du, A tutorial on transverse Ising model, Phys. Rev. E 58, 5355 (1998).
formulating and using QUBO models, ArXiv:1811.11538 [28] A. Ozaeta, W. van Dam, and P. L. McMahon, Expec-
(2019). tation values from the single-layer quantum approximate
[12] S. Butenko and P. M. Pardalos, Ph.D. thesis, USA, 2003, optimization algorithm on Ising problems, Quantum Sci.
aAI3120100. Technol. 7, 045036 (2022).

020327-17
JERNEJ RUDI FINŽGAR et al. PRX QUANTUM 5, 020327 (2024)

[29] P. Hauke, H. G. Katzgraber, W. Lechner, H. Nishimori, and [47] L. Bittel and M. Kliesch, Training variational quantum
W. D. Oliver, Perspectives of quantum annealing: Meth- algorithms is NP-hard, Phys. Rev. Lett. 127, 120502
ods and implementations, Rep. Progr. Phys. 83, 054401 (2021).
(2020). [48] J.-G. Liu, X. Gao, M. Cain, M. D. Lukin, and S.-T. Wang,
[30] A. Abramé and D. Habet, ahmaxsat: Description and eval- Computing solution space properties of combinatorial opti-
uation of a branch and bound Max-SAT solver, J. Satisf. mization problems via generic tensor networks, SIAM J.
Boolean Model. Comput. 9, 89 (2015). Sci. Comput. 45, A1239 (2023).
[31] C. M. Li and F. Many, MaxSAT, hard and soft constraints, [49] J. R. Finžgar, M. J. A. Schuetz, J. K. Brubaker, H. Nishi-
Front. Artif. Intell. Appl. 185, 613 (2021). mori, and H. G. Katzgraber, Designing quantum annealing
[32] B. Borchers and J. Furman, A two-phase exact algorithm schedules using Bayesian optimization, ArXiv:2305.13365
for MAX-SAT and weighted MAX-SAT problems, (2023).
J. Comb. Optim. 2, 299 (1998). [50] S. Kirkpatrick, C. D. Gelatt, and M. P. Vecchi, Optimization
[33] Y. J. Patel, S. Jerbi, T. Bäck, and V. Dunjko, Reinforcement by simulated annealing, Science 220, 671 (1983).
learning assisted recursive QAOA, ArXiv:2207.06294 [51] A. Ignatiev, A. Morgado, and J. Marques-Silva, RC2: An
(2022). efficient MaxSAT solver, J. Satisf. Boolean Model. Com-
[34] K. Hukushima and K. Nemoto, Exchange Monte Carlo put. 11, 53 (2019).
method and application to spin glass simulations, J. Phys. [52] S. Hadfield, Z. Wang, B. O’Gorman, E. Rieffel, D. Ven-
Soc. Jpn. 65, 1604 (1996). turelli, and R. Biswas, From the quantum approximate
[35] D. J. Earl and M. W. Deem, Parallel tempering: Theory, optimization algorithm to a quantum alternating operator
applications, and new perspectives, Phys. Chem. Chem. ansatz, Algorithms 12, 34 (2019).
Phys. 7, 3910 (2005). [53] M. B. Hastings, Classical and quantum bounded depth
[36] F. Romá, S. Risau-Gusman, A. Ramirez-Pastor, F. approximation algorithms, ArXiv:1905.07047 (2019).
Nieto, and E. Vogel, The ground state energy of the [54] C. Huang, et al., Efficient parallelization of tensor net-
Edwards–Anderson spin glass model with a parallel tem- work contraction for simulating quantum computation, Nat.
pering Monte Carlo algorithm, Physica A 388, 2821 Comput. Sci. 1, 578 (2021).
(2009). [55] J. Gray and S. Kourtis, Hyper-optimized tensor network
[37] Z. Zhu, C. Fang, and H. G. Katzgraber, Borealis–A gen- contraction, Quantum 5, 410 (2021).
eralized global update algorithm for Boolean optimization [56] M. W. Johnson, et al., Quantum annealing with manufac-
problems, Optim. Lett. 14, 2495 (2020). tured spins, Nature 473, 194 (2011).
[38] S. Boettcher and A. G. Percus, in Computational Modeling [57] A. D. King, et al., Coherent quantum annealing in a pro-
and Problem Solving in the Networked World (Springer US, grammable 2000 qubit Ising chain, Nat. Phys. 18, 1324
Chandler, Arizona, 2003), p. 61. (2022).
[39] M. Mahajan and F. Slivovsky, in 26th International Con- [58] L. Caha, A. Kliesch, and R. Koenig, Twisted hybrid algo-
ference on Theory and Applications of Satisfiability Testing rithms for combinatorial optimization, Quantum Sci. Tech-
(SAT 2023) (Schloss Dagstuhl-Leibniz-Zentrum für Infor- nol. 7, 045013 (2022).
matik, Alghero, Italy, 2023). [59] F. Wagner, J. Nüßlein, and F. Liers, Enhancing quantum
[40] P. Erdös and A. Rényi, On random graphs I, Publ. Math. algorithms for quadratic unconstrained binary optimization
Debrecen 6, 290 (1959). via integer programming, ArXiv:2302.05493 (2023).
[41] M. Dupont, B. Evert, M. J. Hodson, B. Sundar, S. Jeffrey, [60] [Link]
Y. Yamaguchi, D. Feng, F. B. Maciejewski, S. Hadfield, M. [61] E. Farhi and A. W. Harrow, Quantum supremacy
S. Alam, Z. Wang, S. Grabbe, P. A. Lott, E. G. Rieffel, D. through the quantum approximate optimization algorithm,
Venturelli, and M. J. Reagor, Quantum-enhanced greedy ArXiv:1602.07674 (2019).
combinatorial optimization solver, Sci. Adv. 9, eadi0487 [62] E. Bae and S. Lee, Recursive QAOA outperforms the origi-
(2023). nal QAOA for the MAX-CUT problem on complete graphs,
[42] M. M. Halldórsson and J. Radhakrishnan, Greed is good: ArXiv:2211.15832 (2023).
Approximating independent sets in sparse and bounded- [63] R. Niedermeier and P. Rossmanith, New upper bounds for
degree graphs, Algorithmica 18, 145 (1997). maximum satisfiability, J. Algorithms 36, 63 (2000).
[43] M. Böther, O. Kißig, M. Taraz, S. Cohen, K. Sei- [64] R. S. Andrist, M. J. A. Schuetz, P. Minssen, R. Yalovetzky,
del, and T. Friedrich, in The Tenth International Con- S. Chakrabarti, D. Herman, N. Kumar, G. Salton, R. Shay-
ference on Learning Representations, ICLR 2022, Vir- dulin, Y. Sun, M. Pistoia, and H. G. Katzgraber, Hardness of
tual Event, April 25–29, 2022 ([Link], 2022), the maximum independent set problem on unit-disk graphs
[Link] and prospects for quantum speedups, ArXiv:2307.09442
[44] S. Lamm, P. Sanders, C. Schulz, D. Strash, and R. F. Wer- (2023).
neck, Finding near-optimal independent sets at scale, J. [65] J.-G. Liu, L. Wang, and P. Zhang, Tropical tensor network
Heurist. 23, 207 (2017). for ground states of spin glasses, Phys. Rev. Lett. 126,
[45] D. Hespe, C. Schulz, and D. Strash, Scalable kernelization 090506 (2021).
for maximum independent sets, ACM J. Exp. Algorithmics [66] V. Bergholm, et al., PennyLane: Automatic differentiation
24, 1 (2019). of hybrid quantum-classical computations, ArXiv:1811.
[46] J. Lee, A. B. Magann, H. A. Rabitz, and C. Arenz, Progress 04968 (2022).
toward favorable landscapes in quantum combinatorial [67] H. G. Katzgraber, Introduction to Monte Carlo methods,
optimization, Phys. Rev. A 104, 032401 (2021). ArXiv:0905.1629 (2011).

020327-18