Quantum-Informed Recursive Optimization Algorithms

We propose and implement a family of quantum-informed recursive optimization (QIRO) algorithms for combinatorial optimization problems. Our approach leverages quantum resources to obtain information that is used in problem-specific classical reduction steps that recursively simplify the problem. These reduction steps address the limitations of the quantum component and ensure solution feasibility in constrained optimization problems. Additionally, we use backtracking techniques to further improve the performance of the algorithm without increasing the requirements on the quantum hardware. We demonstrate the capabilities of our approach by informing QIRO with correlations from classical simulations of shallow (depth $p=1$) circuits of the quantum approximate optimization algorithm (QAOA), 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 available online on Amazon Braket to find large independent sets of graphs. In summary, our scheme achieves results comparable to classical heuristics, such as simulated annealing and greedy algorithms, even with relatively weak quantum resources. Furthermore, enhancing the quality of these quantum resources improves the performance of the algorithms, highlighting the potential of QIRO. Notably, the modular nature of QIRO offers various avenues for modifications, positioning our work as a blueprint for designing a broader class of hybrid quantum-classical algorithms for combinatorial optimization.


I. INTRODUCTION
Quantum optimization has been identified as a promising area of research towards practical quantum advantage.On noisy intermediate-scale quantum (NISQ) devices, much effort has been dedicated to studying hybrid quantum-classical algorithms such as the quantum approximate optimization algorithm (QAOA) [1].Importantly, QAOA is a local algorithm.This means that at any constant circuit depth, only qubits that are separated by less than a certain distance in the interaction graph of the optimization problem are able to communicate.Together with the solution space properties of certain optimization problems such as the Z 2 symmetry, or the overlap gap property [2], locality was shown to severely limit the performance of QAOA [3][4][5][6].
Because introducing non-local updates in quantum algorithms comes with additional hardware requirements, applying non-local updates classically has been proposed as an alternative in recursive QAOA (RQAOA) [3,7].Here, values of variables are iteratively frozen by rounding the correlations between variables as measured in the quantum state prepared by QAOA.Hence, as variables are removed from the optimization problem the distances between nodes in the interaction graph are reduced iteratively.As such, the non-local effects introduced via the new connections between previously unconnected nodes counterbalance the locality inherent to QAOA.
FIG. 1. Schematic visualizing the core principles of the quantum-informed recursive optimization (QIRO) algorithm presented in this work.Quantum resources (e.g., QAOA) are used to obtain information, e.g., in the form of one-point and two-point correlations, respectively.This information is then used to simplify the problem through problem-specific classical update rules.Here we exemplify the update rules with simple examples for the maximum independent set and satisfiability problems -details on the update rules can be found in Section 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 been fully simplified.
torial optimization [8], and tailor the classical subroutines to the particular optimization problem of interest, thereby enhancing the algorithm's performance.
Moreover, the inclusion of problem-specific update rules comes with additional benefits.It allows us to broaden the scope of our algorithm beyond gate-based architectures to analog devices, as shown in our experiments with the QuEra Aquila neutral atom analog quantum processor.This approach is reminiscent of previous approaches leveraging spin-freezing schemes on quantum annealers and in various Monte Carlo methods [9].Furthermore, for problems with hard constraints, the update rules can offer options to enforce feasibility by design, even in the presence of noise.
Finally, we propose backtracking as a strategy to further improve the performance of QIRO, by attempting to identify and correct non-ideal decisions made at earlier stages of the algorithm.Backtracking provides a way to enhance algorithmic performance without necessitating an increase in circuit depth, as is commonplace in quantum optimization algorithms.Taken together, the robustness of QIRO, its applicability to analog devices, and the possibility to improve its efficacy without demanding better quantum resources, make QIRO a promising NISQ algorithm.
In this work we provide concrete implementations for two paradigmatic NP-hard combinatorial optimization problems, namely, maximum independent set and maximum satisfiability.However, our work should be seen mainly as a proof-of-concept demonstration of a general template for designing hybrid quantum-classical algorithms for combinatorial optimization.We study the performance of QIRO by means of large-scale simulations of low-depth QAOA circuits with up to two hundred variables.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.In Section II we introduce the two problem classes studied in this work, in Section III we provide details about QIRO and in Section IV we study the performance of QIRO on the maximum independent set and maximum satisfiability problems.We discuss the implications of our findings in Section V and we conclude by suggesting potential extensions to our framework in Section VI.

II. PROBLEM CLASSES
Because problem-specific update rules are at the core of QIRO, we first introduce the specific optimization problems considered in this work.We formalize these problems, and provide mappings to quantum Hamiltonians Ĥc , which can in turn be implemented on different quantum devices [10,11].

A. Maximum Independent Set
The maximum independent set (MIS) problem is a paradigmatic NP-hard combinatorial optimization problem.It is of commercial relevance in areas such as network design or traffic optimization [12].For a graph G = (V, E) with vertex set V and an edge set E, the MIS problem amounts to finding the largest subset S ⊆ V such that no two vertices in S are adjacent [see Fig. 5a below for an example].For a binary vector x ∈ {0, 1} |V | finding the MIS is equivalent to finding the ground state of the following (classical) Hamiltonian: Here, λ > 1 is a penalty term to enforce the independence constraint [10].However, one should be mindful of the pitfalls of enforcing the independence constraint as a soft penalty term, because low-energy states (other than the ground states) of H(x) might not necessarily obey the independence constraint.Performing the mapping z i = 2x i − 1 to Ising variables z i = ±1, and promoting the spin variables to quantum operators z i → Ẑi translates Eq. ( 1) to a quantum Hamiltonian suitable for quantum devices [10].Moreover, the MIS problem on a family of graphs called unit disk graphs (UDGs) has been shown to have a natural mapping to neutral atom quantum devices [13,14].There, the Rydberg blockade mechanism-which prevents simultaneous excitations of nearby atomsenforces the independence constraint.Finding the MIS of a general UDG, despite being a specialized family of graphs, remains an NP-Hard problem [13] and has commercial applications (e.g., in network design [15]).For additional details we refer the reader to Appendix D.

