# On Circuit Depth Scaling For Quantum Approximate Optimization

V. Akshay<EMAIL_ADDRESS>Skolkovo Institute of Science and
Technology, 30 Bolshoy Boulevard, Moscow 121205, Russian Federation  H.
Philathong  Skolkovo Institute of Science and Technology, 30 Bolshoy
Boulevard, Moscow 121205, Russian Federation  E. Campos  Skolkovo Institute of
Science and Technology, 30 Bolshoy Boulevard, Moscow 121205, Russian
Federation  D. Rabinovich  Skolkovo Institute of Science and Technology, 30
Bolshoy Boulevard, Moscow 121205, Russian Federation  I. Zacharov  Skolkovo
Institute of Science and Technology, 30 Bolshoy Boulevard, Moscow 121205,
Russian Federation  Xiao-Ming Zhang  Department of Physics, City University of
Hong Kong, Tat Chee Avenue, Kowloon, Hong Kong SAR, China  J.D. Biamonte
Skolkovo Institute of Science and Technology, 30 Bolshoy Boulevard, Moscow
121205, Russian Federation

(September 3, 2024)

###### Abstract

Variational quantum algorithms are the centerpiece of modern quantum
programming. These algorithms involve training parameterized quantum circuits
using a classical co-processor, an approach adapted partly from classical
machine learning. An important subclass of these algorithms, designed for
combinatorial optimization on currrent quantum hardware, is the quantum
approximate optimization algorithm (QAOA). Despite efforts to realize deeper
circuits, experimental state of the art implementations are limited to fixed
depth. However, it is known that problem density—a problem constraint to
variable ratio—induces under-parametrization in fixed depth QAOA. Density
dependent performance has been reported in the literature, yet the circuit
depth required to achieve fixed performance (henceforth called critical depth)
remained unknown. Here, we propose a predictive model, based on a logistic
saturation conjecture for critical depth scaling with respect to density.
Focusing on random instances of MAX-2-SAT, we test our predictive model
against simulated data with up to 15 qubits. We report the average critical
depth, required to attain a success probability of 0.7, saturates at a value
of 10 for densities beyond 4. We observe the predictive model to describe the
simulated data within a $3\sigma$ confidence interval. Furthermore, based on
the model, a linear trend for the critical depth with respect problem size is
recovered for the range of 5 to 15 qubits.

##  I Introduction

The quantum approximate optimization algorithm (QAOA) was introduced to solve
combinatorial optimization problems by Farhi et al. [1]. In this algorithm,
the candidate solution is prepared using a fixed ansatz structure with tunable
parameters that are variationally adjusted to minimize an objective function.
QAOA was inspired by adiabatic quantum computation [2, 3, 4], where a system
prepared in the ground state of an initial Hamiltonian evolves towards the
ground state of a problem Hamiltonian, given a slow enough evolution. The
circuit structure of QAOA shares resemblance to adiabatic evolution via a
finite depth Trotterization procedure, and thus allows for recovering the
exact ground state of the problem Hamiltoninan in the infinite depth limit. In
contrast, a fixed $p$–depth QAOA circuit, consists of alternating applications
of the complex exponentiation of the problem Hamiltonian and a so called mixer
Hamiltonian, accounting for $2p$ variational parameters.

Although QAOA recovers adiabatic evolution in the infinite depth limit, short
depth circuits have indeed been shown to offer some benefits in approximating
solutions to low density instances of combinatorial optimization problems [1,
5, 6]. Farhi et al. [1] studied the performance of QAOA applied to the MAX-CUT
problem on 3-regular graphs and reported that for $p=1$, QAOA guarantees a cut
that is at least $0.6924$ of the optimal cut. Furthermore, the quality of
approximation is shown to improve with increasing $p$. Beyond combinatorial
optimization problems, Farhi and Harrow [7] reported that the split operator
structure of QAOA might be a suitable candidate for establishing quantum
supremacy, since the output distributions of such circuits may not be
efficiently sampled classically. Moreover, in the works of Jiang et al. [8]
and Morales  et al. [9], the $O(\sqrt{N})$ query complexity of Grover’s
algorithm is recovered by a class of QAOA circuits, thus establishing quantum
advantage within the QAOA framework. These along with universality results
[10, 11, 12] make QAOA a promising model for NISQ era devices.

