Fault-tolerant resource estimation using graph-state compilation on a modular superconducting architecture

Due to the size of a utility-scale input circuit, a decomposition strategy is required to avoid memory bottlenecks on RRE. Fortunately, quantum algorithms often comprise repeating blocks, which facilitate decomposing them.

Our decomposition strategy is to slice $U$ in the time direction into $n_{\text{widgets}}$ “widgets” (or subcircuits) of $n_{\text{input}}$-width with a user-defined maximum allowed depth of $d_{\text{max}}$. We emphasize there are no entangling gates connecting the subcircuits. One could then execute the complete algorithm as a sequence of the subcircuits, $U=[U_{0},\dots,U_{n_{\text{widgets}}-1}]$ where

where we have defined $\ket{\text{state}}_{0}\equiv\ket{\text{input}}$ and $\ket{\text{state}}_{n_{\text{widgets}}-1}\equiv\ket{\text{output}}$. We call this process “circuit widgetization”, pictured in fig. 1, which allows us to avoid the memory overheads associated with decomposing the complete circuit to elementary operations. Notice the particular choice of widgetization above is likely not optimized for reducing space-time costs; however, it was selected due to the ease of implementation.

After the widgetization process, one can construct an approximate version of $S_{\text{est}}$. We propose and verify a widgetization-based set $\tilde{S}_{\text{est}}$ below as a good approximation of $S_{\text{est}}$ for resource estimation purposes,

In this widgetization-based estimation procedure, we first obtain estimation primitives for each widget based on its graph and then define the full circuit resources estimates as,

where $g_{i}$ is the graph-state compiled from the $i$th widget. eq. 6 highlights the critical assumptions we made while generating all the estimations in this work and formally states our “stitching” strategy. We stitch the individual graphs and scheduling lists to approximate the complete $\mathcal{G}$, $S^{\mathcal{G}}_{\text{prep}}$, $S^{\mathcal{G}}_{\text{consump}}$, and $S^{\mathcal{G}}_{\text{meas}}$. Note that while literal stitching of the widgets would be required to execute the circuit on a MBQC, it is not required for resource estimation purposes [25]. The computational requirements required for obtaining resource estimates can thus be reduced by calculating the “widgetized" quantities above on the software.

To prepare for graph compilation each widget must be parsed and transpiled so that it only contains logical Clifford+$T$+$R_{z}(\theta)$ using exact gate equivalences. Here, $\theta$’s are nontrivial small angles. Later on, when we want to perform single qubit logical measurements, at the end during the distillation-consumption stage, each $R_{z}(\theta)$ node needs to be decomposed to a chain of Clifford+$T$ itself – a process sometimes referred to as gate synthesis[28, 29].

The transpilation to logical Clifford+$T$+$R_{z}(\theta)$ is also performed by RRE. Performing gate synthesis at the end allows us to generate, cache, and reuse subgraphs $\{g_{0},\dots,g_{n_{\text{widgets}}-1}\}$ independent from logical and hardware-level configurations, $\{\text{config}\}$.

Once transpilation is complete, we need to generate the unified sets of graph states $\{g\}$ and schedules $\{S^{g}_{\text{prep}}\}$, $\{S^{g}_{\text{consump}}\}$, and $\{S^{g}_{\text{meas}}\}$ (and some other intermediary lists detailed below), which we argue suffice for the resource estimations of $\{U\}$.

RRE employs the stabilizer-tableau-based Jabalizer FT-compiler[30] as an essential intermediate layer to transform individual logical subcircuits into corresponding subgraphs, consumption schedules, MBQC-formatted circuits, and other critical information. In particular, each MBQC-formatted circuit is an optional qasm output of the Jabalizer that corresponds to a completely positive trace-preserving (CPTP) map transforming the larger space of all graph qubits to the smaller space of $\ket{\text{state}}\bra{\text{state}}$ (i.e. equivalent to the execution of $U\in[U_{0},\dots,U_{n_{\text{widgets}}-1}]$). The circuit is extracted from the corresponding graph, Pauli corrections, and the consumption schedule and contains logical $H$, Phase and $CZ$ gates, as well as mid-circuit measurements. We refer readers to [30, 24, 19] for complete details on the different components of Jabalizer’s scalable compilation method.