B. Maximum Satisfiability
The satisfiability problem (SAT) is arguably one of the 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 x i ∈ {0, 1}, SAT is a decision problem which asks for the existence of an assignment of the variables x i that satisfies ϕ.Each clause can then be written as where literals ℓ j k are either a variable x j k or its negation xj k , for j k ∈ {1, . . ., n}, and K j is the length of the j-th clause.For a fixed K j = K the corresponding problem of whether or not the Boolean formula with m clauses ϕ = m j=1 Cj has a satisfying assignment, is referred to as K-SAT.
Here, we focus on the optimization version of satisfiability, referred to as maximum satisfiability (MAX-SAT).In MAX-SAT, the goal is to find an assignment of the variables x i that minimizes the number of violated clauses.The MAX-K-SAT problem deals with propositional logic formulae with at most K literals per clause.An important parameter of a (MAX-)K-SAT problem is the clause-to-variable ratio α := m/n, which strongly influences the algorithmic hardness of random problem instances [17].Research indicates that there is a phase transition in computational complexity [18].While for K-SAT an easy-hard-easy transition is observed [19,20], MAX-K-SAT exhibits a transition from easy underconstrained to hard overconstrained problems at a certain α c [21].In what follows we limit ourselves to K = 2 (clauses contain either one or two literals), as MAX-2-SAT is already NP-hard [22], despite its decision version (2-SAT) being solvable in linear time [23].Here, the phase transition occurs at α c = 1 [24].
We use the fact that for Boolean (binary) variables xi ≡ (1 − x i ), solving MAX-2-SAT for a given formula ϕ is equivalent to finding a ground state of the following (classical) Hamiltonian: ( It is important to keep in mind that ℓ j k are nothing but the Boolean variables x j k (or their negations xj k ).Hence, by using the transformation into Ising variables outlined in the case of the MIS problem in Section II A, we can likewise encode the MAX-2-SAT problem onto quantum devices.

III. ALGORITHM
We now provide details about QIRO.Essentially, in QIRO quantum resources are utilized to inform update rules that recursively simplify the problem at hand (see Fig. 1).The idea is a generalization of previous work, in particular RQAOA [3,7] (also described in Appendix B).However, in contrast to RQAOA, the update steps used in QIRO are problem-specific.This comes with several advantages and facilitates cross-pollination with ideas from classical optimization [8].
At each iteration of the QIRO algorithm, for a given problem Hamiltonian we first prepare a low-energy quantum state, generically in the form of a superposition of low-energy candidate solutions to the optimization problem.Next, we use this quantum state to extract information (e.g., correlations between variables), which is then passed to a classical update step with the goal of simplifying the problem.We design this update step to complement the quantum part of the algorithm by performing operations better suited to classical hardware.This may include non-local steps that address the identified limitations of (local) quantum algorithms, or enforcing hard constraints.As such, the classical update step (and the specific design thereof) is crucial to the performance of QIRO, and is typically problem-specific.As such, it might need to be modified by the end user, depending on their needs.However, templates for the update steps may be borrowed from the literature on classical combinatorial optimization [25].The procedure of quantum state preparation and the subsequent classical update step is repeated until the size of the problem is reduced sufficiently, such that, for example, we are able to solve it exactly by a brute-force search or other means if the reduction does not yield problem sizes amenable to exact solvers.Finally, we propose the use of backtracking to identify and rectify problem reductions that led to an unfavorable outcome, thus providing an additional strategy to obtain enhanced solutions without requiring deeper quantum circuits.
In what follows, we provide details on the individual components of QIRO.We begin by describing the quantum state preparation methods utilized in this work, followed by the proposed update rules for the MIS and MAX-2-SAT problems.Next, we introduce the general strategy of backtracking and describe the particular implementation used in this work.Finally, we combine the building blocks to present the complete QIRO algorithm.

A. Quantum State Preparation
It is desirable that the quantum state which is prepared at the beginning of each iteration has a low expectation value with respect to the cost Hamiltonian Ĥc .Therefore, we hope that the information guiding the classical reduction steps is extracted from (a superposition of) configurations corresponding to good candidate solutions to the optimization problem.A plethora of protocols for preparing low-energy states may be found in the literature; in this manuscript we consider QAOA [1] and adiabatic state preparation [26,27], as outlined next.

QAOA
Let us first introduce QAOA.For both the MAX-2-SAT and MIS problems considered in this manuscript, the cost Hamiltonian can be written in the general (Ising) form Ĥc = i,j J ij Ẑi Ẑj + i h i Ẑi , where we have omitted any constant terms.Using the cost Hamiltonian Ĥc from Section II, and defining the mixer Hamiltonian as Ĥmix = − n i=1 Xi , the QAOA circuit is given by an alternating application of the time evolution generated by these two operators.Specifically, taking the initial state as |ψ 0 ⟩ := |+⟩ ⊗n , the state produced with QAOA at depth p is The 2p parameters (β, γ) are then determined by a classical optimization routine such that ⟨ψ(β, γ)| Ĥc |ψ(β, γ)⟩ is minimized.
In this work, we predominantly focus on the lowest QAOA depth (p = 1), with the goal to establish a lower bound for the performance of QIRO.It is natural to expect that the improved quality of the quantum information extracted from higher depth p > 1 QAOA circuits (in the absence of noise) would enhance the performance of the algorithm.Our p > 1 experiments on small problem instances shown in Appendix E corroborate this intuition.
For the state produced by QAOA at p = 1 the expectation values of the single-and two-body terms comprising Ĥc can be calculated classically in polynomial time, enabling us to perform numerical simulations of systems with several hundred variables [3,7,28].We note that there is no method known (to us) how this could be done efficiently for p > 1.For details about our implementation we refer the reader to Appendix A.

Adiabatic state preparation
We next turn or attention to adiabatic state preparation.We present a general description of the protocol and defer the details pertaining to our implementation of the MIS problem on neutral atom quantum processors to Appendix D.
At the beginning of an adiabatic protocol, we start in an easy-to-prepare ground state of a simple Hamiltonian, e.g., Ĥmix .We then slowly transform the Hamiltonian from Ĥmix to Ĥc .If the rate of change is slow enough compared to the smallest instantaneous spectral gap (between the instantaneous ground and first excited states), the adiabatic theorem guarantees that at the end of the protocol, the ground state of Ĥc is obtained [26,27].Inspired by adiabatic protocols, the quantum annealing computational paradigm has been subject of widespread interest in recent decades [29].In quantum annealing one forgoes the strict requirement of adiabaticity and, in general, uses analog protocols to prepare near-ground states.These can, in turn, be used to inform the classical update steps.

B. Update Rules
We now discuss example update rules for the considered optimization problems.We begin by preparing a low-energy quantum state |ψ⟩ from which information is distilled.In QIRO this information comes in the form of one-point ⟨ψ| Ẑi |ψ⟩ and two-point correlations ⟨ψ| Ẑi Ẑj |ψ⟩, which we store in the matrix M (see Algorithm 1).Potential generalizations, such as using expectation values of other observables, are left for future research.
If the lowest-depth (p = 1) of QAOA is used to generate these correlations, the only nonzero values of M in an ideal numerical simulation as described in Appendix A (i.e., free of shot and environmental noise) correspond to pairs of variables that are connected in the original interaction graph, i.e., to pairs of variables with nonzero coupling coefficients J ij .Therefore, we design our update rules to only consider correlations between such connected variables.Inspired by RQAOA [3,7], we next use these correlations to successively simplify the problem with each update step.However, in contrast to RQAOA, the QIRO update steps are problem-specific in that they exploit the structure of the given optimization problem.
As such, at any stage of the algorithm, the problem that remains to be solved is a valid instance of the same class of optimization problems.This will be crucial when deploying QIRO to solve the MIS problem on neutral atom quantum hardware.
In the two specific cases detailed in the remainder of this section, only the largest correlation in terms of its absolute value is used to inform the update step, with ties broken at random.Intuitively, the largest correlation corresponds to a variable (pair) which appears in a certain configuration with high probability in the lowenergy state |ψ⟩.Thus, it is reasonable to use this information to inform the update step.For example, if a large positive correlation ⟨ Ẑi Ẑj ⟩ ≈ 1 is measured, it is likely that most of the low-energy configurations have x i and x j take the same value.Generalizations of this deterministic scheme to randomized (e.g., evolutionary) schemes will be studied in the future.
As we lay out next, these update steps are carried out in the original formulation of the respective optimization problems, facilitating the adoption of techniques from the extensive literature on classical combinatorial optimization.We call this update step Reduce and describe it in Algorithm 1.
▷ Problem-specific simplification is informed by the correlations stored in M .4: return P ′ , S ′ The problem-specific updates employed for the optimization problems considered in this work constitute the Simplification function in Algorithm 1, and are described next.However, we emphasize that the design choices we make are by no means the only option, with many modifications and extensions possible.

MIS Simplification
At the beginning of the Simplification routine for MIS, we find the entry of M with the largest absolute value.Depending on the sign of the correlation with the largest absolute value and on whether it lies on the diagonal (corresponding to a one-point correlation) or not (corresponding to a two-point correlation), we perform a different simplification.The simplifications are designed such that the correlation-informed reductions are consistent with the independence constraint.The four cases are visualized in Fig. 2 and described next: (a) If M ii ≥ 0 was selected, we set the i-th vertex to be in the independent set (IS).We then remove all vertices connected to the i-th node from the graph, as including them would violate the independence constraint.
(b) If M ii < 0 was selected, we remove the i-th vertex from the graph.
was selected, we remove both nodes from the graph -this is because the only positively correlated assignment of connected vertices consistent with the independence constraint is not to include either node in the IS.
(d) If M ij < 0 [for (i, j) ∈ E] was selected, we remove nodes which are connected to both node i and node j simultaneously, i.e., we remove every node k for which (k, i) ∈ E and (k, j) ∈ E. Intuitively, we know that if the variables are negatively correlated, one of them will be assigned to be in the IS.Thus, we only remove variables which are connected to both of them simultaneously.
If no nodes were removed from the graph using the described simplifications [only possible in the case (d)], we repeat the procedure using the next largest correlation in terms of its absolute value.Finally, we identify all connected components (i.e., connected subgraphs that are not part of any other connected subgraph) of the graph.If the connected component contains less than a certain number of vertices n c (for concreteness, we use n c = 15 here), we compute its maximum independent set by a brute-force search, and remove it from the graph.This step is mostly performed in an effort to save resources when the algorithm is deployed on real quantum hardware.
Crucially, the rules delineated above explicitly ensure that the independence constraint is obeyed.Hence, the solutions obtained by repeated applications of Simplification are definitely feasible.This holds true irrespective of the quality of the (quantum) information used to perform the simplification, thus bestowing the scheme with the resilience required in the NISQ era.However, although the update rules improve the performance of the algorithm by ensuring solution feasibility, the quality of the obtained solutions ultimately hinges on the quality of the supplied correlations.Interestingly, in contrast to the RQAOA update rules the QIRO MIS update rules proposed here are local according to the definition of locality in Section I, as no new connections between variables are created.However, the RQAOA update rules lead to intermediate Hamiltonians that generically do not take the form of an MIS Hamiltonian as in Eq. (1).Conversely, at each stage of the QIRO algorithm the reduced problem is a valid MIS problem (see also Fig. 2).Thus, if one can encode the original MIS problem onto an analog quantum device (e.g., a neutral atom quantum processor), the same can be done with intermediate, simplified problems.As such, one might be able to harness the nonlocal effects generated by the quantum many-body time evolution on a neutral atom quantum processor, circumventing the issues of local quantum algorithms.Previously, such adiabatic protocols have been shown to generate solutions surpassing those obtained by running QAOA, even at depths higher than p = 1 used in this manuscript [14].

MAX-2-SAT Simplification
In the case of MAX-2-SAT the Simplification step 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, we apply a different reduction: (a) If M ii was selected, we set the variable x i → [sign(M ii ) + 1]/2, i.e., a positive (negative) onepoint correlation corresponds to assigning x i to TRUE (FALSE).
(b) If M ij (for i < j) was selected, we replace Intuitively, this update step captures the relationship between variables x i and x j as inferred from the correlation.If there is a positive (negative) correlation, we assign x i the same (opposite) value as x j .
We note that only the two-point correlation reduction step described in (b) constitutes a nonlocal update because it potentially introduces a new link between previously unconnected variables.Because we only consider formulae in conjunctive normal form (i.e., conjunctions of disjunctions), setting a certain variable to the value required by the clause is sufficient to satisfy the clause (see the example update step in Fig. 1).Hence, we can remove the clause from consideration, thus simplifying the problem.Specifically, all clauses in which a given literal evaluates to TRUE can be removed.Conversely, if the literal evaluates to FALSE, only the literal itself can be removed from the clause.If a certain clause becomes empty (i.e., if both literals were removed from it), it corresponds to an unsatisfied clause.
Next, we apply inference rules, which is a standard technique in the field of satisfiability solvers [30,31].Inference rules deduce information from the structure of the current Boolean formula, simplifying it by assigning values to certain variables and thus speeding up further computation.The rules guarantee that the optimal solutions of the current and simplified problem are of equal quality.The chosen inference rules (e.g., the pure literal rule [32]) are conceptually simple rules that incur a small computational overhead.Further details on the specific inference rules we employ are provided in Appendix C. Finally, if the number of remaining variables falls below n c (we use n c = 10 for MAX-2-SAT), the solution of this simplified problem is found by a brute-force search of the remaining solution space.

C. Quantum Informed Recursive Optimization Algorithm
We now possess all of the ingredients for the quantuminformed recursive optimization (QIRO) algorithm.At each step of QIRO, we begin by preparing a quantum state (as described in Section III A).We proceed by using the correlation information by following the problemspecific prescriptions from Section III B. The simple endto-end QIRO procedure is summarized in Algorithm 2.

D. Backtracking
We now turn our attention to backtracking.QIROas specified in Algorithm 2-is a polynomial-time algorithm.Therefore, under common complexity theoretical assumptions, there exist instances of NP-hard problems that QIRO will not be able to solve to optimality.This means that some of the reduction steps performed by QIRO to arrive at a candidate solution might be wrong.A similar observation was made for RQAOA in Ref. [33], prompting the authors to consider reinforcement learning to steer RQAOA update steps.Here, we propose backtracking to derive improved solutions starting from an initial candidate solution obtained by QIRO.Through backtracking, we obtain a method to refine solutions by using the same quantum resources more frequently.This approach stands in stark contrast to most other quantum optimization algorithms, which usually require enhanced quantum resources (e.g., increasing the depth p in QAOA) to deliver improved solutions.
In this work, we only employ this strategy for the MAX-2-SAT problem.However, extending this approach to other optimization problems such as the MIS problem is straightforward.
As visualized in Fig. 3 the simplification process of QIRO can be represented by means of a search tree.Decisions (i.e., update steps) made by the algorithm correspond to arrows in the search tree, which lead to simplified problems represented by nodes.Following the update rules in Section III B (i.e., successively applying Reduce) traces out a path through the search tree (see red arrows and nodes in Fig. 3), which eventually leads to a candidate solution (see Solution (0•••00) in Fig. 3).Here, the superscript (0 • • • 00) indicates that this is the solution obtained in the initial passage through the search tree.
We would like to consider alternative solutions, obtainable by deviating from the initial path through the search tree which led to the initial candidate solution Solution (0•••00) .To this end, we define the Reduce function, which reverses the decision previously made at the revisited node (i.e., partial problem and solution).We emphasize that this initial reversal of the previously made decision requires no additional measurements from the quantum device -we simply perform the opposite reduction to the one we performed initially.For example, if at first the decision was made to assign the value of a variable x i = TRUE, we reverse that decision and set x i = FALSE instead.Thus, we diverge from the original path through the search tree as visualized in Fig. 3.We note that Reduce requires the information about the solution it is deviating from in order to be able to perform the reversed decision.
In principle, exploring the entire search tree essentially corresponds to a brute-force search of the complete solution space.Because we are interested in designing heuristic algorithms that run in polynomial time, we limit the permitted exploration.We find that a simple rule for FIG. 3. Schematic visualization of the backtracking procedure used to obtain improved solutions.Starting from the original problem P (0) and an empty solution S (0) , we perform a series of reductions (red arrows & nodes) that recursively simplify the problem, until the first candidate solution (Solution (0•••00) ) is reached.This corresponds to applying QIRO as outlined in Algorithm 2. Backtracking, as visualized by the dashed teal arrows, then amounts to revisiting reduction steps made along the path.At the chosen ancestor node we have backtracked to, we can then make an alternative reduction, here denoted by Reduce (see black arrows in Fig. 3).Concretely, Reduce amounts to reversing the initially made decision, not needing extra calls to the quantum device.In this work, we limit ourselves to making a binary decision at each node.As such, we can label the reduced problem and partial solution at each node with a binary string encoding the path connecting the node to the original problem node.Each 0 (1) in the superscript labelling P and S corresponds to one application of Reduce (Reduce) starting from the original problem.
which the number of QAOA circuit executions scales quadratically in the problem size is sufficient to produce satisfactory results for the considered problem sizes of MAX-2-SAT.However, more (or less) sophisticated approaches can be designed and tailored to the particular needs and resources of the end user.Next, we describe the specific backtracking strategy used in our simulations and refer the reader to Algorithm 3 for the pseudocode.Append (P, S) to L.
First, we generate the initial candidate solution Solution (0•••00) by successively applying Reduce as visualized in Fig. 3. Put differently, we first apply QIRO without backtracking.Next, we backtrack to each of the ancestor nodes of the initial candidate solution (as visualized by dashed teal lines in Fig. 3).Each of the revisited parent nodes corresponds to a different simplified problem and a corresponding partial solution.The decision made in the initial passage through the search tree is then reversed (by applying Reduce), whereby a new node in the search tree (e.g., P (01) , S (01) in Fig. 3) is created.Here, each 0 (1) in the superscript labelling P and S corresponds to one application of Reduce (Reduce) starting from the original problem.Subsequently, if the problem has not been fully reduced yet, from this new node a new candidate solution is generated by repeatedly applying Reduce.In contrast to the initial application of Reduce the subsequent applications of Reduce require additional calls to the quantum device -see the definition of Reduce in Algorithm 1.In other words, each revisited node leads to a new candidate solution that is identical to the initial candidate solution in the values of the variables fixed before the revisited node.From the revisited node onward the new candidate solution differs from the initial solution.Finally, the output of the algorithm is the best candidate solution found in this way.We note that the order in which the parent nodes are revisited is not important.

IV. RESULTS
Next, we present the results obtained by applying QIRO to random instances of the MIS and MAX-2-SAT problems.In numerical simulations we compare the performance of quantum algorithms (QIRO, RQAOA) and two prominent classical heuristics (simulated annealing, greedy algorithms).Of course more complex tailored classical approaches can be used, such as the paral-lel tempering [34][35][36], borealis [37], and extremal optimization [38] heuristics, as well as a plethora of exact solvers benchmarked at the annual (MAX-)SAT competition [39].
Efficient classical simulations of p = 1 QAOA allow us to study the system size scaling on problem instances with various connectivities of the underlying interaction graph.Results using the lowest depth of QAOA should mainly be interpreted as a lower bound on the performance of QIRO.Furthermore, we investigate the robustness of QIRO and RQAOA against perturbations of the optimal parameters of the QAOA circuits generating the correlations.This allows us to emulate failures of the classical optimizer to find the optimal parameters.Notably, we also deploy QIRO on the QuEra Aquila neutral atom device (as available on Amazon Braket), to (approximately) solve the MIS problem on unit disk graphs with more than a hundred nodes.

A. MIS
We begin by presenting the results of QIRO and RQAOA as applied to the MIS problem.First, we consider the algorithmic performance as a function of the system size on the well-studied family of Erdős-Rényi graphs, with fixed average degrees [40].As a classical comparison, we use two variants of a greedy algorithm.This choice is motivated by recent findings that information from quantum devices can be used to enhance the performance of greedy solvers [41].First, we consider a simple random greedy algorithm, which iteratively assigns random nodes to the independent set, and removes nodes adjacent to the selected node such that no violations are possible.As the second comparison we use the minimal degree greedy algorithm, which iteratively assigns a node with the minimal degree to the independent set and likewise removes adjacent nodes to ensure the validity of the final solution [42].We note that the minimal degree greedy algorithm was designed chiefly for sparse graphs (i.e., of low average degree).This is indicated by the performance guarantee of the greedy algorithm, which is inversely proportional to the graph's average degree [42].However, we note that even for the sparse graphs considered throughout this section, stronger algorithms exist in the literature -see Ref. [43] for benchmarks.
Importantly, because all of the algorithms considered here operate by iteratively simplifying the problem, their performance is indicative of the fitness of the information used to perform these simplifications.
For each system size and average degree we generated a set of fifty random Erdős-Rényi graphs and solved these with the ReduMis heuristic solver from the KaMIS project [44,45].We use the solution obtained by ReduMis as (a proxy for) the optimal solution, and report the percentage of graphs for which each of the algorithms considered performs on par with the ReduMis reference.FIG. 4. System size scaling for the success rate of the random greedy (Rnd.greedy), minimal degree greedy (Min.greedy), RQAOA, and QIRO for Erdős-Rényi graphs with different expectation values of the average degree.At each system size and average degree fifty graphs were generated, and ten runs of each algorithm were performed.We report the median percentage of graphs for which we find independent sets of the size found by the classical ReduMis benchmark algorithm.
The error bars denote the minimum and maximum performance observed across the ten runs.Note that the results for QIRO and RQAOA are obtained using simulations of QAOA at depth p = 1.
The results for graph sizes with n = 40 -200 nodes and average degrees of 3, 5, and 12 are shown in Fig. 4. We used λ = 1.1 in the MIS Hamiltonian (see Eq. ( 1)) for RQAOA and QIRO.The performance of the two quantum algorithms (QIRO and RQAOA) is found to lie between the random and minimal degree greedy algorithms.This highlights the limitations of the quantum correlations obtained from QAOA at p = 1, which only "sees" the immediate neighbors of every node.As such, it is unsurprising that the minimal degree greedy algorithm -which makes the optimal greedy move based on information from nearest neighbors-outperforms QIRO and RQAOA informed by p = 1 QAOA.However, our results in Appendix E indicate that (at least on small graphs), QIRO with QAOA at p > 1 outperforms the greedy benchmark.
Furthermore, RQAOA (slightly, but consistently) outperforms QIRO.This is likely due to the fact that the  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 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).
QIRO simplification rules as defined in Sec.III B 1 are manifestly local.The fact that the gap between QIRO and RQAOA seems to diminish for denser graphs (i.e., for increasing average degree) corroborates this intuitive explanation, as we expect locality to be a weaker limitation if the underlying interaction graph has a higher connectivity.
A major advantage of QIRO (as compared to, e.g., RQAOA) is that the problem-specific classical subroutines can be used to enforce feasibility of the solutions generated (i.e., all obtained candidate solutions are guaranteed to be independent sets).To showcase this feature, we performed numerical simulations of QIRO and RQAOA with p = 1 QAOA with suboptimal parameters.The goal of these experiments is to consider the robustness of the algorithms against lower-quality quantum information.The parameters are generated by sampling a uniform 30 × 30 grid in the (β, γ)-parameter plane, and choosing the parameters corresponding to q-th quantile in the energy values, among the 30 × 30 energies obtained.Given their relevance for neutral atom quantum devices [13,14], we use unit disk graphs (UDGs) (see Appendix D) for these experiments.
The results for an ensemble of fifty UDGs with 137 nodes and two different values of the penalty term λ are shown in Fig. 5b.For every graph and parameter quality q ten runs of each algorithm were performed.We show the mean approximation ratio of valid solutions produced by RQAOA and QIRO as plotted against the quality quantile q of the parameters used to generate the correlations by the QAOA.The parameters at q = 0 correspond to the fully optimized parameters.As expected, QIRO only produces valid solutions.Remarkably, if perfect (i.e., optimized) parameters are supplied to RQAOA, it likewise produces exclusively valid solutions.However, for both values of λ, there is a sharp drop-off in the number of valid solutions for RQAOA.Because finding optimal parameters in an actual experimental setting is challenging [46,47], the robustness of QIRO is a valuable feature, in particular for realistic experimental scenarios.Additionally, results in Figure 5b suggest that increasing the quality of the quantum correlations has a positive impact on the performance of both QIRO as well as RQAOA.This is in agreement with the outcome of our experiments at p > 1 in Appendix E, where a clear of improved QIRO performance with increasing p is observed.Surprisingly, results shown in Fig. 5b indicate that the 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 penalty term λ should favor valid solutions.In general, we observe complex behaviour in the (λ, q) parameter plane (see Fig. 10 and the associated discussion in Appendix F).We note that determining the magnitude of penalty terms 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 remains to be solved is a valid MIS problem.Moreover, if the starting graph is a UDG (as in Fig. 5a) then also all successive (simplified) graphs are -see also Fig. 2. Be-cause UDGs can be naturally embedded on neutral atom quantum processors [13,14], we can use such devices to produce the required quantum correlations.In what follows, we present results from the QuEra Aquila neutral atom quantum processor, as accessed through Amazon Braket.The details on the (quasi-)adiabatic protocols used to generate the correlations on the QuEra Aquila device can be found in Appendix D. Due to the limited accessibility of the device, we report the results on a perinstance basis.We sort the instances in increasing order of their 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.
In Fig. 6 we display the approximation ratios |IS|/|MIS| obtained by QIRO with correlations from the QuEra device (QIRO-QuEra) and classical simulations of p = 1 QAOA (QIRO-QAOA), respectively.On the QuEra device, we generate the correlations using the protocols outlined in Appendix D. We compare these to the sizes of the independent sets found by the minimal degree greedy algorithm.The reference maximum independent set is computed using the tensor network algorithm from Ref. [48].For clarity, results of RQAOA are omitted from the plot, because the approximation ratios found for UDGs were generally comparable to those found by QIRO-QAOA (as seen in Fig. 5b).
Our results indicate that the quantum correlations from the QuEra machine are better at guiding the QIRO procedure than numerical simulations of QAOA at p = 1.This is in line with previous results for solving the MIS problem using analog protocols on such hardware [14,49].A larger sample size would be required to solidify such claims.However, initial results appear to suggest that QIRO-QuEra performs comparably to the minimal degree greedy algorithm, reinforcing our belief that improved quantum correlations boost the performance of QIRO.Strikingly, the minimal degree greedy algorithm fails dramatically on some graphs (e.g., instance 1), where it achieves an approximation ratio of around ∼ 0.6, while QIRO-QuEra is able to solve the same problems (nearly) optimally.This preliminary evidence showcases the potential utility of quantum correlations for guiding optimization algorithms.