Although several advantages have been reported [1, 5, 6, 8, 9, 13, 14, 7, 15],
effects limiting the performance of QAOA have also been discovered [16, 17,
18]. Considering the energy error in approximation as a performance metric,
Akshay et al. [16, 17] demonstrated reachability deficits, an effect limiting
performance for fixed depth QAOA on high density problem instances. Deeper
circuits therefore become a necessity and hence come with overheads in
classical outerloop optimization. While approaches such as layerwise training
aim to reduce this classical computational cost, Campos et al. [18] showed
that such a strategy would not improve performance beyond some threshold
depth. Other approaches aimed at reducing the classical computation cost
exploit what are known as concentration effects [19, 20, 21, 22]. Although
these effects allows one to heuristically guess near-optimal parameters at
fixed depth and for increasing number of qubits, they fail to address the
depth required to guarantee fixed performance.

Here we propose a logistic saturation conjecture which states that the minimum
depth required for QAOA to achieve fixed performance up to some tolerance
(henceforth called critical depth) scales logistically with problem density.
We numerically test this conjecture for the case of MAX-2-SAT with up to $15$
qubits and maximum density of $4$, and observe (i) the predictive model
describes the simulated data within a $3\sigma$ confidence interval and (ii)
the critical depth scales linearly with problem size in this dataset.

##  II Background

###  II.1 Quantum Approximate Optimization Algorithm

The objective function of a combinatorial optimization problem is prescribed
as a map from the assignment of variables to non-negative real numbers
$\mathcal{C}:\\{0,1\\}^{\times n}\rightarrow\mathbb{R_{+}}$. Under a vector
space embedding [23, 24], the objective function then induces a corresponding
problem Hamiltonian $\mathcal{H}_{c}$, which encodes the solutions to the
problem instance in its ground state space,

$\mathcal{H}_{c}=\sum_{z\in\\{0,1\\}^{\times
n}}\mathcal{C}(z)\outerproduct{z}{z},$ (1)

where $\ket{z}$ are the computational basis vectors in the $2^{n}$ dimensional
Hilbert space of $n$ qubits.

The optimization problem aims to find assignments such that the objective
function $\mathcal{C}(z)$ is minimized. This is equivalent to finding a state
$\ket{z^{*}}$ such that

$\ket{z^{*}}\in\text{span}\\{\underset{z\in\\{0,1\\}^{\times
n}}{\operatorname*{arg\,min}}\bra{z}\mathcal{H}_{c}\ket{z}\\}.$ (2)

Measuring the state $\ket{z^{*}}$ in the computational basis outputs a
bitstring which is a solution to the optimization problem, $\mathcal{C}$.

To approximate solutions to an optimization problem instance using a $p$-depth
QAOA circuit, we generate ansatz states
$\ket{\psi(\boldsymbol{\gamma},\boldsymbol{\beta})}$ on a quantum computer
with tunable real parameters
$\boldsymbol{\gamma},\boldsymbol{\beta}\in[0,2\pi)^{\times p},[0,\pi)^{\times
p}$. These states are prepared by applying an alternating sequence of $2p$
parameterized unitary gates on the initial state
$\ket{+}^{\otimes{n}}=\frac{1}{\sqrt{2^{n}}}\sum_{z\in\\{0,1\\}^{\times
n}}\ket{z}$ as,

$\ket{\psi(\boldsymbol{\gamma},\boldsymbol{\beta})}=\prod_{k=1}^{p}e^{-i\beta_{k}\mathcal{H}_{x}}\cdot
e^{-i\gamma_{k}\mathcal{H}_{c}}\ket{+}^{\otimes{n}},$ (3)

with the mixing Hamiltonian, $\mathcal{H}_{x}=\sum_{l=1}^{n}X_{l}$. Here
$X_{l}$ is the Pauli matrix applied to the $l^{th}$ qubit.

The ansatz state in Eq. (3) is iteratively adjusted via an outer-loop
optimization routine to obtain an energy approximation as:

$\min_{\boldsymbol{\gamma},\boldsymbol{\beta}}\leavevmode\nobreak\
\bra{\psi(\boldsymbol{\gamma},\boldsymbol{\beta})}\mathcal{H}_{c}\ket{\psi(\boldsymbol{\gamma},\boldsymbol{\beta})}\geq\min\left(\mathcal{H}_{c}\right).$
(4)