Given an individual logical subcircuit $U_{i}$, containing only Clifford+$T$+$R_{z}(\theta)$, Jabalizer can generate the following set of quantities,

Here, $S_{\text{frames}}^{g_{i}}$ is a set tracking in which Pauli corrections (frames) must be applied to the measurements at the end. $S_{\text{output}}^{g_{i}}$ denotes the set or map of output nodes in the graph state. Because we decompose the complete algorithm in the time with a uniform width covering all qubits, the number of intermediate input and output nodes are always the same, so knowing $\{S^{g_{i}}_{\text{input}}\}$ are not required for resource estimation purposes. $n^{g_{i}}_{\text{logical}}$ is the maximum quantum memory needed to process $U_{i}$. $S_{\text{consump}}^{g_{i}}$ is the graph consumption schedule to perform individal logical qubit measurements. $S_{\text{meas}}^{g_{i}}$ is the list identifying which rotational basis must be used for each of the logical measurements. Finally, $\mathcal{E}^{\text{Jabalizer}}_{\text{MBQC}}$ is the MBQC-formatted qasm extracted from the previous outputs, which corresponds to an CPCT map from graph to output qubits executing the corresponding portion of the logical algorithm. For the widget in step $i$, this correspondence can be written as $\mathcal{E}^{\text{Jabalizer}}_{\text{MBQC,i}}(\ket{\text{state}_{i},00\dots 0}\bra{\text{state}_{i},00\dots 0})=\ket{\text{state}_{i+1},00\dots 0}\bra{\text{state}_{i+1},00\dots 0}$.

The process of eq. 7 is repeated for all $n_{\text{widgets}}$ subcircuits, providing the user with the information required for the approximate estimation set, eq. 5, excluding the preparation schedules, $\{S^{g}_{\text{prep}}\}$, requiring another tool detailed below.

RRE employs the Substrate Scheduler [31] to map the graph states $\{g\}$, calculated by Jabilizer in eq. 7, to the set of graph initialization schedules $\{S_{\text{prep}}^{g}\}$. These initialization schedules are optimized for the quantum bus in the superconducting hardware layout detailed in section I.3. We assume the bus contains linear arrangements of surface code data qubits to store the graph alongside a linear ancillary system. The bus is laid out in comb-like patterns in the cryogenic modules, facilitating the measurements. The Substrate Scheduler currently finds an optimum schedule based on an input of individual $g$ mapping to this layout. The operation of the scheduler can be improved by including additional information on node types and adding support for customized bus architecture in the future.

The robust verification of our compilation and estimation methodology at different layers was a priority when designing RRE. To this end, we have developed a verification protocol, fig. 2, to precisely check the correctness of individual subgraphs $\{g\}$ as long as the widgets originated from small circuits with 10’s of logical qubits.

Following the process shown in the fig. 2, the input logical unitary, $U$ (e.g. a small algorithm or a subcircuit of a bigger one in the qasm format), is provided to the test module. As detailed above, RRE facilitates running Jabalizer on the input to output a unified JSON containing the graph state, consumption schedule, and Pauli tracking information. The MBQC-formatted qasm for $U$ is particularly relevant to the verification protocol, which we label with $\mathcal{E}^{\text{Jabalizer}}_{\text{MBQC}}$ in the figure to highlight the correspondence discussed earlier. Now, independent from the operations of the compiler, one can build the inverse of the original algorithm for small-size cases, i.e. $U^{\dagger}$. For example, if $U$ is unitary for a logical QFT4, then $U^{\dagger}$ is InverseQFT4 manually generated. We then use a quantum simulator software to run noiseless simulations of the staggered circuit instructions shown in the lowest panel of fig. 2 with the input state of the verification protocol set, as an example, to $\ket{s_{0}=0,s_{1}=0,\cdots,s_{n_{\text{input}}-1}=0,00\cdots 0}$. By the unitarity principles, if the $\mathcal{E}^{\text{Jabalizer}}_{\text{MBQC}}$ was valid (which implies the validity of the graph state and the rest of Jabalizer outputs), then simulation output counts must exactly correspond to $\ket{00\cdots 0}$ for the example input.