B. MAX-2-SAT
Next, we analyze the performance of QIRO and QIRO with backtracking (QIRO + BT) on random MAX-2-SAT instances with clauses of length two.We compare our results to those obtained with RQAOA and simulated annealing (SA) [50], an established and relatively simple classical heuristic.More detailed simulations with other, more powerful heuristics are left for future studies.Details about our implementation of SA can be found in Appendix G.
We first analyze the system size scaling of the algorithms' performance for random MAX-2-SAT instances 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 parallel tempering (PT) [34] to compute (approximately) optimal solutions for instances where the RC2 solver timedout.Details of our implementation of PT can be found in Appendix G.We then report the number of times the algorithms perform at least as well as PT or the RC2 solver.
Figure 7 showcases the results of our experiments.QIRO performs better than RQAOA for all clause-tovariable ratios α and over the entire range of the considered problem sizes.This shows that performance gains can be made by including very simple problem-specific update rules that require only modest classical resources.We note, however, that the improvement is especially prominent for lower values of α, and diminishes for problems with higher α (connectivity).Moreover, due to the inclusion of inference rules, the Boolean formula is simplified in fewer update steps, thus requiring fewer calls to the quantum device (and fewer optimizations of the QAOA circuit) than in RQAOA.
As expected, we notice a decrease in performance with increasing problem sizes, with both QIRO and RQAOA often struggling to find optimal solutions at the largest problem sizes considered.However, with the addition of backtracking, QIRO+BT is able to optimally solve more than half of the instances for almost all system sizes and all values of α considered here.For α = 2 it even slightly outperforms our SA implementation.It should be noted that while the number of calls to the quantum device increases significantly (i.e., from linear to quadratic in the problem size) if backtracking is included, the requirements on the quantum hardware in terms of the gate and qubit count are kept constant.This fundamentally differs from QAOA, where an improved performance can only be obtained by increasing the circuit depth p, and thus the quantum resources required.