Algorithmic performance of QAOA can be assessed under performance metrics such
as (i) energy error in approximation and (ii) ground state overlap [16, 17,
18]. Both these quantifiers feature shortcomings as energy approximations fail
to indicate closeness to optimal solutions while computing ground state
overlap requires prior knowledge of optimal solutions. However, it can be
demonstrated that both metrics inter-relate and therefore have some advantages
as a theoretical tool for studying algorithmic performance.

Let $\mathcal{H}^{\dagger}=\mathcal{H}\geq 0$ be an optimization problem
instance with $d$-degenerate ground state energy $\lambda_{0}$. In its
eigenbasis $\ket{\phi_{j}}$, we can write
$\mathcal{H}=\lambda_{0}\sum_{j=1}^{d}\outerproduct{\phi_{j}}{\phi_{j}}+\sum_{k>d}\lambda_{k}\outerproduct{\phi_{k}}{\phi_{k}}$,
with subsequent excited state energies,
$\lambda_{0}<\lambda_{k}\leq\lambda_{max}$ for $k>d$, separated by a gap
$\Delta=\min_{k>d}|\lambda_{0}-\lambda_{k}|$. For any arbitrary state
$\ket{\psi}$, the ground state overlap is calculated as:

$g(\psi)=\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}.$ (5)

It is straightforward now to derive (see Appendix B) the bounds relating the
energy error in approximations $\bra{\psi}\mathcal{H}\ket{\psi}-\lambda_{0}$,
and the ground state overlap $\text{g}(\psi)$:

$1-\frac{\bra{\psi}\mathcal{H}\ket{\psi}-\lambda_{0}}{\Delta}\leq g(\psi)\leq
1-\frac{\bra{\psi}\mathcal{H}\ket{\psi}-\lambda_{0}}{\lambda_{max}}.$ (6)

Note that the relation works when the error in approximation is less than the
spectral gap $\Delta$. These bounds are the extended version of the
variational stability lemma in [10] to the degenerate case and allows one to
recover an estimate on the ground state overlap from the measured energy.

###  II.2 Satisfiability

Boolean satisfiability or SAT, is the problem of determining satisfiability of
Boolean formulae. More specifically, one aims to decide whether a given
formula can be satisfied by assigning truth values to the variables. Limited
to $k$-variables or literals per clause, and expressed in the conjunctive
normal form (CNF), $k$-SAT is NP-complete for $k\geq 3$ [25]. $k$-SAT clauses
are randomly generated to form instances by uniformly selecting unique
$k$-tuples from the union of a variable set (cardinality $n>k$) and its
element wise negation. Generated random instances admit an order parameter
called clause or problem density given by,

$\alpha=m/n,$ (7)

where $m$ is the number of clauses and $n$ is the number of Boolean variables
in the $k$-SAT instance.

The clause density plays an important role in the study of satisfiability
phase transitions. It has been observed that for random $k$-SAT instances, an
abrupt change in satisfiability occurs at some critical clause density,
$\alpha_{c}$. Though this criticality has been studied only empirically for
$3$-SAT ($\alpha_{c}\approx 4.27$), it is proven for $2$-SAT [26].
Furthermore, it is observed that the computational resources required for SAT
solvers increase at the transition point with an easy-hard-easy pattern. This
behaviour suggest that the hard instances concentrate around the critical
clause density [27, 28, 29, 30, 31, 26].

In this work we focus on the optimization version of $k$-SAT, otherwise called
MAX-$k$-SAT. MAX-$k$-SAT is the canonical NP-hard optimization problem (for
$k\geq 2$) [32] where one seeks to determine variable assignments that
maximize the number of satisfied clauses in a given instance. Similar to
decision $k$-SAT, MAX-$k$-SAT also features a phase transition with an easy-
hard pattern [33]. Although the point of criticality differs for the
optimization and decision variants of these problems, they coincide for $k=2$
with $\alpha_{c}=1$ [34].

With standard Ising embedding techniques [23, 24], MAX-2-SAT instances can be
mapped onto $2$-local Ising Hamiltonians of the form:

$\mathcal{H}=\sum_{i<j}\mathcal{J}_{ij}Z_{i}Z_{j}+\sum_{i}h_{i}Z_{i},$ (8)

with appropriate coefficients and scale such that: (i) for satisfiable
assignments, $\ket{z}$, $\bra{z}\mathcal{H}\ket{z}=0$ and (ii) for un-
satisfiable assignments $\bra{z}\mathcal{H}\ket{z}=m^{\prime}$, where
$m^{\prime}$ is the number of clauses violated by the assignment. A detailed
construction of such Hamiltonians is shown in Appendix A. By construction, the
ground state space of $\mathcal{H}$ is spanned by assignments that violate the
least number of clauses in a given MAX-$2$-SAT instance.