We have successfully simulated and verified graphs produced in RRE as per fig. 2 for several small-size input circuits containing QFT, single TOFFOLI, and Clifford instructions (where some were widgetized from larger algorithms) for all possible input stats as

In this section, we formulate the space-time volume requirements for graph state processing parts based on $\{\text{config}\}$, $\tilde{S}_{\text{est}}$, and the architecture detailed in section I.3.

We consider an architecture with spatially optimal quantum bus and time-optimal compilation strategies. The architecture contains two fridges based on the assumption that superconducting platforms have a low quantum cycle duration and would benefit from minimizing physical footprint and interleaving graph preparation and consumption. This means we always create $n_{\text{logical}}$ surface code patches per linear register of the bi-linear bus (one line for data and another one for ancilla qubits). We arrange patches in a comb-like arrangement of the bi-linear bus for processing the subgraphs and hosting ancillae in the fridges. Any remaining space is always used to fit as many $T$-factories fridge to service the bus. There will be an equal number of $n_{\text{distilleries}}$ per each fridge – further details below.

The space-time volume of the complete graph creation on the bi-linear quantum bus before each single qubit logical measurement,

where $\mathbb{L}_{\text{prep}}(\tilde{S}_{\text{prep}})=\sum_{\text{sub}\in\tilde{S}_{\text{prep}}}L_{\text{sub}}$ is the sum of the sub-schedule lengths, $L_{\text{sub}}$ for set $\tilde{S}_{\text{prep}}$. $L_{\text{sub}}$ is the number of distinct measurement steps for the graph initialization part. For notational convenience we will extend and use this notation of $\mathbb{L}$ for other sets as well.

We now detail how we select the class of $T$-state distillation widgets placed in each modules, as described in fig. 4. We always choose the type of distilleries from widgets detailed by Litinski in Table 1 of [32]. Furthermore, we consistently assume a superconducting homogeneous physical error rate of $p=10^{-3}$ for all quantum operations or, at least, for the dominant error rates, whether from measurements or two qubit physical gates. Therefore, from this table, we have extracted the widgets matching and summarized those in table 1 forming a “Widget Lookup Table”, which we recursively go through to find a suitable T-factory type. table 1 is part of $\{\text{config}\}$ that the user provides and can customize for RRE.

table 1 details a selection of distillation widgets numerically simulated in [32]. The following are important for the “lookup table”: The output probability, $p_{out}$, sets the error rate for the ancilla state and bounds the error rate on any T-gate implemented using the ancilla. The number of qubits, 2D patch dimensions, and physical cycle sweeps depend on the widget construction. We have only considered the selection of distillation widgets that produce $T$-states (hence ignoring any entries for $CCZ$ in the original table). In our architecture, we assume there are

distilleries per fridge servicing the widgets with purified T-states for each node as measured.

$N_{\text{avail-qubits}}$ is the total number of physical qubits available for each QPU fridge, and $Q$ is the number of physical qubits required for a distillery. The above ignores the geometrical considerations required to match fridge and widget structures. Such considerations are planned as improvements in future RRE versions.

Next, we must include the resource cost for the gate-synthesis decomposition to Clifford+$T$ executed at the logical measurement step of each graph node, consuming resource states of the distillation factories, at the end of each subgraph processing step. We assume the graph nodes are either of the type requiring arbitrary $Z$-basis measurements or a single $T$-basis one. The measurements that require a $Z$-axis rotation are decomposed into $T$-gate sequences of length $L_{\varepsilon}=\lceil c_{0}\log(1/\varepsilon)+c_{1}\rceil$ as instructed by the gate synthesis method available to RRE. $\varepsilon$ is known as the diamond-norm error for decompositions to Clifford+$T$ – see below and [29] for details.

For each distilled $T$-state, $C$ cycles are needed to prepare the state, and $d$ is required to interact the state with the graph node and measure it in the $X$ or $Z$-basis. Hence for each node measurement, assuming an $R_{z}$-basis measurements, the volume increases as $2n_{\text{logical}}\mathbb{L}_{\text{prep}}d\rightarrow 2n_{\text{logical}}\mathbb{L}_{\text{prep}}d+2L_{\varepsilon}(C+d)$. Hence, for all nodes of both types $T$ and $R_{z}$ that need to be measured sequentially, the space-time volume is