V. DISCUSSION
In this work we have extended previous efforts on hybrid quantum-classical algorithms for combinatorial optimization, by introducing a family of quantum-informed recursive optimization (QIRO) algorithms.We propose the use of quantum-informed classical update steps to recursively simplify the original optimization problem, allowing for the incorporation of problem-specific classical subroutines (e.g., satisfiability inference rules) from the rich literature on classical combinatorial optimization [25].We have shown via numerical simulations for random instances of two NP-hard combinatorial optimization problems that such classical subroutines can significantly strengthen the algorithm's performance.Additionally, we compare the performance of QIRO to RQAOA, providing, to the best of our knowledge, the first systematic numerical study of RQAOA on problems that are not Z 2 -symmetric.Importantly, we think that QIRO offers a promising approach for solving optimization problems with hard constraints, which have previously proven challenging for quantum algorithms [10,11,52].Using the MIS problem as a guiding example, we have shown how classical update steps can be used to enforce hard constraints and only produce feasible solutions.This is true regardless of the quality of the quantum resources, endowing QIRO with the robustness essential for algorithms within the NISQ era.
Interestingly, we observe that keeping the quantum resources equal, RQAOA outperforms QIRO on the MIS problem.We believe that this is a consequence of the local nature of the QIRO update rules for MIS.This is corroborated by our finding that improved quantum correlations from the neutral atom device enhance QIRO's performance compared to when QAOA at p = 1 is used to generate the correlations.The intuitive explanation is that the quantum many-body evolution on a neutral atom quantum processor provides sufficient non-local effects to counteract the locality of the update steps.In fact, adiabatic protocols analogous to those used to generate correlations here have previously been shown to produce candidate solutions surpassing those produced by QAOA at depths higher than p = 1 [14].Taken together, our results may serve as yet another indication that locality is a limiting factor for (quantum) optimization algorithms [53].
Furthermore, we have confirmed the intuitive expectation that higher-quality quantum correlations should enhance the performance of the algorithm, as evidenced by means of extensive numerical simulations, as well as experiments on an actual neutral atom quantum device.This result suggests that advances in quantum hardware should lead to further performance gains.
We note that we performed no optimization of the protocols used to generate the correlations in our experi-ments on the QuEra Aquila neutral atom device.Previously, dramatic improvements in the quality of candidate solutions has been obtained by first optimizing the schedules [14,49].Moreover, the (unoptimized) protocols and device used in this study (and specified in Appendix D) have been shown to produce candidate solutions typically corresponding to either very small independent sets, or possessing multiple violations of the independence constraint [49].Thus, it is reasonable to expect that such schedule optimization schemes could greatly benefit the quantum state preparation step in QIRO.Moreover, the ability of QIRO to find good (feasible) solutions despite the limitations of the quantum hardware serves as yet another demonstration of its robustness.
In addition, we have introduced the use of backtracking to obtain improved results via additional calls to existing quantum hardware, rather than requiring deeper circuits, as is usually the case in quantum (optimization) algorithms.This scheme has led to promising results for the MAX-2-SAT problem, where QIRO enhanced with backtracking could optimally solve most random problem instances, even at comparatively large instance sizes.Because large instances are accessible in numerical simulations (at p = 1 of QAOA), we were able to sensibly compare the performance of QIRO(+BT) with widely-used classical heuristic solvers such as simulated annealing.
Finally, we would like to stress that our work should be seen as a first demonstration of a general template for a larger class of hybrid quantum-classical algorithms for combinatorial optimization.We firmly believe that hybrid techniques will be necessary to unlock the advantages of quantum devices.The modular nature of QIRO offers a plethora of ways in which the presented toolbox can be modified; we touch upon some potential extensions in the next section.