##  III Logistic Saturation Conjecture

The state space which can be prepared with a $p$-depth QAOA circuit is defined
as,

$\Omega=\underset{\boldsymbol{\gamma},\boldsymbol{\beta}}{\bigcup}\ket{\psi(\boldsymbol{\gamma},\boldsymbol{\beta})},$
(9)

with $\ket{\psi(\boldsymbol{\gamma},\boldsymbol{\beta})}$ given by Eq. (3). In
general, for a given instance size $n$, density $\alpha$ and depth $p$, the
variational state space may not cover the whole Hilbert space. Therefore,
optimization over $\Omega$ may not converge to the exact ground state energy.
The difference, or the energy error in approximation,

$\text{f}(p,\alpha,n)=\min_{\psi\in\Omega\subseteq\mathbb{C}_{2}^{\otimes
n}}\bra{\psi}\mathcal{H}\ket{\psi}-\min_{\phi\in\mathbb{C}_{2}^{\otimes
n}}\bra{\phi}\mathcal{H}\ket{\phi},$ (10)

quantifies performance of $p$-depth QAOA. Whenever $f(p,\alpha,n)>0$, QAOA
suffers a reachability deficit. This behaviour has been studied in [16, 17],
and it is observed that the onset of these deficits correspond to the phase
transition criticality of the considered optimization problem. For the case of
MAX-2-SAT, it was shown that $f(p,\alpha,n)$ grows with increasing clause
density beyond $\alpha_{c}\approx 1$. Though better performance is achieved at
the cost of increasing QAOA depth, a similar trend for deficits is reflected.

Motivated by these findings, we propose a logistic saturation conjecture:

Given some fixed tolerance on performance, $\epsilon>0$. For QAOA on instances
characterized by density $\alpha$, the critical depth $p^{*}$, required to
guarantee approximations that fall within an $\epsilon$ tolerance is given by:

$p^{*}(\alpha)=\min\\{p\leavevmode\nobreak\ |\leavevmode\nobreak\
\text{f}(p,\alpha,n)\leq\epsilon\\}.$ (11)

We conjecture the critical depth to scale with density as,

$p^{*}(\alpha)\approx\frac{p_{max}}{1+e^{-\kappa(\alpha-\alpha_{c})}},$ (12)

where $\alpha_{c}$ is the critical density, $\kappa$ the logistic growth rate,
and $p_{max}$ the saturation value.

##  IV Results

We first recover the findings of [16, 17] in Figure 1 by studying the
performance of fixed depth QAOA on $300$ randomly generated MAX-$2$-SAT
instances on 15 variables and clauses ranging from 1 to 60
($\alpha\in\left[0.06,4.0\right]$).

Figure 1: Average energy error in approximation, $f$ in Eq. (10), vs clause
density for MAX-$2$-SAT instances on $15$ variables. Data points represent the
average value obtained over statistics of $300$ random instances with error
bars representing $3\sigma$ for the estimated mean. Colors indicate varying
QAOA depths $p=\\{5,10,15\\}$, illustrating improved performance at higher
depths.

Numerically, it is seen that the average performance worsens monotonically
with increasing clause density. This is also evident (see Fig. 2) for the case
of ground state overlap as described by stability bounds in Eq. (6).

Figure 2:  Ground state overlap vs clause density for depth $p=10$ QAOA on
MAX-2-SAT instances with 15 variables. Stability lower bound is calculated by
considering average energy error in approximation obtained over statistics of
$300$ random instances. Points in orange indicate the corresponding average
ground state overlap.

It becomes apparent that for guaranteeing fixed performance, the QAOA depth
must necessarily increase with increasing clause density.

To test the logistic saturation conjecture, we numerically calculate the depth
required to attain performance within an $\epsilon$ tolerance. This is done by
incrementing depth, until the energy error in approximations falls within
acceptance, $\leavevmode\nobreak\ \text{f}(p,\alpha,n)\leq\epsilon$. For each
instance this process is repeated and the average critical depth is
calculated. Figure 3 illustrates scaling of the average critical depth with
respect to clause density and compares the numerical data against the
predictive model in Eq. (12). We see that the model is able to match the
numerical data well within a $3\sigma$ confidence interval for the estimated
mean.