where we introduced $N^{T}_{tot}=N^{R_{z}}_{\text{init}}L_{\varepsilon}+N^{T}_{\text{init}}$, which is our proxy for the total number of $T$ gates, $T$-count, and a quantity typically much larger than $N_{\mathcal{G}}$. Note that we do not have to include the space-time volume of the actual distillation widgets in the calculations as their errors are subsumed into the $p_{out}$ of the distilled state, specified in [32, Fig. 1]. In order to write the total space-time volume and simplicity of the estimation of $d$, we assume the worst-case scenario of no parallelization coming from $n_{\text{distilleries}}$ distilleries, which is equivalent to finding the largest volume and $d$ satisfying rate conditions with a single T-factory (we will revert to exploiting all possible distilleries for hardware time calculations). However, when widgets from table 1 are chosen to have a second (20-to-4) layer (the second and third widgets), we still replace $N_{\mathcal{G}}\rightarrow N_{\mathcal{G}}/4$ and $N^{T}_{tot}\rightarrow N^{T}_{tot}/4$ affecting the space-time volume, as the distillation widget produces 4-outputs simultaneously.

A single cycle time of the surface code has a failure probability governed by the following physical-to-logical rates scaling (e.g. see [6, 33]):

where the coefficients, $\kappa$ and $p_{\text{thresh}}$, are directly informed by Minimum Weight Perfect Matching (MWPM) decoder simulations and device modeling. In particular, decoder simulations provide logical and physical probabilities while we sweep through $d$ and consume one qubit and two qubit error rates and channels offered by device modeling – see section V for details. Consequently, a single unit cell containing $d$ cycles fails with probability

and failure rate of a volume $V^{\text{sequential}}_{\text{total}}$ of these cells is

Assuming that $p_{C}\ll 1$ then $(1-p_{C})^{d}\approx 1-p_{C}d$, and thus $p_{\text{log-cell}}\approx p_{C}d$ and,

Subtracting each side from 1, and taking the logarithm of each side, this becomes,

We remind the reader that the total algorithm failure rate, $p_{\text{algo-fail}}$, represents the failure probability of the entire space-time volume and, hence, the computation. This number is budgeted by the user. We now need to solve eq. 14, which can be rewritten,

where $J=\log(1-p_{\text{algo-fail}})$. The minimum code distance required is the smallest $d$ integer where the LHS of the above equation is negative. The only free variables above are the distillation widget’s cycle time, $C$, and $\varepsilon$. Note $\mathbb{L}_{\text{prep}}$ and $n_{\text{logical}}$ were previously extracted from Jabalizer, Pauli Tracker, and the Substrate Scheduler.

To determine $C$ and $\varepsilon$, we perform nested loops to iterate through widgets and then gate synthesis’ $\varepsilon$. The inner loop feeds back into the Clifford+$T$ decomposition as the failure rate of a unit cell, $p_{\text{log-cell}}$, sets the approximation error as $\varepsilon<p_{\text{log-cell}}$. By choosing a particular widget, we fix $C$, which also fixes the distilled state’s output error rate, $p_{out}$.

In a compiled quantum algorithm, we wish to ensure that the implemented $T$-gates have precisely the same error rate as any other gates in the compiled circuit. A Pauli initialization or measurement takes a single unit cell in the space-time volume. Hence, its failure rate is governed by eq. 13. This error rate must be dominant. Consequently, the error output, $p_{out}$, from the distillation circuit should be less than the failure rate of this cell (i.e., we must ensure that the output of distillation is at least as good as a unit cell in the space-time volume). Hence, we first choose the smallest value of $C$ from table 1 and solve eq. 15 for $d$. We then check that, for a physical error rate of $p=10^{-3}$ and our solution for $d$ if

is satisfied for the corresponding $p_{\text{out}}$ entry. If this is not satisfied, we must choose a better-performing state distillation widget in the outer loop. We move to the next entry on the list and repeat until we find a distillation widget that performs at least as well as the failure of a logical unit cell.

Once the particular distillation widget is chosen, our total number of required physical qubits and total intra-modular runtime, including first nonparallel graph preparation and all distillation-injection-consumption operations, can be calculated respectively as,