VI. OUTLOOK
Several modifications to the QIRO algorithm as presented here can be envisioned.Exploring update rules for optimization problems beyond those examined in this manuscript would be a compelling avenue for future research.Likewise, it would be intriguing to devise the update rules for weighted versions of the problems addressed here.Moreover, one could attempt to extend the information extracted from the quantum state beyond just the maximal correlation.This could be accomplished by a randomized strategy, with the selection probabilities determined as some function of the correlation matrix.Furthermore, devising update rules incorporating correlations beyond those between nearest neighbors in the interaction graph might be useful, especially in combination with deeper QAOA circuits or on analog devices with long-range interactions (e.g., Rydberg atom arrays).Additionally, one could exploit the locality of QAOA, and simultaneously round multiple correlations at each update step, if the corresponding variables are outside of each other's sphere of influence.Thus, one could save resources by requiring fewer calls to the quantum routines, and thus mitigate one of the limitations of the presented scheme.On the other hand, if deeper quantum circuits (e.g., QAOA at p > 1) are available, extracting correlations beyond those between variables adjacent in the interaction graph might be beneficial.Finally, considering problem-tailored alternatives to a brute-force search of the remaining solution space at a threshold problem size n c may significantly improve the performance and resource-efficiency of QIRO.
It is unreasonable to expect that the simple backtracking strategy presented here would be equally successful for larger problem instances, or instances from different problem classes.Hence, one might need to devise improved backtracking strategies, both in terms of the resources required, as well as in terms of the performance.One possibility is the use of reinforcement learning to guide the reduction process, as proposed by Patel et.al [33].Resources might be saved by efficient pruning of the search tree using existing methods [30], identifying and ignoring unfavorable parts of the solution space.Therefore, fewer quantum state preparations would be required.Finally, one could use backtracking to enforce constraints, by identifying steps where a constraint was violated, backtracking to it, and altering the decision made at that step.
It would likewise be interesting to consider alternative methods of generating the quantum correlations.This could include efficient numerical methods to obtain correlations with the QAOA at p > 1, by means of, e.g., tensor network simulations [54,55].Naturally, it would be interesting to acquire correlations from higher depth QAOA circuits executed on real digital quantum computers, or by generating them using quantum annealers based on superconducting flux qubits [56,57].One could consider optimizing analog protocols (using e.g., Bayesian optimization [49]) to obtain better correlations, and enhance the effectiveness of the algorithm.Moreover, modifying the cost function with respect to which the quantum resources are optimized might be helpful, similarly to the approach employed by Caha et al. [58].Finally, one could attempt to generate correlations using classical approximation algorithms, as proposed by Wagner et al. [59].Long story short, there is a lot to be done, but this paper is already long enough.
Conveniently, the modular nature of QIRO allows us to analyze the promise of quantum computers for combinatorial optimization, e.g., by comparing QIRO with different quantum and classical methods of computing the correlations.Should shortcomings of quantum approaches (e.g., locality) be identified in this manner, the classical subroutines in QIRO can be designed to rectify them.As such, we hope that our work can bring us closer to determining whether quantum devices can provide value for combinatorial optimization problems.
One may wonder why this should be a quantum algorithm of interest, if it is amenable to efficient classical simulations.The answer to this lies in the fact that the situation changes dramatically at higher depths of QAOA with p > 1 where no method for analytically computing the one-point and two-point correlations classically in polynomial time is known (to us).The same holds true at p = 1 when calculating the expectation values of k-body interaction terms, with k ≥ 3. Therefore, these cases could lead to an algorithm which provides quantum advantage.In general, we expect an increasing performance of QIRO for higher depths of QAOA, as indicated by the MIS experiments in Section IV A. We note that similar observations were previously made by Bravyi et al. [3].
To overcome the limitations of the QAOA due to its local structure [3,53], Bravyi et al. [3,7] introduced recursive QAOA (RQAOA).The algorithm optimizes a QAOA state based on the current problem Hamiltonian and uses this information to reduce the problem via the elimination of a single variable.Recursive application of this step reduces the problem until a specified threshold of remaining variables n c is reached.The remaining problem can be solved by a classical solver (e.g., by a brute-force search).The elimination step is motivated by the rounding of fractional correlations in relaxations of linear programs [3].While hardly a fair comparison in terms of resources required, RQAOA has been shown to significantly outperform QAOA [7,62].
It was first applied to the Z 2 -symmetric Max-Cut problem, which comprises exclusively 2-body terms.Here, we (trivially) extend the approach to also include singlebody terms.We additionally note that a generalization to Hamiltonians with an arbitrary degree of interaction can be envisaged.However, for Hamiltonians with Kbody terms with K ≥ 3 one faces an additional obstacle, as the degree of the interactions can increase when an elimination step is performed.This comes with additional requirements on the quantum hardware; as we are mostly interested in NISQ-era algorithms, we do not consider such extensions here.

