A fault-tolerant continuous-variable measurement-based quantum computation architecture

Continuous variable measurement-based quantum computation on cluster states has in recent years shown great potential for scalable, universal, and fault-tolerant quantum computation when combined with the Gottesman-Kitaev-Preskill (GKP) code and quantum error correction. However, no complete fault-tolerant architecture exists that includes everything from cluster state generation with finite squeezing to gate implementations with realistic noise and error correction. In this work, we propose a simple architecture for the preparation of a cluster state in three dimensions in which gates by gate teleportation can be efficiently implemented. To accommodate scalability, we propose architectures that allow for both spatial and temporal multiplexing, with the temporal encoded version requiring as little as two squeezed light sources. Due to its three-dimensional structure, the architecture supports topological qubit error correction, while GKP error correction is efficiently realized within the architecture by teleportation. To validate fault-tolerance, the architecture is simulated using surface-GKP codes, including noise from GKP-states as well as gate noise caused by finite squeezing in the cluster state. We find a fault-tolerant squeezing threshold of 12.7 dB with room for further improvement.


I. INTRODUCTION
In measurement-based quantum computation (MBQC), gates are implemented by projective measurements on a multi-mode entangled cluster state, circumventing the complex coherent unitary dynamics required in conventional gate-based quantum computation [1]. As such, the cluster state is a critical resource for MBQC, and its number of modes and structural design defines the size of a potential measurement-induced algorithm. A particularly promising platform for scaling and controlling the structure of a cluster state is the optical continuous variable (CV) platform [2,3], where large cluster states can be deterministically generated and controlled, and efficiently measured by homodyne detection. This has been proven by the realizations of large-scale CV cluster states in both one dimension [4][5][6] and two dimensions [7,8]. Moreover, the versatility of the CV optical platform has been further corroborated by the recent demonstrations of single-and multi-mode gates using high-efficiency projective measurements on one-dimensional [9] and two-dimensional cluster states [10].
MBQC based on CV is however inherently noisy due to the impossibility of generating maximally entangled CV cluster states: The generation of maximal CV entanglement requires squeezed states of infinite squeezing and thereby infinite energy, which is not feasible. Therefore, inevitably, Gaussian noise will be added to the quantum information during computation. To combat this additive noise, information is encoded as special qubits in CV bosonic modes of infinite dimension. By encoding such qubits into the bosonic modes, using e.g. a cat-code [11], a binomial code [12], or the Gottesman-Kitaev-Preskill (GKP) code [13], the Gaussian noise can be corrected at the cost of being converted into Pauli errors on the encoded qubit. These Pauli qubit errors must then be corrected by some qubit quantum error correction scheme. Implementing qubit error correction efficiently in MBQC puts stringent requirements on the underlying cluster state. As an example, the local connectivity in the cluster states support only coupling between nearest-neighbor modes, so topological error correction is a natural choice for qubit error correction [14]. This, in turn, requires a three-dimensional (3D) cluster state for MBQC [15][16][17].
Different proposals on 3D cluster state generation and topological MBQC exist. Fukui et al. [18] suggested a scheme for fault-tolerant MBQC based on topological error correction, but their scheme assumes the availability of a highly complex 3D cluster state of encoded qubits. Wu et al. [19] proposed an optical setup for the generation of a 3D cluster state using time and frequency multiplexing. However, in their proposal, gates are implemented by gate teleportation through four-mode square cluster states leading to increased gate noise. In another work of Fukui et al. [20], an all-temporally encoded 3D cluster state is proposed, but this scheme is experimentally highly challenging as it requires the construction of 12 squeezing sources and real-time feed-forward opera-tions. Moreover, no schemes for qubit encoding and qubit error correction was put forward. The most complete work on CV MBQC to date is carried out by Bourassa et al. in Ref. [21] in which a computation architecture for the generation of a 3D cluster state combined with topological MBQC is proposed. However, the suggested architecture is based on spatial encoding, rendering the number of spatial resources very large (as this number scales linearly with the computation size). Moreover, their scheme relies on a very large number of experimentally challenging on-line swap and sum gates which they assume to be ideal. Experimental work towards topological quantum computation has been demonstrated in other platforms, including a 9-qubit code in a photonic platform with polarization encoded qubits [22], and a 7-qubit code in an ion-trap platform [23] and a superconducting platform [24]. Still, thousands of qubits are required for large fault-tolerant codes [25].
In our work, we propose a simple, scalable, and complete architecture for topological MBQC and validate the fault-tolerance of the computation scheme. It is based on gate teleportation on parallel one-dimensional (1D) cluster states, or wires, arranged in a 3D lattice and coupled by variable beam-splitters for two-mode gates. As such, the setup is a variation of the well-demonstrated 1D cluster state generation [4,7,9] with added variable beam-splitters. Combined with GKP-encoded qubits [13], the scheme allows for universal computation, while fault-tolerance is achievable by encoding logical qubits in the topological surface code [25][26][27]. Furthermore, the scheme, being based on gate teleportation, is compatible with a recently proposed GKP correction protocol that dispenses with demanding coupling to ancillary GKP-qubits [28]. We validate the fault-tolerance of the full scheme by a thorough simulation that includes both noise in the GKP-qubits and-unlike previous worksgate noise caused by finite squeezing in the cluster state. As a result, when combining the topological surface code with GKP error correction in the surface-GKP code [29], we find a squeezing threshold of 17.3 dB. We continue to propose a variation of the surface-GKP code-the surface-4-GKP code with four GKP corrections during the surface code syndrome measurements-by which we upgrade the squeezing threshold to 12.7 dB while leaving room for further improvements. Related schemes in Refs. [18,21] have better thresholds, but they assume an ideal cluster state, i.e. without gate noise. We obtain a comparable threshold of 10.2 dB when ignoring gate noise.
The paper is organized as follows. In section II we present the computation scheme and describe the implementation of the required gates. In section III we focus on GKP error correction within the computation scheme, and in section IV we implement the surface code for qubit error correction and validate the fault-tolerance properties by performing simulations. In section VI we discuss the results and conclude the paper. n×m wires of two-mode entangled states are prepared in n×m resource preparation gadgets (res. prep.). States with input information for computation, |ψin GKP , are switched into each wire for computation using an optical switch. With a beamsplitter (marked with a red arrow) and two detectors for each wire, single-mode gates are implemented in each wire when feeding the resource preparation gadget with squeezed states, |0 sq [9,30]. Feeding the resource preparation gadget with special GKP states, namely qunaught states, |∅ GKP , allows for GKP quadrature correction of encoded GKP qubits [28]. Neighbouring wire setups are connected using variable beamsplitters (VBS) allowing tunable coupling for implementing two-mode gates and thereby enabling multi-mode computation.