where $t$ is the characteristic gate time (which sets the quantum tock of all intra-modular operations for the architecture), $t^{\text{intra}}_{\text{tock}}=8td$ and $N_{\text{tocks}}=L^{0}_{\text{prep}}+N^{T}_{tot}(1+\frac{C}{d})$. Note we only need to include the temporal cost of initializing the first subgraph $g_{0}$ with $L^{0}_{\text{prep}}(S^{0}_{\text{prep}})$ since the remaining graph preparations will be performed in parallel with distillation-consumption operations as demonstrated in the following section.

To validate the assumptions of the two-fridge architecture, detailed below, we must also verify that for each $n_{\text{widgets}}$ of scheduling steps, the preparation time for the upcoming subgraph is always greater than the cumulative time to perform intra-modular distillation-consumption operations. For every widget $i$, these include distillation and transferring (injecting) $T$-states, $t^{i}_{\text{$T$-inject-intra}}$, and all surface code operations and parity measurements to consume graph nodes, $t^{i}_{\text{consump-intra}}$. Formally,

In the above, we exploit the parallelization coming from multiple distilleries. $n^{i}_{\text{reps-consump}}=\left\lceil\frac{n^{i}_{T}}{n_{\text{distilleries}}}\right\rceil+L_{\varepsilon}\left\lceil\frac{n^{i}_{R_{z}}}{n_{\text{distilleries}}}\right\rceil$ and $n^{i}_{\text{reps-inject}}=\left\lceil\frac{n^{i}_{T}+L_{\varepsilon}n^{i}_{R_{z}}}{n_{\text{distilleries}}}\right\rceil$ indicate the number of nonparallel surface code cycles at each stage. $n^{i}_{T}$ and $n^{i}_{R_{z}}$ denote the number of $T$ and $R_{z}$-type nodes in the sub-step $i$ of the consumption schedule, respectively. Notice also $t^{\text{tot}}_{\text{$T$-inject-intra}}:=\sum_{i=0,\dots,n_{\text{widgets}}-2}t^{i}_{\text{$T$-inject-intra}}=8t(C+d)N^{T}_{tot}$.

The total FT-hardware time then becomes

where,

In the above, $t^{\text{handover}}_{\text{tock}}=8t^{\text{handover}}d$ and $t^{\text{handover}}$ is the characteristic time for the inter-modular operations to teleport Bell states needed to stitch and prepare the next subgraphs on the real FTQC. $n_{\text{inter-pipes}}$ is the number or bandwidth of inter-modular pipes connecting fridges, each with a width of $d$. Overall, $t^{\text{tot}}_{\text{handover-inter}}$ is the cumulative time for the graph handovers, which scales by the number and the speed of inter-modular pipes. $t^{\text{decode-delay}}_{tot}$ is the delay added to perform decoding, fully classically, for all QEC cycles. We suggest a method to estimate $t^{\text{decode-delay}}_{tot}$ in section A.1.2 based on tensor network decoding[34, 35].

We have summarized the workflow as algorithm 1 below.

As shown in Figure 4, the system for which resources are estimated consists of two identical modules each of which handles for graph state creation and consumption as well as magic state distillation for local T-states. In each module graph creation is accmplished by laying out a comb of logical qubits and acillae fig. 7 to approximate a folded bi-linear chain for the bus to apply measurements to each logical qubit mediated by the ancilla bus. This bi-linear bus pattern takes up as much of the module as required, to host the graph state processing for the desired algorithm. It is likely that there are more efficient layouts for embedding this.

T-states are provided by magic state distillation factories that are laid out among the remaining physical qubits. According to the total algorithm and widget accuracy the suitable magic state distillation factory will be selected, from the options outlined in table 1[32]. The ancilla bus will link all T-factories and logical qubits and connect along any inter-module connections to create a continuous bus across both modules.

For executing algorithms that may require more resources than a single module, we propose to couple two modules using qubits only connected along the ancilla bus. This is shown in fig. 4. Here the widgetization of the logical algorithm allows for some degree of parallelization between modules. As shown in fig. 3 each module will prepare the graph state for a given widget while the other module consumes the graph for the previous widget. Then the output nodes of the graph will be teleported to the nodes of the newly prepared widget before consumption of that graph. This cycle repeats until the all of the widget sub-units are processed and the entire algorithm is executed.