We further investigate the reliability of the predictive model by treating the
critical density as a fitting parameter in Eq. (12) and compare against the
theoretical value of $\alpha_{c}\approx 1$. In Figure 4 fitting parameters are
recovered for $n=5$ up to $15$ and indeed we see the critical density
recovered from the predictive model matches the theoretical value up to finite
size variations.

Figure 3: Average critical depth $p^{*}$ vs clause density. Data points
represent the average $p^{*}$ calculated according to Eq. (11) across problem
density for MAX-$2$-SAT with tolerance $\epsilon=0.3$. Colors correspond to
different problem sizes $n=\\{5,10,15\\}$ with errorbars representing
$3\sigma$ for the estimated mean. Solid curves represent the least squares fit
of the data to the predictive model as described in Eq. (12).  Figure 4:
Fitting parameters vs number of qubits for QAOA on MAX-$2$-SAT with tolerance
$\epsilon=0.3$. Note that the recovered critical density $\alpha_{c}\approx
1$, and growth rate $\kappa$, show very little variability with respect to
number of qubits. In contrast, the the saturation value $p_{max}$ illustrates
a linear trend in the considered range of $5$ to $15$ qubits.

##  V Numerical Details

The computationally intensive calculations were performed on the KunLun
supercomputer with a maximum CPU frequency of $2.1$ GHz. The code was written
in Python and C++. Python code is responsible for providing easy-to-use
interfaces and C++ code is responsible for the actual simulation tasks as the
underlying simulation code. The code executed in parallel with a resource of
1500 cores. The calculations were performed for 300 random instances for each
respective problem density.

We used a heuristic optimisation strategy, motivated by layerwise training, to
find good initial parameterization for the ansatz such that gradient based
classical optimizers may avoid most local minima.

In order to find optimal parameters for depth $p$, we start by optimizing a
depth one ansatz (single layer). Then keeping the previous parameters fixed, a
new layer is added and trained. This is followed by a simultaneous
optimization of all parameters and the process is repeated until the desired
depth $p$ is attained.

In our implementation of this heuristic, at each optimization step, the best
result out of 25 random seed runs of L-BFGS-B optimizer is chosen.

##  VI Conclusion

QAOA is among the most studied gate based approaches for combinatorial
optimization in NISQ era devices. Although analytical techniques address QAOA
either at low depth or in the adiabatic limit, $p\to\infty$, understanding
performance at intermediate depths remains largely open. Furthermore, problem
density induces underparameterization and limits the performance of fixed
depth QAOA. An open question in this regard is the required depth for QAOA to
succeed. Additionally, by treating the circuit depth of QAOA as a
computational resource, it becomes paramount in complexity studies to know the
scaling of depth with respect to problem size. Understanding these open
problems paves way towards identifying quantum advantage in QAOA.

In this paper we introduce a methodology based on a predictive model that aims
to address these open problems. We first relate QAOA energy approximations and
success probability (or ground state overlap) via Eq. (6). This relation then
determines success/failure of the algorithm up to some $\epsilon$ tolerance.
Based on this acceptance criterion, we numerically recovered the critical
depth needed for QAOA to achieve success probability of at least 0.7 for the
case of MAX-$2$-SAT. Considering up to $15$ variables and clause density up to
4 we indeed observe that our proposed model describes the data within a
$3\sigma$ confidence interval. Furthermore, we recovered for the first time,
the scaling of critical depth with respect to problem size which illustrates a
linear trend in our considered range of problem sizes. Though this finite
range is insufficient to assert quantum advantage, we anticipate future works
to test the accuracy of the presented model over comprehensive ranges of
densities and problem sizes.

###### Acknowledgements.

Author contributions. All authors conceived and developed the theory and
design of this study and verified the methods. All authors wrote the paper.
Data and code availability. The data that supports the findings of this study
are available within the article. The code for generating the data will be
available on reasonable request. Competing interests. The authors declare no
competing interests.