8:
if i = j then ▷ One-point correlation.Here we provide the details on the update rules employed in the QIRO update step for MAX-2-SAT.Inference rules are a commonly used building block of many classical satisfiability solvers [31].They are used to simplify the formula by fixing the assignments of specific variables and therefore accelerate the remaining computation.Here, we use only very simple inference rules, that have low computational (and conceptual) complexity.The selected rules ensure that the optimal solutions of both the original and the simplified problem are equally good in quality.Despite their simplicity, the inference rules improve the performance of QIRO, both in terms of quality (see Section IV B), as well as in terms of resource efficiency.This is because fewer calls to the quantum device are required, due to the simplifications performed by the inference rules.The selected rules, listed next, were proposed in Refs.[30,31].
• The pure literal rule [32]: If a literal exclusively occurs with positive (negative) polarity, then the value of the corresponding variable is set to TRUE (FALSE).
• The almost common clause rule [32]: If a MAX-2-SAT instance includes the clauses x i ∨x j and xi ∨x j , then both clauses can be replaced with a singleliteral clause (i.e., unit clause) x j .If neither x i nor xi appears in the formula anymore, the variable can be set arbitrarily either to TRUE or FALSE.
• The complementary unit clause rule [63]: If the formula contains single-literal clauses x i and xi then these two clauses are removed.Again, if neither x i nor xi appears in the formula thereafter, the assignment of the variable can be set arbitrarily either to TRUE or FALSE.
• The dominating unit clause rule [63]: If the total count of clauses, regardless of their length, that contain a variable x i (x i ) is not greater than the count of unit clauses that contain xi (x i ), then the variable x i is set to FALSE (TRUE).
The complete QIRO update step Reduce, as introduced in Section III B 2 is as follows: Algorithm 5 Reduce (MAX-2-SAT) Input: Formula ϕ, partial solution S. Output Simpl.formula ϕ, updated partial solution S. 1: Prepare low-energy quantum state |ψ⟩.