We can estimate the time required to execute the schedule of FT operations across this system of resources accounting for local and inter-module interactions. To do this, we abstract the teleportation interaction between modules as a single time described as an intermodule handover time, which we are assuming to be 1 $\mu$s, and we apply this timing for any schedule operations that cross a cryostat boundary. For a fault tolerant teleporation, as mentioned in section I.1.6, the total handover time is dependent on the required code distance for the algorithm.

In addition to algorithmic execution times and logical qubit resources, RRE also can compute the hardware resources required to execute a given algorithm, including signal line counts, functional chip area and number of amplifiers. For superconducting qubits that operate at milliKelvin temperatures, heatloads are a primary concern and power consumption to operate the computer are estimated.

Cryogenic energy consumption calculations for this architecture are based on per-line thermal loads. These are calculated assume a signaling solution with all-superconducting input and output lines below 4K, with conventional dispersive readout using a Traveling Wave Parametric Amplifier (TWPA) and High-Electron-Mobility Transistors (HEMT). Each qubit is assumed to require a microwave line for XY control and a flux bias line for tuning qubit frequency. Each pair of qubits has a tunable coupler [55] that requires a single flux line. Dispersive readout is assumed to read out 10 qubits multiplexed in parallel. As heat loads are dominated by HEMT dissipation and static heat loads, only those are considered currently. Example numbers are provided in the open source repository, based on estimates from [56]. Power consumed by the overall system is scaled by rough estimates for the power required to operate 4K cryo-plants [57] (500 kW of power for 500 W of cooling at 4 K) and dilution refrigerators (1 kW of power for 1 $\mu$W of cooling at 20 mK) as in the published work[57] and manufacturer specifications. We expect these estimates to change over time as they are dependent on detailed designs of the signal chain.

The Los Alamos National Laboratory (LANL) has created a public repository containing circuit generation and analysis of problems that may be amenable to quantum computing, and ostensibly are of interest to the scientific community, [59]. Drawing inspiration from their investigation of magnetic lattices, we perform resource estimation for simulation of a time-invariant Transverse-Ising Hamiltonian,

on a triangular lattice, shown in Figure 8.

Here, $X_{\textbf{i}}$ and and $Z_{\textbf{i}}$ represent the Pauli X and Z operators, respectively, acting on lattice site i, and angle brackets denote lattice sites connected by an edge.

Based on the studies documented by LANL, we carried out resource estimates of simulation of the Hamiltonian in eq. 22 for total time of $t=1$ on a 2D triangular lattice of size $N\times N$ where $N$ is an integer between $2$ and $10$.

We also investigated Fermi-Hubbard simulation as a core computational capability towards solving problems related to room temperature superconductivity. These problems typically involve studying the ground state energy, state configurations, and other observable properties of the system, including time-dependent quantities.

The Fermi-Hubbard Hamiltonian we consider is given by,

where single and double angled brackets denote all nearest neighbor and next-nearest-neighbor lattice site pairs respectively. $c^{\dagger}_{\textbf{j},\sigma}(c_{\textbf{j},\sigma}$) creation (annihilation) operators of a spin-$\sigma$ Fermion on the jth lattice site where $\sigma\in\{\uparrow,\downarrow\}$. $J$ (primed or unprimed) and $U$ capture the magnitudes of the hopping and repulsion energies respectively. $n_{\textbf{k},\sigma}=c^{\dagger}_{\textbf{k},\sigma}c_{\textbf{k},\sigma}$ is the number operator, $h_{z}$ is an external applied magnetic field, and $\mu$ is the chemical potential.

To motivate the specific Fermi-Hubbard simulations we will obtain resource estimates for, we consult the literature to understand configurations that are of interest to the community. Broadly speaking, methods for simulating a Fermi-Hubbard system can be broken into two families; stochastic and deterministic. As mentioned in [60], deterministic methods do not “suffer from numerical issues seen in stochastic simulations when dealing with systems with non-zero chemical potential”. More over this instability due to the “sign problem” becomes prohibitive for moderate $\mu\sim 0.3$. This suggests separating instances into those with $\mu=0$, and those with greater than $\mu>0.3$.