## References

  * [1] Edward Farhi, Jeffrey Goldstone, and Sam Gutmann.  A quantum approximate optimization algorithm.  arXiv preprint arXiv:1411.4028, 2014.
  * [2] Tadashi Kadowaki and Hidetoshi Nishimori.  Quantum annealing in the transverse ising model.  Physical Review E, 58(5):5355, 1998.
  * [3] Edward Farhi, Jeffrey Goldstone, Sam Gutmann, Joshua Lapan, Andrew Lundgren, and Daniel Preda.  A quantum adiabatic evolution algorithm applied to random instances of an np-complete problem.  Science, 292(5516):472–475, 2001.
  * [4] Sergio Boixo, Troels F Rønnow, Sergei V Isakov, Zhihui Wang, David Wecker, Daniel A Lidar, John M Martinis, and Matthias Troyer.  Evidence for quantum annealing with more than one hundred qubits.  Nature physics, 10(3):218–224, 2014.
  * [5] Edward Farhi, Jeffrey Goldstone, and Sam Gutmann.  A quantum approximate optimization algorithm applied to a bounded occurrence constraint problem.  arXiv preprint arXiv:1412.6062, 2014.
  * [6] Zhihui Wang, Stuart Hadfield, Zhang Jiang, and Eleanor G Rieffel.  Quantum approximate optimization algorithm for maxcut: A fermionic view.  Physical Review A, 97(2):022304, 2018.
  * [7] Edward Farhi and Aram W Harrow.  Quantum supremacy through the quantum approximate optimization algorithm.  arXiv preprint arXiv:1602.07674, 2016.
  * [8] Zhang Jiang, Eleanor G Rieffel, and Zhihui Wang.  Near-optimal quantum circuit for grover’s unstructured search using a transverse field.  Physical Review A, 95(6):062317, 2017.
  * [9] Mauro ES Morales, Timur Tlyachev, and Jacob Biamonte.  Variational learning of grover’s quantum search algorithm.  Physical Review A, 98(6):062333, 2018.
  * [10] Jacob Biamonte.  Universal variational quantum computation.  Physical Review A, 103(3):L030401, 2021.
  * [11] Mauro ES Morales, JD Biamonte, and Zoltán Zimborás.  On the universality of the quantum approximate optimization algorithm.  Quantum Information Processing, 19(9):1–26, 2020.
  * [12] Seth Lloyd.  Quantum approximate optimization is computationally universal.  arXiv preprint arXiv:1812.11075, 2018.
  * [13] Eric Anschuetz, Jonathan Olson, Alán Aspuru-Guzik, and Yudong Cao.  Variational quantum factoring.  In International Workshop on Quantum Technology and Optimization Problems, pages 74–85. Springer, 2019.
  * [14] Wen Wei Ho and Timothy H Hsieh.  Efficient variational simulation of non-trivial quantum states.  SciPost Phys, 6:29, 2019.
  * [15] Zhi-Cheng Yang, Armin Rahmani, Alireza Shabani, Hartmut Neven, and Claudio Chamon.  Optimizing variational quantum algorithms using pontryagin’s minimum principle.  Physical Review X, 7(2):021027, 2017.
  * [16] V Akshay, H Philathong, Mauro ES Morales, and Jacob D Biamonte.  Reachability deficits in quantum approximate optimization.  Physical Review Letters, 124(9):090504, 2020.
  * [17] V Akshay, H Philathong, I Zacharov, and J Biamonte.  Reachability deficits in quantum approximate optimization of graph problems.  Quantum, 5:532, 2021.
  * [18] E Campos, D Rabinovich, V Akshay, and J Biamonte.  Training saturation in layerwise quantum approximate optimization.  Physical Review A, 104(3):L030401, 2021.
  * [19] Edward Farhi, Jeffrey Goldstone, Sam Gutmann, and Leo Zhou.  The quantum approximate optimization algorithm and the sherrington-kirkpatrick model at infinite size.  arXiv preprint arXiv:1910.08187, Oct 2019.
  * [20] Michael Streif and Martin Leib.  Training the quantum approximate optimization algorithm without access to a quantum processing unit.  Quantum Science and Technology, 5(3):34008, Jul 2020.
  * [21] V. Akshay, D. Rabinovich, E. Campos, and J. Biamonte.  Parameter concentrations in quantum approximate optimization.  Physical Review A, 104:L010401, Jul 2021.
  * [22] Dave Wecker, Matthew B Hastings, and Matthias Troyer.  Training a quantum optimizer.  Physical Review A, 94(2):022309, 2016.
  * [23] JD Biamonte.  Nonperturbative k-body to two-body commuting conversion hamiltonians and embedding problem instances into ising spins.  Physical Review A, 77(5):052331, 2008.
  * [24] Andrew Lucas.  Ising formulations of many np problems.  Frontiers in Physics, 2:5, 2014.
  * [25] Stephen A Cook.  The complexity of theorem-proving procedures.  In Proceedings of the third annual ACM symposium on Theory of computing, pages 151–158, 1971.
  * [26] Andreas Goerdt.  A threshold for unsatisfiability.  Journal of Computer and System Sciences, 53(3):469–486, 1996.
  * [27] James M Crawford and Larry D Auton.  Experimental results on the crossover point in random 3-SAT.  In Artificial Intelligence, volume 81, pages 31–57, 1996.
  * [28] Ehud Friedgut, Jean Bourgain, et al.  Sharp thresholds of graph properties, and the k-SAT problem.  Journal of the American mathematical Society, 12(4):1017–1054, 1999\.
  * [29] Bart Selman and Scott Kirkpatrick.  Critical behavior in the computational cost of satisfiability testing.  Artificial Intelligence, 81(1-2):273–295, 1996.
  * [30] Bart Selman, David G Mitchell, and Hector J Levesque.  Generating hard satisfiability problems.  Artificial intelligence, 81(1-2):17–29, 1996.
  * [31] Vasek Chvátal and Bruce Reed.  Mick gets some (the odds are on his side)(satisfiability).  In Proceedings., 33rd Annual Symposium on Foundations of Computer Science, pages 620–627. IEEE Computer Society, 1992.
  * [32] Francisco Barahona.  On the computational complexity of ising spin glass models.  Journal of Physics A: Mathematical and General, 15(10):3241, 1982\.
  * [33] Weixiong Zhang.  Phase transitions and backbones of 3-sat and maximum 3-sat.  In International Conference on Principles and Practice of Constraint Programming, pages 153–167. Springer, 2001.
  * [34] Don Coppersmith, David Gamarnik, MohammadTaghi Hajiaghayi, and Gregory B Sorkin.  Random max sat, random max cut, and their phase transitions.  Random Structures & Algorithms, 24(4):502–545, 2004.