II. COMPUTATION SCHEME
The concept of our computation scheme is illustrated in Fig. 1. The scheme consists of parallel wires, each corresponding to temporally encoded one-dimensional cluster states [4,6] on which single-mode gates can be implemented by projective measurements using a beamsplitter (marked in Fig. 1 with a red arrow) and two detectors [9,30]. Input states, |ψ in , can be swapped into each wire for computation using an optical switch. To enable multi-mode computation, the setups of neighbouring wires are connected with variable beam-splitters (VBS), which allow for a tunable coupling of wires for implementing two-mode gates. This architecture may be implemented spatially as depicted in Fig. 1 with a large grid of wire setups forming a 3D cluster state encoded in (space) 2 × time. Such spatial encoding requires spatially scalable resources and may be possible with integrated photonics. As an alternative, in Fig. 2 we propose an all-temporal encoded version of the computation scheme which allows for a simple experimental implementation and easy scalability. In the following, while focusing on the temporally encoded architecture when describing the computation scheme in detail, all the presented methods, results, and conclusions also hold true for the spatially encoded architecture-even a combination of the spatial and temporal encoding architectures might lead to a similar computation scheme with identical conclusions.
The temporally encoded scheme in Fig. 2 consists of three parts: the preparation of resource states, the in- (a) Temporally encoded computational setup of our scheme. The scheme consists of three parts: The resource preparation gadget (res. prep.); the computational level where the computation takes place; and the temporally delocalized measurement device (TDMD) for gate implementation by projective measurements. This scheme utilizes temporal multiplexing of two spatial modes, A and B, marked in the computational level. (b) Wires of two-mode entanglement at the computational level shown in the time domain, here for the simple case of nm = 9. Bold lines represent two-mode entanglement, while thin lines indicate the temporal overlap of A and B. The wires begin with nm input states in temporal modes 0 to 8, switched in using an optical switch in A of the computational level. The colors of the wires have no physical meaning and are merely used to indicate different wires. (c) Wires rearranged into a 3D time lattice where the input is encoded onto the n × m end surface while gates are implemented by teleportation along the wires in the third dimension. Here, the red arrows represent the first beam-splitter of the TDMD, while the dotted blue and green arrows represent the variable beam-splitters of the TDMD. The first 10 temporal modes of (b) from 0 to 9 are labeled in (c).
jection of input states at the computational level, and the measurements, enabled by a temporally delocalized measurement device (TDMD). Note, the term 'computational level' refers to the location in the setup at which information is encoded and computation takes place. In some of our previous works this computational level is referred to as the 'logic level' [10,31]. However, in this work we reserve the term 'logic' for qubit error correction in section IV. As ancillary input for the resource preparation, we switch between squeezed vacuum states, |0 sq , when implementing gates by projective measurements, and GKP qunaught states [28], |∅ GKP , when performing GKP error correction. In section III, GKP error correction with ancillary |∅ GKP -states is described, while throughout this section we focus on gate implementation with ancillary |0 sq -states.
At the resource preparation stage, the spatial modes A and B are initially occupied by squeezed vacuum states, |0 sq , which are squeezed along the orthogonal quadratures (q−p)/ √ 2 and (q+p)/ √ 2, respectively. Here,q and p are the electric field amplitude and phase (or position and momentum) quadratures for which we use the = 1 convention, corresponding to a vacuum variance of 1/2. Each pair of squeezed states is then interfered on a balanced beam-splitter, leading to two-mode entanglement withqp-correlations. This is an approximate cluster state equivalent to a conventional two-mode squeezed state (withqq-andpp-correlations) that is phase-rotated by π/4 in both modes [32,33]. As a unitary operator for the balanced beam-splitter, we useB = e −iπ(qi⊗pj −pi⊗qj )/4 with corresponding symplectic matrix acting on (q i ,q j ,p i ,p j ) T quadrature vectors, and represented graphically with an arrow pointing from mode i to j. Note, in this work we prepare two-mode cluster states, however, we could as well have considered preparation of conventional two-mode squeezed states (witĥ qq-andpp-correlations) which are equivalent to cluster states under phase-rotation that may be absorbed into the measurement bases.
After the interference at the beam splitter, the modes of A are delayed by nm temporal modes, leading to synchronization of the modes A and B of two-mode entangled states that initially are separated by nm temporal modes [34]. The result is nm decoupled wires of twomode entangled states, illustrated in the time domain in Fig. 2(b) for nm = 9. Here, each color indicates different wires with the bold lines indicating two-mode entanglement, while the thin lines indicate temporal overlap of A and B. The nm wires constitute the computational level in which computation is performed. Using an optical switch, nm input modes to the computation can be switched into the computational level at A-optical switching has been demonstrated in a continuous variable quantum setting in [35,36]. In Fig. 2(b) nine input modes have been switched into A in temporal modes 0 to 8.
The nm wires can be arranged in a 2D grid such that they form a 3D square lattice, as shown in Fig. 2(c). The third dimension is in principle arbitrarily deep. As such, information is encoded on a surface while computation proceeds along the third dimension by teleportation using the TDMD. The TDMD consists of a balanced beamsplitter, two VBSs, two delays of 1 and n − 1 temporal modes, and two homodyne detectors (HD) measuring A and B in basesq(θ) =q cos θ +p sin θ. The arrangement is illustrated in Fig. 2(a). Each VBS can vary between two settings: when implementing single-mode gates, the VBSs are left 'open' such that the modes A and B do not interfere, corresponding toÎ A ⊗Î B ; when implementing two-mode gates, one of the two VBSs are 'enabled' to be functioning as a balanced beam-splitter with the symplectic matrix in Eq. (1) interfering A and B. Such a variable beam-splitter may be implemented in various ways, for instance as a Mach-Zehnder interferometer with a controllable phase in one arm, or by polarization control combined with polarization-dependent beam-splitters [36][37][38].
When the VBSs are left open, the TDMD simply implements a two-mode joint Bell measurement which enacts a single-mode gate teleportation through the twomode entangled resource state [9,30]. The state in computation is teleported from temporal mode k in A, (A, k), to mode (A, k + nm). In this process, depending on the HD basis settings, θ A,k and θ B,k , the gate operation is implemented on the teleported state, where θ ± = (±θ A,k + θ B,k )/2 [39]. HereR(θ) = e −iθ(q 2 +p 2 )/2 and S(s) = e i ln(s)(qp+pq)/2 are the rotation and squeezing operators. Note, k labels the temporal modes at the computational level, while at the HDs, modes in A are delayed by n temporal modes relative to modes in B. All single-mode Gaussian gates can be implemented with two iterations of Eq. (2) [40]. Enabling one of the two VBSs, two-mode gates can be implemented between nearest neighbours in the 3D time lattice. Two-mode gates between (A, k) and (A, k + 1) are implemented by enabling the first VBS, while enabling the second VBS allows two-mode gates between (A, k) and (A, k + n). In the 3D time lattice of Fig. 2(c), the VBSs are represented by dotted arrows. To encode the surface code described in section IV, we implement two different symmetric two-mode gates:Ĉ Z (g) = e igq⊗q andĈ X (g) = e −igp⊗p . They are controlled-phase gates that displace one mode inp (orq) by an amount gq (or gp) controlled by the other mode. We note thatĈ X (g) does not correspond to a controlled-not gate.Ĉ Z (g), orĈ X (g), constitutes together with Eq. (2) a universal Gaussian gate set. In practice,Ĉ Z (g) andĈ X (g) cannot be implemented in a single computation step without some Fourier by-products of π/2 phase-rotations, F =R(π/2). To implement the surface code in section IV with a minimum number of computation steps, we make use of 4 variations ofĈ Z (g) andĈ X (g) with different by-products, each listed in table I with their re-Two-mode gate Basis setting, (θ A,k , θ B,k , θ A,k+j , θ B,k+j ) Two-mode gates with input and output in modes (A, k) ⊗ (A, k + j) and (A, k + nm) ⊗ (A, k + nm + j), respectively, and their required basis settings. Here, j = 1 when enabling the first VBS (marked by blue in Fig. 2), and j = n when enabling the second VBS (marked by green in Fig. 2). The order of the tensor products are arranged with earlier temporal modes first. Note, as apparent from the basis settings, (F † ⊗F )ĈZ (g) =ĈX (g)(F † ⊗F ) and CZ (g)(F ⊗F † ) =ĈX (g)(F † ⊗F ). However, when implementing the gates in the surface code in section IV, it is useful to consider them individually, as we will have the Fourier by-products to cancel.
quired basis settings for implementation. These are with by-products ofF ⊗F † orF † ⊗F when implemented on modes (A, k) ⊗ (A, k + j) where j = 1 or n depending on which VBS is enabled. When implementing the surface code, the gates are arranged such that the Fourier by-products cancel. Finally, as the resource squeezed states |0 sq are finitely squeezed, all gate implementations will inevitably produce excess noise which accumulates on the computational modes throughout the computation. Due to the Gaussian nature of the quadrature distribution of |0 sq , this gate noise leads to a Gaussian convolution of the quadratures of all computational modes [10,30,31]. Assuming the variance of the squeezed quadrature of |0 sq to be σ 2 = e −2r /2 (where r is the squeezing parameter), the variance of the uncorrelated gate noise is which will be added symmetrically in each quadrature of the computational modes. In addition to gate noise, an implemented gate also results in a displacement of the computational modes depending on the projective measurement outcomes. Since the measurement outcomes are known, this displacement can be compensated for by another cancelling displacement operation. However, in practice, these ubiquitous displacement operations need not be executed directly onto the output modes; they can simply be accounted for in post-processing of the measurement outcomes. Therefore, in this work, we will ignore these displacements, while for practical implementation, one has to keep these in mind when analysing the measurement outcomes. For a derivation of implemented gates considered in this section, together with their resulting gate noise and displacements, see appendix A.