6:
if sign(Mij) = 1 then 7: Replace xi → xj in ϕ and extend S accordingly.Replace xi → xj in ϕ and extend S accordingly.Importantly, the MAX-2-SAT QIRO update rules presented here preserve the problem structure in that only clauses with 2 or less variables appear in the formula.Consequently, the Hamiltonian remains in the same form as in Eq. (A1).

Appendix D: MIS on a neutral atom quantum device
Here, we provide details about the adiabatic protocols on a neutral atom quantum processor that we used to prepare the quantum states used to guide the QIRO algorithm for the MIS problem.The MIS problem on unit disk graphs (UDGs) naturally emerges in the context of neutral atom quantum processors based on Rydberg atom arrays [13,14].UDGs are a family of graphs, where nodes are connected only if the (Euclidian) distance between is smaller than some threshold radius R d .
Here, we generate such graphs by positioning nodes on an underlying square lattice with lattice constant a.We set √ 2a < R d < 2a, such that we get nearest neigbhor, and diagonal connectivity.We chose a lattice constant a = 5.3 µm, inspired by a previous study [49].The graphs presented here were generating by positioning 137 nodes on such an underlying lattice with 14×14 sites.This corresponds to a filling factor of roughly 70%.An example graph is shown in Figure 5a.
Furthermore, we classified the generated instances in terms of their classical hardness parameter where N M denotes the number of independent sets of size M .The hardness parameter was shown to be related to the performance of quantum and classical algorithms alike [14].However, we note that evidence has recently been put forward that HP only influences the performance of Markov chain Monte Carlo based algorithms (e.g., SA) [64].
Here, we consider graphs spanning several orders of magnitude in HP, which we computed using the tropical tensor network algorithm provided within the GenericTensorNetworks library [48,65].Specifically, the graphs considered in Fig. 6 have been sorted in order of increasing hardness, with the exact values of HP given in Table I.We now briefly summarize how neutral atom arrays can be used to solve the MIS problem on UDGs [13,14].The dynamics of neutral atom quantum devices are governed by the Hamiltonian Ĥ(t) = Ĥdr + Ĥcost , where the individual terms read: Here, Ω is the Rabi frequency, ∆ is the laser detuning, and the interaction term V ij = C 6 ∥x i − x j ∥ −6 2 , with C 6 the van der Waals coefficient.The particular value C 6 depends on the atom species used in the experimental setup.The number operator ni := 1 2 (1 2 − Ẑi ) counts the number of Rydberg excitations on the i-th site.
An essential property of the Hamiltonian in Eq. (D2) is the so-called Rydberg blockade phenomenon, in which two atoms cannot simultaneously be in the (excited) state |1⟩ if the distance between them is smaller than the Rydberg blockade radius R b ≡ (C 6 /ℏΩ) obey the independence constraint on UDGs if one assigns atoms that are in the |1⟩ to the MIS [13].Additionally, if we set 0 < ∆ < V ij , the ground state of H cost maximizes the number of excitations, while not violating independence, hence corresponding to solutions to the original MIS problem.We thus start the adiabatic protocol with ∆ i < 0, such that the ground state is |0⟩ ⊗|V | and then (slowly) increase the detuning until some final positive value, as shown in Fig. 8, such that the ground states of the final Hamiltonian correspond to the solutions of the MIS problem.
We ran such protocols on the QuEra Aquila device, which is available through the Amazon Braket service.The detailed parameters can be found in the caption of Fig. 8, and are partially motivated by a previous study [49].At each stage of the QIRO algorithm we estimate the correlations using a hundred measurements of the quantum state prepared by the adiabatic protocol described in Figure 8.