##  Appendix A Ising Embedding of SAT instances

We start by mapping logical bits $\\{0,1\\}$ to states
$\\{\ket{0},\ket{1}\\}\in\mathbb{C}_{2}$ respectively. The objective function
to be minimized is equivalent to minimizing a Hamiltonian $\mathcal{H}_{SAT}$,
constructed by the real linear extension of $\\{P_{0},P_{1},\openone\\}$ as:

$x_{j}\longrightarrow P_{j}^{0},\qquad\neg x_{j}\longrightarrow P_{j}^{1},$
(13)

and

$\wedge\longrightarrow+,\qquad\vee\longrightarrow\otimes.$ (14)

Here $P_{j}^{1}=\outerproduct{1}{1}_{j}$ and
$P_{j}^{0}=\outerproduct{0}{0}_{j}$ are rank one projectors acting on the
$j^{th}$ qubit. The construction ensures each clause in a given SAT instance
is mapped onto projectors, $P_{ij\cdots
k}^{\alpha\beta\cdots\gamma}=P_{i}^{\alpha}\otimes
P_{j}^{\beta}\otimes\cdots\otimes P_{k}^{\gamma}$, which penalize
unsatisfiable assignments $\ket{x}\in\mathbb{C}_{2}^{\otimes n}$ with a
penalty $\bra{x}P_{ij\cdots k}^{\alpha\beta\cdots\gamma}\ket{x}=1$.

The full Hamiltonian $\mathcal{H}_{SAT}$ is then constructed by summing over
all the clauses in an instance,

$\mathcal{H}_{SAT}=\sum_{l=1}^{m}\mathcal{C}_{l}\\{P_{ij...k}^{\alpha\beta...\gamma}\\},$
(15)

where $\mathcal{C}_{l}$ assigns the value of
$\alpha,\beta,\cdots,\gamma\in\\{0,1\\}$, corresponding to the negation of the
literals indexed by $i,j,\cdots,k\in\\{1,2,\cdots,n\\}$ in the $l^{th}$
clause. It is straightforward to see that the ground state space of this SAT-
Hamiltonian is spanned by assignments that violate minimal number of clauses.

Decomposing the projectors onto the $Z$ Pauli basis as,

$P^{\alpha}_{j}=\frac{1}{2}(\openone+(-1)^{\alpha}Z_{j}),$ (16)

Eq. (15) can be re-written in the form of a generalized Ising Hamiltonian with
at most k-body interaction. In case of 2-SAT, the Hamiltonian
$\mathcal{H}_{2SAT}$, contain at most 2-body interactions as,