Here we reference a few exemplar studies in the literature for both stochastic and deterministic methods.

Some examples of stochastic studies include a study using Auxiliary-field quantum Monte Carlo (AFQMC) to perform ground state estimation on lattice sizes of $4\times 4$ to $16\times 16$ and $U/J=4,8$, $J^{\prime}=0$. [61]. In addition, [62] used Determinant Quantum Monte Carlo (DQMC) to study doped ground state estimation using $U/J=6$, $J^{\prime}/J=-.25$, and lattice sizes of $8\times 8$, $12\times 8$, $12\times 12$, $16\times 8$. In contrast, examples of deterministic studies included work on numerical linked cluster expansions [63]. In this particular study, researchers selected random repulsion strengths from an uncertainty "box" about $U_{0}$. Simulations were for a periodic $10\times 10$ lattice with $U_{0}/J=8$, $J^{\prime}=0$. Others researchers used Projected Density Matrix Embedding (PDME) to perform a 2D Fermi-Hubbard study with $U/J=2,4,6,8$, $J^{\prime}=0$ on a $40\times 40$ lattice [64]. Finally, in [60] the authors use Projected Entangled Pair States (PEPS) with a bond dimension of $D=10$, and claim that accuracy increases with $D$. Most of their results are from a $3\times 4$ honeycomb lattice with $U/J\in[0,6]$ and $J^{\prime}=0$, but the authors were able to calculate the ground state via imaginary time evolution on a $15\times 15$ bipartite honeycomb lattice.

Of note is that in [60] it was determined that the time required to run these studies on a CPU varies like,

where $\chi$ is the truncated bond dimension and $d$ the problem dimension. Memory scaled as,

The authors of the study empirically found $\chi=3D$ was sufficiently accurate for their studies, with error about $1\%$. These estimates provide an interesting guidepost for the tipping point at which FTQCs quantum computers may outperform classical computers for these types of problems.

Taking these studies as examples of the types of systems of interest to the community it suggests obtaining resource estimates for Fermi-Hubbard Hamiltonian simulations consisting of a 2D square lattice for total time $t=1$, with $J=1$, $J^{\prime}=0$, $h_{z}=0$, $\mu=0.5$, $U=8$, and lattice sides of length $N=2,3,4,5,6,7,8$.

Owing to its optimal scaling characteristics, in terms of queries to the block encoding of the Hamiltonian [10], we have chosen Quantum Signal Processing (QSP) as the Hamiltonian simulation algorithm for which we obtain resource estimate for. QSP is briefly described below, but more details can be found in [10], [65], [66], and [67].

Given a Hamiltonian $H$ and simulation time $t$, QSP involves performing a sequence of single qubit rotations, parameterized by a sequence of phase angles, interlaced with a quantum walk-like operator. This allows one to apply a polynomial approximation of $e^{-iHt}$ to the input target state. The sequence of phase angles are responsible for the polynomial approximation, and are computed classically offline, according to the error tolerance of the approximation error desired, captured by the user defined parameter $\epsilon_{\text{qsp}}$. $\epsilon_{\text{qsp}}$ is also inversely related to the probability of success of the QSP algorithm. Thus, for an accurate simulation, with a high probability of success, $\epsilon_{\text{qsp}}$ may have to be set quite small.

For implementation details, we follow the examples accompanying the pyLIQTR [68] package (v0.3.1), as well as the work done in [69, 70].

For demonstration purposes, we set the desired probability of success of our algorithm, $p_{\rm algo}$, to be equal to 0.95. In real world applications this is typically fixed to meet outcome requirements (such as accuracy, wall clock time, etc). Denoting the probability of success of a QSP circuit, $p_{\rm qsp}$, we therfore require,

Since circuit size is inversely related to $\epsilon_{\text{qsp}}$, it is prudent to choose $\epsilon_{\text{qsp}}$ to be as large as possible while still guaranteeing that $p_{\rm qsp}$ exceeds the algorithm target of 0.95. According to [10] $(1-2\epsilon_{\text{qsp}})\leq p_{\rm qsp}$. Thus, to ensure eq. 24 is satisfied we set,

and so in this work we choose $\epsilon_{\text{qsp}}$ to be equal to 0.025.