Appendix E: Higher depth QAOA experiments
To provide additional support to the claim that a higher quality of correlations yields better performance of QIRO, we performed numerical simulations of QIRO using correlations from QAOA at p > 1.In contrast to p = 1, no analytical formulae for the required correlators are known at p > 1, forcing us to resort to statevector simulations of QAOA circuits.Consequently, we had to limit our experiments to solving the MIS problem on smaller graph instances with 12 nodes, randomly chosen from the Erdős-Rényi ensemble as in Section IV A. For each fixed expected degree, we generated a hundred random instances, and found their MIS using a brute force search of the solution space.We then ran QIRO informed by correlations from the QAOA implementation in the PennyLane library [66], using the gradient descent optimizer with 300 iterations and 15 restarts of the optimizer at different random initial values of the variational parameters.For these experiments we set n c = 1 because of the comparatively smaller problem instance sizes.
The results in Fig. 9 show the fraction of optimal in- dependent 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 with increased QAOA depth, for both considered graph densities.Together with the results in Section IV A, this provides yet another indication that enhanced correlations lead to better performance of QIRO.Moreover, these experiments indicate that correlations from deeper circuits enable QIRO to surpass the performance of the greedy benchmark.Although this matches the intuitive expectations, further experiments on larger instances and quantum hardware are required to strengthen these claims.

Appendix F: Further RQAOA validity experiments
In Section IV A we assessed the performance of QIRO and RQAOA at different qualities of QAOA parameters.As the results presented in Fig. 5b exhibited interesting behavior, we here report an extended set of experiments, where RQAOA validities for additional magnitudes of the penalty term λ [cf.Eq. ( 1)] are shown.
Interestingly, the results partially defy the naïve expectation that a larger penalty term should favor enforcing constraints.This is especially prominent in the region between λ = 1.2 and λ = 5, where the robustness of RQAOA seems to decrease.However, the trend reverses again for λ = 10.
While a thorough analysis of the artefacts observed in Fig. 10 is beyond the scope of this paper, we offer some insights into what might underpin the peculiar behavior at intermediate values of λ.We suspect that this is caused by an intricate interplay between the deformations of the spectrum of the cost Hamiltonian Eq. ( 1) when λ is varied, with the particularities of QAOA states at p = 1.As λ is varied, only the energy levels of states that correspond to configurations with violations of the indepen-  5b, at each (q, λ) the same fifty UDGs of size 137 were generated and ten runs of RQAOA were performed.We here report the ratio valid (i.e., independent) sets found by RQAOA.Black stars indicate values of (q, λ) where all found solutions were feasible.
dence constraint are changed.Thus, the nature of the excited states can change dramatically, i.e., states corresponding to configurations with violations might shift relative to states corresponding to independent sets.It is difficult to predict the amplitudes of such excited states in the optimal QAOA states at p = 1, let alone at perturbed parameters.As such, a pattern as the one seen is not inconceivable.However, more experiments are required to pin down the exact underlying mechanism.For a MAX-2-SAT instance with n variables, starting from a random initial configuration, we perform simulated annealing with 600•n attempted variable flips.The annealing schedule consists of evenly spaced inverse temperatures β = T −1 , with an initial value β = 0 and a final inverse temperature of β = 6.
For parallel tempering (PT) [34][35][36] we initialize 12 replicas of the problem in random configurations with different temperatures.We chose the following fine-tuned the temperatures: These were chosen in accordance with the prescription that replica exchanges should be accepted with probabilities between 0.2 and 0.8 [67].We tested this empirically for a collection of representative MAX-2-SAT instances.We then used the following PT scheme.In each run of PT we carried out 15000 cycles, with one simulation cycle consisting of one Monte Carlo sweep (n single spin flips) per replica followed by an attempted replica exchange between all replica pairs with neighboring temperatures, starting from the lowest temperature.The configuration with the lowest energy found in any of the replicas during the course of the algorithm was finally returned as the solution.

1 FIG. 2 .
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 reduction rules for the case of a positive and negative two-point correlation.

FIG. 5 .
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 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).

FIG. 6 . 6 FIG. 7 .
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.

6 FIG. 9 .
FIG.9.The fraction of optimally solved instances by QIRO using QAOA at depths p ∈ {1, 2, 3}, and by the minimal degree greedy algorithm.The height of the bars show the mean fraction across ten runs of the algorithms, with the error bars indicating the best and worst performance across those runs.

FIG. 10 .
FIG.10.Ratio of valid solutions obtained by RQAOA for different qualities of the QAOA parameters q (at p = 1) and different penalty factors λ.As in the results in Fig.5b, at each (q, λ) the same fifty UDGs of size 137 were generated and ten runs of RQAOA were performed.We here report the ratio valid (i.e., independent) sets found by RQAOA.Black stars indicate values of (q, λ) where all found solutions were feasible.
Appendix G: Parallel tempering and simulated annealing setup Sn c ← Bruteforce( Ĥ) ▷ Solve the remaining problem.15: S ← Reconstruct(Sn c , L) ▷ Extract the solution from L.

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