$\begin{split}\mathcal{H}_{2SAT}&=\sum_{l=1}^{m}\mathcal{C}_{l}\\{P^{\alpha}_{j}\otimes
P^{\beta}_{k}\\}\\\
&=\frac{1}{4}\left(m\openone+\sum_{j}h_{j}Z{j}+\sum_{j<k}\mathcal{J}_{jk}Z_{j}Z_{k}\right),\end{split}$
(17)

with appropriate coefficients $h_{j}$ and $\mathcal{J}_{jk}$.

##  Appendix B Derivation of Energy-Overlap Bounds.

Once a variatinal state $\ket{\psi(\boldsymbol{\theta})}$ has been prepared on
a quantum computer as a function of a set of tunable parameters
$\boldsymbol{\theta}$, we would like to know how close this state is to the
solution space of the problem being solved. We will adapt and extend the Lemma
1 (Variational Stability) from [10] to the degenerate case.

Let $\mathcal{H}\geq 0$ be some Hamiltonian with $d$–degenerate ground state
energy $\lambda_{0}\geq 0$. The Hamiltonian in its eigenbasis can be expressed
as,

$\mathcal{H}=\lambda_{0}\sum_{j=1}^{d}\outerproduct{\phi_{j}}{\phi_{j}}+\sum_{k>d}\lambda_{k}\outerproduct{\phi_{k}}{\phi_{k}}.$
(18)

Here we reorder the eigenbasis so that
$\lambda_{0}<\lambda_{1}\leq\cdots\lambda_{k}\cdots\leq\lambda_{max}$.
Consider the expectation value of the Hamiltonian on an arbitrary state
$\ket{\psi}$,

$\bra{\psi}\mathcal{H}\ket{\psi}=\lambda_{0}\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}+\sum_{k>d}\lambda_{k}\absolutevalue{\innerproduct{\phi_{k}}{\psi}}^{2}.$
(19)

Assuming the gap $\Delta=\lambda_{1}-\lambda_{0}$, is known,

$\begin{split}\bra{\psi}\mathcal{H}\ket{\psi}&=\lambda_{0}\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}+\sum_{k>d}\lambda_{k}\absolutevalue{\innerproduct{\phi_{k}}{\psi}}^{2}\\\
&\geq\lambda_{0}(\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2})+(\lambda_{0}+\Delta)(\sum_{k>d}\absolutevalue{\innerproduct{\phi_{k}}{\psi}}^{2})\\\
&=\lambda_{0}(\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2})+(\lambda_{0}+\Delta)(1-\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2})\\\
&=\lambda_{0}+\Delta(1-\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}).\end{split}$

Rearranging the terms then gives the required lower bound on the ground state
as:

$\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}\geq
1-\frac{\bra{\psi}\mathcal{H}\ket{\psi}-\lambda_{0}}{\Delta}.$ (20)

For the upper bound, we substitute for each excited state energies
$\lambda_{k}$ in Eq. (19) by $\lambda_{0}+\lambda_{max}$,

$\begin{split}\bra{\psi}\mathcal{H}\ket{\psi}&\leq\lambda_{0}\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}+(\lambda_{0}+\lambda_{max})\sum_{k>d}\absolutevalue{\innerproduct{\phi_{k}}{\psi}}^{2}\\\
&\leq\lambda_{0}(\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2})\\\
&\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\
+(\lambda_{0}+\lambda_{max})(1-\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}).\\\
\end{split}$ (21)

Rearranging the terms now gives the required upper bound as:

$\sum_{j=1}^{d}\absolutevalue{\innerproduct{\phi_{j}}{\psi}}^{2}\leq
1-\frac{\bra{\psi}\mathcal{H}\ket{\psi}-\lambda_{0}}{\lambda_{max}}.$ (22)

By combining both bounds, we obtain the stability theorem that relates energy
error in approximation with the success probability (or the ground state
overlap) for the degenerated case:

$1-\frac{f}{\Delta}\leq\sum_{i=1}^{d}\absolutevalue{\innerproduct{\phi_{i}}{\psi}}^{2}\leq
1-\frac{f}{\lambda_{max}},$ (23)

where $f=\bra{\psi}\mathcal{H}\ket{\psi}-\lambda_{0}$.

Notice that for the case of MAX-$2$-SAT instances, $\lambda_{max}\leq m$ where
$m$ the number of clauses in the instance and $\Delta\geq 1$.