III. GKP QUADRATURE CORRECTION
As mentioned above, noise will be added to the computation modes at each gate implementation due to the finite amount of squeezing of the resource states. To correct for this noise and thus prevent noise accumulation, we consider a quadrature noise correction scheme that relies on bosonic qubit encoding in the infinite-dimensional Hilbert space. This noise correction scheme, however, comes with the cost of introducing qubit errors which must be subsequently corrected by a qubit error correction scheme. The first correction layer, the quadrature correction scheme, will be discussed in this section while the second correction layer, the qubit error correction, will be the subject of section IV.
Several schemes for encoding qubits into bosonic harmonic oscillators of infinite Hilbert space dimension exist, including cat-codes [11], binomial codes [12], and the Gottesman-Kitaev-Preskill (GKP) code [13]. Since the gate noise of our computation scheme is additive quadrature noise, GKP-encoding where a qubit is encoded in the mode quadratures as Dirac combs is most suitable. The GKP code is also suitable for correcting excitation loss errors since excitation loss can be converted via quantumlimited amplification into additive quadrature noise [41][42][43][44]. Furthermore, as the gate noise (with variance given in Eq. (3)) is added symmetrically in phase space, we consider GKP qubits encoded on square grids in phase space with a 2 √ π × 2 √ π unit cell. For such encoded qubits, a universal Clifford gate set is realized by the Gaussian gates {R(π/2),P (1),Ĉ Z (1)} together with √ π displacements in phase space. For a comprehensive review of the GKP code, see [13,45,46].
Information encoded in GKP qubit states, |ψ in GKP , is launched into the computation scheme at the computational level as shown in Fig. 2(a). These states are not ideal as they are subjected to the finite energy constraints (similar to the squeezed states). This means that the uncertainties of the individual spikes of the quadrature comb of the GKP state are not zero but have a finite value. Mathematically, the delta functions of the Dirac comb in the GKP state quadrature wave function are replaced by finitely squeezed Gaussian functions, each with a variance of σ 2 GKP such that theq-quadrature wave functions of the approximate GKP Pauli-Z eigenstates, |j GKP (where j = 0, 1), are Here E(q) is an overall envelope that can be chosen to satisfy the Fourier relations between the orthogonalq and p quadratures with equally squeezed spikes-different finitely squeezed approximations of GKP states exist with different envelopes [47]. In the following, for the sake of simplifying the simulation of the fault-tolerance squeezing threshold in section IV, we ignore the overall envelope, i.e. we set E(q) = 1. Doing so corresponds to a noisier GKP state that is an incoherent mixture of ideal GKP states. This is therefore a conservative assumption that does not lead to false positive results in the noise model [29]. Moreover, for the small σ GKP values considered in this work, the envelopes are correspondingly broad, and we expect that ignoring these will have little effect on the simulated error thresholds presented here. We further assume that the squeezing of the GKP spikes in bothq andp quadratures is the same as that of the |0 sq states in the resource preparation, The Gaussian noise accompanying gate implementation results in the variance of the GKP spikes increasing by σ 2 gate in both quadratures for every single gate. To prevent this, GKP quadrature correction is performed, preferably after every gate. Traditionally, this is done by coupling each quadrature to ancillary GKP states, which are then measured, and the result is fed forward to displacements of the computational qubit (or compensated for in following measurement outcomes) [13]. However, this on-demand coupling of desired modes with encoded GKP qubits to ancillary GKP states requires either active squeezing, which is experimentally hard to realize, or projective measurements, which add noise. Instead, we use the new approach by Walshe et al. [28] where GKP quadrature correction is realized by qubit teleportation using ancillary GKP qunaught states and is directly compatible with our computation scheme. The GKP qunaught state, |∅ GKP , is the 1-level version of the generalized GKP qud it state with a √ 2π spacing between the spikes in the quadrature wave functions [28,48]: where the overall envelope is ignored as well, and we assume the spike-squeezing in both quadratures to equal that of the GKP qubit states with variance σ 2 GKP . As such, |∅ GKP holds no information, but interfering two |∅ GKP states on a beam-splitter results in a two-mode GKP-qubit Bell state-for more information, see [28]. This state can then be used for GKP-qubit teleportation with support only on the GKP grid in phase-space so that a noisy GKP qubit is projected into a purified GKP qubit by the teleportation.
The implementation of the GKP quadrature correction in [28] is shown in Fig. 3. In the resource preparation gadget, we switch from |0 sq to |∅ GKP states. After interference on the first beam-splitter, a GKP Bell state is prepared at the computational level instead of a twomode CV cluster state. For teleportation of a noisy GKP qubit through the GKP Bell state, a Bell-measurement of the noisy GKP qubit and one mode of the Bell state should be carried out by the TDMD. This is done by leaving the two VBSs open and measuring in theq andp basis in spatial modes A and B, respectively. The corresponding graph in a small section of the 3D time lattice is 3. Implementation of GKP quadrature correction by qubit teleportation. Qunaught states, |∅ GKP , are injected in the resource preparation gadget, thereby preparing a twomode GKP-qubit Bell state shown as two connected rectangles before the nm-delay. In the computational level after the nm-delay, one part of the Bell state overlaps in time with the GKP-qubit state to be corrected, |ψ GKP , shown as a circle. Setting the TDMD to perform a Bell measurement (the two VBSs are open and left out in the figure) |ψ GKP is teleported through the Bell state and projected into a purified GKP qubit state by the Kraus operator in Eq. (6). Below is the corresponding graph as it will appear in the 3D time lattice of Fig. 2(c) as well as the corresponding circuit diagram.
shown in Fig. 3 together with the corresponding circuit. The resulting Kraus operator, projects the noisy input state into a purified GKP qubit state. Here, √ 2q are displacements in theq andp quadratures, respectively, depending on the measurement outcomes m A and m B , N is a normalization factor, also depending on the measurement outcomes, and is a noisy GKP projector (here |j GKP are the approximate GKP Pauli-Z eigenstates with the quadrature wave function in Eq. (4), assuming E(q) = 1 and equal squeezing of spikes inq andp quadratures).
The output values, m A and m B , are integer multiples of π/2 plus some noise associated with the finite squeezing of the GKP qubit and qunaught states. As such, theX(−m A √ 2) andẐ(−m B √ 2) displacements in Eq. (6) mainly corresponds to Pauli-X and Pauli-Z operations on the encoded qubit, and is a natural result of the teleportation similar to regular qubit teleportation. These displacements may be compensated for by unitarily displacing the teleported state back inq(p) by m A(B) √ 2 rounded to the nearest integer of √ π, or simply by shifting the final measurement outcomes. However, due to the inevitable noise in m A and m B , occasionally m A(B) √ 2 will be rounded to the wrong integer of √ π which then results in a faulty displacement operation. This induces a qubit error. The probability for this error to occur is [29] p σ (z) = n∈Z exp −(z − (2n + 1) where the residual analogue information, z = R(m A(B) √ 2), when rounding is given by In Eq. (7), being the spike variance of the GKP qubit before teleportation. For example, if the GKP qubit to be corrected has gone through one gate, then σ 2 in = σ 2 GKP +σ 2 gate where σ 2 GKP was the GKP qubit spike variance before the gate and σ 2 gate is the gate noise variance in Eq. (3). In [49] it was proposed to use the analogue information from the GKP quadrature correction to improve the concatenated qubit error correction. Here, similar to [29], we use this analogue information through the probability in Eq. (7) to improve the second layer of error correction, the surface code, which is the subject of the next section.

IV. SURFACE CODE
In section III, we showed how to project the continuous variable noise from the finite squeezing in section II into qubit Pauli errors by GKP quadrature correction. However, in order to perform fault-tolerant quantum computation, such Pauli errors must then be corrected using an additional quantum error correcting code operating at the qubit level. Given the nearest neighbour interactions of the computation scheme in section II, topological qubit error correction is a natural choice to correct the Pauli errors. With information encoded on a surface of the computation scheme's 3D time lattice, and gates implemented in the third dimension, we consider the surface code [25][26][27]. Specifically, to compute logical X or Z error rates, we implement the simulation methods of [29] applied to the rotated surface code [50,51]. Such simulation methods are adapted to the computation scheme as described in appendix B and C. We note that the rotated surface code may not be the most resource-efficient code for our computation scheme since it is rotated 45 • with respect to the 3D time lattice, and thereby, computation modes located in the corner of the 3D time lattice may not be utilized-the rotated surface code was chosen in order to easily adapt the simulation method of [29]. Below, in section IV A, we first describe the implementation of the surface code, and then consider it combined with GKP quadrature correction. In section IV B, we then present simulation results of logical error rates and provide a squeezing threshold. In (e), anX Pauli error occurs on a measure-Z ancilla after the first two-mode gate of a surface code syndrome measurement round. The error then propagates to neighbouring data and measure-X ancillas via the subsequent two-mode gates used to measure the surface code stabilizers. In (f), aẐ Pauli error occurs on a data qubit after the first two-mode gate and propagates to two measure-X ancillas. See appendix C for an examination of all possible Pauli errors, how they propagate, and the corresponding edge in the matching graphs for decoding.

A. Implementation of the rotated surface code
A logical qubit is shown in Fig. 4(a) for a distance d = 5 rotated surface code. Information is encoded in d 2 data qubits (white and gray circles). The stabilizers of the code are measured using (d 2 − 1)/2 ancilla qubits prepared in |+ GKP (green circles) and (d 2 − 1)/2 ancilla qubits prepared in |0 GKP (red circles). In what follows, we refer to green and red ancillas as measure-Z and measure-X ancillas.
One round of Z and X-type stabilizer measurements is shown in Fig. 4(b,c). Each stabilizer measurement consists of four two-qubit gates and is thus implemented in four time steps along the third dimension of the 3D time lattice in which the surface code is implemented as shown in Fig. 4(d). Using the optical input-switch in spatial mode A of the setup described in section II, the ancilla qubits, initialized beforehand in the |0 GKP and |+ GKP ∝ |0 GKP + |1 GKP states, are switched into the computational level in the temporal modes corresponding to ancillary modes of the surface code. The measure-Z and -X ancillas are then coupled to neighbouring data qubits usingĈ Z (1) = e iq⊗q andĈ X (±1) = e ∓ip⊗p gates, before being measured in thep andq basis, respectively. To measure such ancillas using the TDMD, the VBSs are left open while the same basis is chosen in spatial modes A and B in which case the measurements commute with the beam-splitter of the TDMD. Note that the state initialization and measurement basis for the ancillas are opposite of what is traditionally used in the surface code since they are coupled to data qubits viaĈ Z andĈ X gates instead of sum-gates,Ĉ NOT = e −iq⊗p . The reason for not using sum-gates is that such gates cannot be implemented in the MBQC scheme considered in this work in a single set of projective measurements. As such, using sum-gates would lead to larger gate error rates compared to the error rates of theĈ Z andĈ X gates.
While the measure-Z ancillas are coupled to data qubits with a constant coupling rate throughĈ Z (1), the measure-X ancillas are coupled to data qubits witĥ C X (1) in step 1 and 4 andĈ X (−1) in step 2 and 3. This is to prevent the propagation of finite squeezing noise among measure-qubits [29] (though this does not matter in the case of GKP quadrature correction during the stabilizer measurements as discussed later). Furthermore, since theĈ Z (1) andĈ X (±1) gates cannot be im-plemented in a single computation step without Fourier by-products as described in section II, the surface code is implemented with the two-mode gates listed in table I, and so, for the different two-mode gates in Fig. 4 we use where the first term in the tensor products is the earlier temporal mode in the computational level. In this way, the Fourier by-products of step 1(3) and 2(4) cancel aŝ FF † =F †F =Î on measure-Z and odd data qubits, and becomesFF =F †F † = −Î on measure-X and even data qubits. Hence such terms have no influence on the encoded information and do not propagate errors. For CV noise, −Î on even data qubits cancels with −Î on measure-X qubits when phase-space displacements propagate in-between measure qubits. We proceed to combine the surface code with GKP quadrature correction, the so-called surface-GKP code. Commonly, in the surface-GKP code, each round of syndrome measurements consists of correction of the GKP data qubits followed by measurements of the surface code stabilizers. In this way, qubit errors induced in the GKP quadrature correction is corrected by the surface code [18,29,52]. However, in the usual surface-GKP code, gate noise accumulates during all four gates of the stabilizer measurements in Fig. 4(b,c). We propose to modify the scheme to perform GKP quadrature correction of each mode after every implemented gate. In other words, for each Z-and X-stabilizer measurement, GKP quadrature correction is performed four times, and we refer to this as the surface-4-GKP code. Unfortunately, when doing so, qubit errors are induced during the surface code stabilizer measurements with a large impact on the fault-tolerant error threshold [25]. Two examples of induced qubits errors, and how they propagate during the stabilizer measurements, are shown in Fig. 4(e,f). Here, a qubitX error on a measure-Z qubit, induced in the GKP quadrature correction after the first gate of the stabilizer measurement, propagates to three data qubits asẐ errors through theĈ Z gates (while an initialẐ error will not propagate throughĈ Z ). From there it further propagates to two measure-X qubits throughĈ X gates. Similarly, aẐ error on a data qubit after the first twomode gate propagates asX errors to measure-X qubits throughĈ X . These errors may lead to faulty syndrome measurements, and can therefore lead to wrong error recovery inducing logic errors, but even then, we will find a significant improvement of the surface-4-GKP code over the surface-GKP code. All possible Pauli errors induced by GKP quadrature correction, and their effect on the stabilizer measurements, are described in appendix C. Note, in the case here with GKP quadrature correction after every gate, having −1 coupling rate inĈ X (−1) of step 2 and 3 is unnecessary, as all CV noise is immediately corrected. However, since a −1 coupling rate requires no extra resources and is solely controlled by the basis settings in table I, we keep it like this to compare with the surface-GKP code.
Finally, the surface code Z-and X-stabilizer measurement outcomes from d rounds of syndrome measurements are recorded in the vertices of 3D Z and X matching graphs with edges corresponding to possible Pauli errors as described in appendix C. For simplicity, here we consider only edges corresponding to single uncorrelated Pauli errors, and ignore possible, but less likely, two-qubit Pauli errors correlated by two-mode gates. Minimum-weight perfect matching (MWPM) [53,54] on these matching graphs is then used as the decoding algorithm to determine data qubit errors and the resulting error recovery. In practice, the error recovery is simply handled by using and updating a Pauli frame [55][56][57][58], similar to how feed-forward can be handled in MBQC by compensating for by-products in the following measurement outcomes [2,3]. For the MWPM to find error paths of highest probability, the edges of the matching graphs are dynamically weighted using Eq. (7) with the residual analogue information from each GKP quadrature correction. In this way, we can infer the probabilities of having induced the Pauli errors represented by each edge (described in appendix C). With each edge representing multiple Pauli errors induced in different GKP quadrature corrections, multiple Pauli error probabilities are combined in each edge weight as where p i is the probability given by Eq. (7) for one GKP quadrature correction taking values between 0 (no error) and 1/2 (minimal error information). For the combined probabilities of multiple edges to correctly add up in an error path determined by the MWPM, the edge weights in the matching graphs are finally taken to be log 2 (p tot ).

B. Simulation results
To establish a fault-tolerant error threshold, we numerically simulate the complete scheme. The GKPencoded data and measure qubits and the qunaught states, |∅ GKP , are all initialized with σ 2 GKP = σ 2 = e −2r /2 variance of the wave functions' GKP spikes as described in section III. The ancillary squeezed vacuum states for gate implementation, |0 sq , are as well squeezed by σ 2 = e −2r /2 leading to quadrature-symmetric gate noise of variance σ 2 gate = 2σ 2 = e −2r as described in section II. Using the Monte Carlo method, logical qubit error rates are simulated as a function of squeezing using Squeezing threshold of the surface-4-GKP code as a function of the probability with which a |∅ GKP -state is replaced with a |0 sq -state in the resource preparation. Here, the threshold is estimated as the crossing point of the d = 7 and d = 9 logic error rates. For zero replacement probability, the threshold is that of (a).
up to 100 000 simulation samples with a stopping condition at the occurrence of 500 combined logicẐ andX qubit error events. The resulting logicalẐ orX error rate (they are equal) is shown in Fig. 5(a) for different code distances d as a function of squeezing level, while the logicalŶ error rate is smaller. The decibel scale is defined relative to the vacuum variance, 10 log 10 [σ 2 /(1/2)]. The resulting squeezing threshold from where the logic error rate decrease with increasing code distance is found to be 12.7 dB of squeezing.
For comparison, in appendix B we also simulate the error rates of other scenarios with the simulation results shown in Fig. 6. For the surface-GKP code with a single GKP quadrature correction before the surface code sta-bilizer measurements, the squeezing threshold increases to 17.3 dB. This is significantly higher than the 12.7 dB squeezing threshold of the surface-4-GKP code due to accumulation of gate noise during the stabilizer measurements. To compare with other MBQC schemes with topological error correction where gate noise is typically not taken into account, and so only includes finite squeezing noise from GKP states, we simulate the surface-4-GKP code with σ 2 gate = 0. The resulting threshold is 10.2 dB of squeezing which agrees well with the 10 dB reported in [18] and the 10.5 dB reported in [21]. Finally, to see the impact of using the residual analogue information of the GKP quadrature correction in the weighting of the matching graphs for MWPM decoding, we simulate the surface-4-GKP code with fixed weighting based on variances of each mode at each point in the code similar to [29]. As expected, the result is a slightly larger squeezing threshold of 13.6 dB.
While the GKP code has been experimentally realized in trapped-ion and circuit QED systems [59][60][61], GKP state generation in optical platforms is yet to be demonstrated, although there are many recent proposals [45,[62][63][64][65][66][67]. At first, GKP state generation is most likely going to be probabilistic. In Ref. [21] it is proposed to combine multiple GKP state generators with optical switches, and then switch between generators with a successful preparation of a GKP state. In this way, the success probability of the GKP state generation, p ∅ , can in principle be brought arbitrarily close to 1. Since the surface-4-GKP code requires a large supply of |∅ GKP states, we consider as our final analysis the multi-GKP state generation scheme of Ref. [21] for |∅ GKP resource state preparation. If all the generators fail to prepare a |∅ GKP state in a given temporal mode for GKP correction, a deterministically generated squeezed vacuum state, |0 sq , is used instead. In this case, if |∅ GKP is replaced by |0 sq in spatial mode A(B), only theq(p) quadrature is corrected in the GKP quadrature correction, while the other quadrature accumulates gate noise of variance σ 2 gate during the correction [28]. The resulting fault-tolerance squeezing threshold is shown in Fig. 5(b) as a function of the probability 1−p ∅ of replacing |∅ GKP states by |0 sq states. Here, a squeezing threshold is seen to exist in a large range of 1 − p ∅ > 0, allowing for a probabilistic supply of |∅ GKP states, while for increasing replacement probability, the squeezing threshold level increases as expected-above 28 dB squeezing it is hard for us to simulate the squeezing threshold due to low error rates of less than 10 −7 . Note, here we still assume successful encoding of the surface code. I.e., the data and measure qubits switched into the setup as |ψ in GKP in Fig. 2(a) are successfully prepared as GKP qubit states. With probabilistic optical GKP state generation, this may be possible using state storage of a probabilistically prepared GKP state until it is switched into the computation scheme [68][69][70][71].

V. DISCUSSION
The squeezing thresholds in this work are derived by assuming a particular noise model in which all resource states are finitely squeezed while all optical propagation and detection losses are set to zero. In practice, however, losses cannot be neglected. Let us denote the transmission of the setup by η. For Gaussian states, 0 < η < 1 leads to the formation of mixed states with reduced effective squeezing. This can be reformulated as an ideal, lossless setup (η = 1) with mixed squeezed vacuum states as input having a lower effective squeezing, and some excess anti-squeezing that does not affect the measurementbased computation [72]. As a result, for η < 1 the 12.7 dB squeezing threshold corresponds to the effectively measured squeezing. Now, for the GKP states, besides a Gaussian convolution in the quadratures, η < 1 leads to a "shrinking" of a GKP state in phase-space. To see this, consider the Heisenberg picture with η modeled as a beam-splitter of η transmission. In this case, an amount 1 − η of vacuum is mixed into the state, adding noise to the quadratures, while a share 1 − η of the state is lost, "shrinking" the state in the quadratures by √ η. The quadrature shrinking is more detrimental to GKP spikes far from the phase-space origin, which are naturally delimited in GKP states of finite squeezing due to the overall envelope in the quadrature wave function. For GKPstates with 12.7 dB of squeezing, we assume this effect to be negligible on the qubit error probabilities for reasonably high efficiencies-we estimate η 0.95 to be doable on optical platforms. We also note that the shrinking effect can be counteracted by linear amplification which on the other hand will further reduce the amount of squeezing [41][42][43][44], effectively resetting η to unity at the cost of lowering the effective squeezing and the purity of the GKP state. Again, the estimated threshold of 12.7 dB refers to the required squeezing after such actions have been implemented. Another detrimental effect that has not been directly accounted for is interferometric phase fluctuations. Similar to optical loss, phase fluctuations lead to mixed squeezed states of reduced squeezing and excess anti-squeezing as well as mixed GKP states with an impact that increases with the quadrature value. Finally, we comment on the scalability of the computation scheme. For the temporal encoding in Fig. 2, the number of modes in which GKP qubits can be encoded for computation, i.e. the size of the encoding plane in the 3D time lattice, depends on the nm-delay in the resource preparation gadget. Increasing the delay length increases the number of encoding modes. However, doing so also increases the optical propagation loss, which puts a limit on the useful delay length. Thus, to continue scaling up, nm must be increased by shortening the temporal modes, in turn increasing the demands on the squeezing and detection bandwidth. In [73], squeezed light with a bandwidth of 2.5 THz was demonstrated, limited by the phase-matching condition of the non-linear down conversion process, while in [74], detection of squeezing up to 3 THz sideband frequency was demonstrated. Assuming proper squeezing, experimental control, and detection in a 2.5 THz bandwidth defining temporal modes of ∼ 1/2.5 THz duration, and assuming a propagation efficiency above 0.95 (0.23 dB attenuation) in a low-loss optical fiber with low optical attenuation of 0.15 dB/km, up to nm ≈ 10 7 computation modes may be realized in the temporally encoded computation scheme. For the spatial architecture of Fig. 1, scalability is similar to other schemes based on spatial encoding. It relies on the availability of resources, and is suitable for integrated photonics [75]. Finally, temporal and spatial encoding may be combined: Consider multiple temporally encoded computational devices, each as in Fig. 2. Using the optical switch at the setup computational level, computation modes can be switched in and out between different devices. Since the setup is optical, the devices are simply connected by optical fibers between the switches of each device without the need of quantum transducers. Furthermore, with the switch being mode selective, each mode of an encoded logical qubit in the surface code can be transferred without the need of decoding and re-encoding the logical quantum state, while measurement of the surface code stabilizers after transfer may be used for error-correcting the transfer line. This is not only suitable for combining temporal and spatial encoding for up-scaling, but is also useful in a quantum internet scheme [76,77], and is made possible by the optical architecture combined with temporal multiplexing on the transfer lines.

VI. CONCLUSION
In this work, we have proposed a simple but complete and scalable architecture for optical CV MBQC that includes quadrature noise correction and qubit error correction using topological codes. The setup consists of simple optical devices such as beam-splitters, delays, optical switches, and variable beam-splitters, where the latter two can be decomposed into beam-splitters and optical phase shifters. The scheme allows for both spatial and temporal encoding, with the temporally encoded version requiring just two squeezing sources. A universal Gaussian gate set is directly implementable, while universal qubit computation is made possible by feeding the setup with GKP states, thereby supplying the required non-Gaussianity [78][79][80]. As the computation scheme is based on gate teleportation on wires of two-mode entangled states, the setup naturally supports the new GKP quadrature correction scheme in Ref. [28], circumventing the need for on-line two-mode gates coupling to ancillary GKP states. Finally, by arranging the GKP qubits in a 2D plane of the cluster state that allows for nearestneighbor interactions, topological codes can be realized. By encoding a variation of the surface-GKP code-the surface-4-GKP code-we show fault-tolerant computation to be possible above a certain squeezing threshold by simulating a logic qubit memory, or an identity gate, of the surface-4-GKP code. In the surface code, Clifford gates can be implemented by braiding [25] or lattice surgery [81] implemented by regulating the surface code syndrome measurements, while non-Clifford gates may be realized using magic states distilled from GKP qubits prepared in a magic state [78,79] and injected into the surface code as input states [25,82].
The fault-tolerant squeezing threshold is found to be 12.7 dB. The estimation of this number takes into account the finite squeezing values of GKP states as well as the gate noise stemming from the finite squeezing values of the generated cluster state on which gates are implemented by projective measurements. However, this squeezing threshold leaves room for improvements: In the matching graph of the decoding algorithm, we have only considered single uncorrelated Pauli errors. As an improvement, considering the surface code to consist of two-mode gates, we may consider a matching graph taking two-qubit error events into account [83][84][85][86]. Furthermore, when weighting the matching graph edges with analogue information from GKP quadrature correction, we simply consider uncorrelated noise in the individual GKP corrections. However, since each mode is subjected to a two-mode gate prior to GKP correction leading to correlated noise on neighbouring GKP qubits, we may improve the estimation of Pauli error probabilities used in the matching graph weighting by jointly considering the analogue information from neighbouring GKP corrections as recently proposed in Ref. [87]. Note, that such modifications are solely implemented at the software level of the error correction decoder and thus requires no modifications to the setup. Another improvement may be found in the gate implementation: Due to the similarity of the GKP quadrature correction and gate implementation, it might be possible to combine the two transformations in one step, that is, implementing a gate while correcting the quadratures. Although the quadrature correction is only considered on single wires [28], it might be possible to generalize it to the two-wire case by which two-mode gates could be implemented during GKP quadrature correction, thereby eliminating gate noise coursed by finite squeezing. If this is possible while maintaining the GKP quadrature correction quality of Ref. [28], the resulting squeezing threshold reduces to 10.2 dB as shown in appendix B.
To derive the gates implemented by projective measurements in section II, consider one computation step on two parallel wires, where the red arrows represent the first beam-splitter in the TDMD, and the gray arrow represents the blue or green VBS in the TDMD for j = 1 or j = n, respectively. Here, a two-mode input state (separable or not) is encoded in modes (A, k), (A, k+j), while the two-mode entangled states are prepared in the resource preparation gadget, and can be written as where |0 q−p and |0 q+p are squeezed vacuum states, squeezed along (q−p)/ √ 2 and (q+p)/ √ 2 quadratures, respectively, and similarly, |0 q and |0 p are squeezed along theq andp, respectively. As such, the two-mode entangled states correspond to two-mode squeezed states rotated in phase-space by π/4, turning them into approximate cluster states with tanh 2r edge weight and isech2r self-loops, where r is the squeezing parameter of the initial squeezed vacuum states [19,33]. Alternatively, we can consider the two-mode entangled states more generally as cluster-type states [32], here with edge weight 1, for which the implemented gate is independent on the squeezing, r, which then only affects the gate noise [10]. The two situations are equivalent: One can change from the former to the latter by normalizing the edge weight [30]. Here, we will consider cluster-type states, since implementing a desired gate in practice (without considering the resulting gate noise) then requires no prior knowledge of the squeezing level.
For single-mode gates the VBSs of the TDMD are left open, and the dashed arrow in Eq. (A1) representŝ I A,k ⊗Î B,k+j . In this case, we can ignore the second wire, and focus on a joint projective measurement of the input mode (A, k) and one mode of the two-mode entangled state, (B, k), resulting in gate teleportation to the output mode (A, k + nm)-exactly the same derivation can be made on the second wire of mode (A, k +j), (B, k +j), and (B, k + nm + j). The corresponding circuit is where |0 q in (A, k + nm) is replaced byR(π/2) |0 p only to follow the traditional convention of cluster states with initial squeezing inp quadratures. Using the method of appendix A in Ref. [31] with (θ A,k , θ B,k ) being the basis setting determining the implemented gate, the corresponding quadrature transformation in the Heisenberg picture can be derived to be Here, G is the symplectic matrix corresponding to the desired single-mode gate operation in Eq. (2), is a gate noise matrix, and is a displacement matrix. Since (θ A,k , θ B,k ), m A,k and m B,k are known, D(m A,k , m B,k ) T can be compensated for by displacing the teleported state back by −D(m A,k , m B,k ) T , or simply by taking this displacement into account in the following measurement outcomes. With finite squeezing in the ancillary modes such that Var{p B,k } = Var{p A,k+nm } = σ 2 = e −2r /2, the noise term N(p B,k ,p A,k+nm ) T leads to quadrature-symmetric gate noise inq A,k+nm andp A,k+nm of In the Wigner function picture, this gate noise corresponds to convolutions in both quadratures by a Gaussian function of variance σ 2 gate , each followed by the application of a corresponding Gaussian envelope due to the Fourier relation betweenq andp [30,31].
When implementing two-mode gates by enabling the first or second VBS of the TDMD, the corresponding circuit is We do not derive a general expression for the implemented gate as a function of the basis setting (θ A,k , θ B,k , θ A,k+j , θ B,k+j ). Instead, we use the method described in Ref. [31]: A cost function is defined based on the implemented gate, a desired target gate, and the gate noise, which is then used in a global search to find the basis setting that implements a desired gate with the minimum gate noise. The resulting basis settings for the gates required to implement the surface code are shown in Table I. We note that those settings only implement two different gates since (F † ⊗F )Ĉ Z (g) =Ĉ X (g)(F † ⊗F ) andĈ Z (g)(F ⊗F † ) =Ĉ X (g)(F † ⊗F ). The reason for considering them as four different gates is to make the implementation of the surface code more intuitive. The basis settings are not unique: other settings exist that implement the same gates with equal gate noise.
The quadrature transformation when applying the basis settings for two-mode gates is

Again,
G is the symplectic matrix corresponding the implemented two-mode gate. To simulate the logic qubit error rate of the surface code, we adopted and modified the simulation in Ref. [29] to the computation scheme of this work. The simulation method is well-described in appendix B of Ref. [29] and is summarized here with focus on the modifications. In the simulation, quadrature noise is simulated as stochastic normally-distributed variables for each quadrature of each mode i, ξ i q and ξ i p . For GKP-states, ξ i q and ξ i p are initialized with random samples from N (0, σ GKP ), where N (0, σ) is a normal distribution of zero mean and σ 2 variance. After each gate, independent random samples from N (0, σ gate ) are added to ξ i q and ξ i p as gate noise. As for homodyne measurements, ξ i q or ξ i p is read out, and the logic value is determined from the closest integer multiple of √ π. Note that, unlike in Ref. [29], we do not consider measure noise or idle noise. In optical platforms, homodyne measurements are carried out with near-unity efficiency (any loss is assumed to just degrade the squeezing as discussed in section V). Furthermore, in MBQC no modes are idle since modes not performing any tasks still have to teleport through the computation step and thereby acquire gate noise instead of idle noise.
For the two-mode gates in the surface code, the simulation here differs from Ref. [29] by usingĈ Z (1) and C X (±1) gates instead of sum-gates. For a two-mode gate between modes i and j, the quadrature noise variables are updated aŝ , C X (±1) : where randG(σ 2 ) returns a random value from N (0, σ). For GKP quadrature correction, instead of coupling to ancillary GKP qubits through sum gates as in Ref. [29], the mode to be corrected is teleported through a twomode GKP qubit Bell state as described in section III. The Bell state is prepared by interfering two GKP qunaught states, denoted ∅1 and ∅2, To teleport, the mode to be corrected, i, and ∅1 are interfered on a beam-splitter and measured inq andp, respectively, with outcomes Finally, to compensate for the Pauli by-products of the qubit teleportation (displacements by √ π), m A √ 2 and m B √ 2 are rounded to the nearest integer multiple of √ π, which is then used to displace the teleportation output mode, ∅2, back, For the sake of simulation, we pass the corrected output mode to the input mode, ξ i q ← ξ ∅2 q and ξ i p ← ξ ∅2 p , such that mode i can be reused in the following simulation. The probability of having induced a qubit error by rounding to a wrong integer of √ π in Eq. (B1) due to input noise in ξ i q and ξ i p , together with initialization noise of ξ ∅1 q , ξ ∅1 p , ξ ∅2 q , and ξ ∅2 p , is inferred using the residual analogue information, 2), in Eq. (8) through the probability in Eq. (7). Finally, it is used for weighting the matching graphs of stabilizer measurement outcomes for the MWPM decoding as described in appendix C. Here, σ 2 in in σ 2 = σ 2 in + σ 2 GKP of Eq. (7) is the quadrature variance of the input mode, and is carefully kept track of in the simulation based on previous gates and corrections.
In one simulation, d + 1 rounds of surface code stabilizer measurements are carried out. Data GKP qubits are initialized in round 1 with σ 2 GKP variance. To stabilize the data qubits, measure GKP qubits and qunaught states are initialized in round 1 to d with σ 2 GKP variance, followed by noisy gates and measurements to build up the matching graphs. In the last round, d + 1, measure GKP qubits and qunaught states are initialized with zero variance to carry out ideal syndrome measurements for determining logic qubit errors induced in round 1 to d. To build up statistics, for each squeezing level and code distance, d, this process is repeated 100 000 times, or until a total of 500 logicX andẐ errors are detected.
We have simulated four different cases, shown in Fig. 6. In three cases, GKP states (qubits and qunaught states) and squeezed vacuum states are initialized with equal variance, σ 2 GKP = σ 2 = e −2r /2, which from Eq. (3) and (5) leads to σ 2 gate = 2σ 2 GKP . In this way, the surface-GKP code with GKP quadrature correction before the surface code stabilizer measurements, and the surface-4-GKP code with four GKP quadrature corrections during the stabilizer measurements, was simulated. To see the impact of using the analogue information from the GKP correction in the weighting of the matching graphs, the surface-4-GKP code was simulated using instead of Eq. (7) [29]. By integrating the wave function marginal distribution in the odd GKP bins, p err (σ) infers the qubit error probability only based on variances without taking the projective measurement outcome into account. Finally, to compare with other MBQC schemes supporting topological error correction, but only taking noise from GKP-states into account, the surface-4-GKP code is simulated using σ 2 gate = 0.

Appendix C: Decoding graphs
In this appendix, we describe the matching graphs used in the MWPM decoding. A section of the Z and X matching graphs is shown in Fig. 7. Each vertex corresponds to a syndrome measurement and is highlighted when the measurement outcome change relative to the previous measurement of the same measure qubit, indicating an error event. The edges correspond to possible errors, and the job of the decoder is to match pairs of highlighted vertices with the most likely error path, which is then translated into corrections of data qubits.
Each horizontal plane of the matching graphs in Fig. 7 corresponds to one round, i, of the syndrome measurements in Fig. 4(b,c). We distinguish between 4 different types of edges: horizontal edges, h; vertical edges, v; diagonal edges, d; and cross edges, c. In the case where errors only occur in between syndrome measurements (i.e. before and after the syndrome measurement circuits in Fig. 4(b,c)), errors on data qubits correspond to h-edges, while errors in the syndrome measurement readout correspond to v-edges. These are often the only edges included in the surface code matching graphs when errors during the syndrome measurements are not considered. However, for the surface-4-GKP code with GKP quadrature correction during the syndrome measurements, most qubit errors are induced during the syndrome measurements, requiring the additional d-and c-edges (often referred to as space-time edges) in the matching graphs for optimal decoding [83].
Below we go through all possible qubits errors which may be induced by each GKP quadrature correction in each of the four steps in the syndrome measurements and describe the corresponding edge in the matching graphs. For simplicity, in this work, we consider only uncorrelated single-qubit errors. Since the surface code consists of two-mode gates, two-qubit errors are possible as well, but we assume them to be negligible in the considered squeezing range. The effect of two-qubit errors on the matching graphs is shown in [83]. Finally, in section C 5 we comment on the graph boundaries, edge weighting, and qubit correction.

Step 1
After the first set of two-mode gates in the syndrome measurements (step 1), we perform GKP quadrature correction. At this point, below we consider singlequbit Pauli-X errors,X = e ±i √ πp , and Pauli-Z errors, Z = e ±i √ πq , (corresponding to √ π displacement inq and p, respectively) in each data and measure qubit. In all the following we will refer to the j'th data, measure-Z, and measure-X qubit as D j , Z j , and X j , and we will use the labeling in Fig. 7 for the distance d = 5 code as an example.
For an odd data qubit, say D 13 , anX error is detected by Z 5 in step 4 of the current round i, while it will be detected by Z 8 in step 1 of the following round i + 1. As a result, this error corresponds to the edge d i 13 of the Z-graph. AẐ error will be detected in X 8 and X 5 in For an even data qubit, say D 8 , anX error is detected by Z 2 and Z 5 in steps 2 and 3 of the current round i, and corresponds to h i 8 of the Z-graph. AẐ error is detected by X 7 in step 4 of the current round i, and by X 5 in step 1 of the following round i+1. Thus this error corresponds to d i 8 of the X-graph. For the measure-Z qubit Z 5 , anX error propagates to D 8,13,14 throughĈ Z (1) in steps 2, 3, and 4 where they cause aẐ error, which eventually will highlight X 7 and X 11 in step 4 of this round i, indicating a false error on D 9 instead of on D 8,13,14 . However, applying aẐ correction on D 9 constitutes together with theẐ errors on D 8,13,14 a code stabilizer, and the error is successfully corrected. Thus, anX error on Z 5 corresponds to h i 9 of the X-graph. AẐ error on Z 5 does not propagate through theĈ Z (1) gates, and becomes instead a detection error at the end of the syndrome measurement, corresponding to v i 5 of the Z-graph.
For the measure-X qubit X 5 , anX error does not propagate through theĈ X (±1) gates, and becomes a detection error at the end of the syndrome measurement, corresponding to v i 5 of the X-graph. AẐ error propagates asX errors to D 7,12,13 throughĈ X (±1) gates, which eventually highlights Z 2 and Z 5 in step 4 of this round, indicating a false error on D 8 . Applying anX correction on D 8 constitutes together with theX errors on D 7,12,13 a code stabilizer, and the error is successfully corrected. Thus aẐ error on X 5 corresponds to h i 8 of the Z-graph.

Step 2
Below we consider single-qubit errors induced in the GKP quadrature correction after the second set of twomode gates in the syndrome measurements (step 2).
AnX error on an odd data qubit is detected similarly to anX error induced in step 1. AẐ error on an odd data qubit, say D 13 , is detected by X 5 in step 3 of this round i, while detected by X 8 in step 2 of round i + 1. Thus the corresponding edge is d i 13 of the X-graph. AnX error on an even data qubit, say D 8 , is detected by Z 5 in step 3 of this round i, while detected by Z 2 in step 2 of round i + 1, and the corresponding edge is d i 8 of the Z-graph. AẐ error on an even data qubit is detected in the same way as aẐ error induced in step 1.
AnX error on the measure-Z qubit Z 5 propagates throughĈ Z (1) to D 8 and D 13 in steps 3 and 4 asẐ error. TheẐ error on D 8 will be detected by X 7 in step 4 of this round i, while theẐ error on D 13 will be detected by X 8 in step 2 of round i + 1. As a result, the corresponding edge is c i 5 of the X-graph. Similar to step 1, aẐ error corresponds to a detection error, i.e. v i 5 of the Z-graph.
AnX error on measure qubit X 5 , similar to step 1, corresponds to a detection error, i.e. v i 5 of the X-graph. AẐ error on X 5 propagates throughĈ X (±1) to D 12 and D 13 in steps 3 and 4 asX error. TheX on D 13 will be detected at Z 5 in step 4 of this round i, while theX error on D 12 will be detected by Z 4 in step 2 of the following round i + 1. As a result, the corresponding edge is c i 5 of the Z-graph.

Step 3
Below we consider single-qubit errors induced in the GKP quadrature correction after the third set of twomode gates in the syndrome measurements (step 3).
AnX error on an odd data qubit is detected similarly to anX error induced in steps 1 and 2. AẐ error on an odd data qubit, say D 13 is detected by X 8 and X 5 in steps 2 and 3 in the following round i + 1, and the corresponding edge is h i+1 13 of the X-graph. AnX error on an even data qubit, say D 8 , is detected by Z 2 and Z 5 in steps 2 and 3 of the following round i+1, and the corresponding edge is h i+1 8 of the Z-graph. AẐ error on an even data qubit is detected in the same way as aẐ error induced in steps 1 and 2.
AnX error on measure-Z qubit Z 5 propagates througĥ C Z (1) to D 13 in step 4 as aẐ error and is detected by X 8 and X 5 in steps 2 and 3 of the following round i+1. Thus the corresponding edge is h i+1 13 of the X-graph. Similar to steps 1 and 2, aẐ error corresponds to a detection error, i.e. v i 5 of the Z-graph. AnX error on measure qubit X 5 , similar to steps 1 and 2, corresponds to a measure, i.e. v i 5 of the X-graph. AẐ error on X 5 propagates throughĈ X (1) as anX error to D 12 in step 4, where it will be detected by Z 4 and Z 8 in steps 2 and 3 of the following round i + 1. Thus the corresponding edge is h i+1 12 of the Z-graph.

Step 4
Errors on data qubits, D j , induced in the GKP quadrature correction after the final fourth set of two-mode gates in the syndrome measurements (step 4) will only be detected in the following round i + 1, and the corresponding edges are h i+1 j of the X-and Z-graph. Since measure qubits are now to be measured after this fourth step of the syndrome measurements, they do not go through GKP quadrature correction, and we do not induce any qubit error. Measure-Z qubits are measured in thep quadrature, and a wrong measurement outcome (caused by finite squeezing of the GKP spikes in the measured quadrature) corresponds to v i j in the Zgraph. Similar, measure-X qubits are measured in thê q quadrature, and a wrong measurement outcome corre-sponds to v i j in the X-graph.

Graph boundaries, weighting, and correction
As described in section III, IV A, and appendix B, from the residual analogue information in the GKP quadrature correction, we can infer the probability of having induced a qubit error by the GKP quadrature correction, which we then use for weighting of the corresponding abovedescribed edges in the matching graphs. With each edge corresponding to several different errors, for each edge, the probabilities from multiple GKP quadrature corrections are combined into p tot by Eq. (10), and the corresponding edge is finally weighted by log 2 (p tot ).
At the boundary of the code, not all diagonal edges are used. Instead, the corresponding horizontal edge in the following round i + 1 is used. As an example, consider anX error induced on D 3 in the distance d = 5 code (see Fig. 4(a) for labeling). This error will be detected by Z 2 in the following round i + 1, but there are no other measure-Z qubits to detect the error in the current round i. Thus, instead of weighting d i 3 , in this case, the error probability is included in the weight of h i+1 3 . In the final round i = d + 1 where an ideal syndrome measurement is performed, all data qubits are first corrected with perfect GKP quadrature corrections (using |∅ GKP -states of infinite squeezing) before measuring the stabilizers of the surface code using measure qubits of GKP-states with infinite squeezing. As a result, this final round only requires horizontal edges.
After the MWPM decoding, vertices in the matching graphs are matched with error paths following edges with the most likely errors based on the edge weighting. A vertical edge in an error path, corresponding to a detection error, requires no correction on data qubits. A horizontal or diagonal edge, h i j or d i j , in an error path in the Z or X graph requires anX orẐ correction on data qubit D j , respectively. Finally, a cross edge, c i j , in an error path in the Z or X graph requires twoX orẐ corrections on neighboring data qubits, respectively: As an example, including c i 5 of Fig. 7(a) in an error path requireŝ X corrections on D 7 and D 8 , or on D 12 and D 13 . Similarly, including c i 5 of Fig. 7(b) in an error path requireŝ Z corrections on D 8 and D 13 , or on D 9 and D 14 .