{'fraction_non_alphanumeric': 0.04891101431238332, 'fraction_numerical': 0.013011257566329128, 'mean_word_length': 4.02535819877189, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 42, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Entanglement Purification with Quantum LDPC Codes and Iterative Decoding

Recent constructions of quantum low-density parity-check (QLDPC) codes provide optimal scaling of the number of logical qubits and the minimum distance in terms of the code length, thereby opening the door to fault-tolerant quantum systems with minimal resource overhead. However, the hardware path from nearest-neighbor-connection-based topological codes to long-range-interaction-demanding QLDPC codes is likely a challenging one. Given the practical difficulty in building a monolithic architecture for quantum systems, such as computers, based on optimal QLDPC codes, it is worth considering a distributed implementation of such codes over a network of interconnected medium-sized quantum processors. In such a setting, all syndrome measurements and logical operations must be performed through the use of high-fidelity shared entangled states between the processing nodes. Since probabilistic many-to-1 distillation schemes for purifying entanglement are inefficient, we investigate quantum error correction based entanglement purification in this work. Specifically, we employ QLDPC codes to distill GHZ states, as the resulting high-fidelity logical GHZ states can interact directly with the code used to perform distributed quantum computing (DQC), e.g. for fault-tolerant Steane syndrome extraction. This protocol is applicable beyond the application of DQC since entanglement distribution and purification is a quintessential task of any quantum network. We use the min-sum algorithm (MSA) based iterative decoder with a sequential schedule for distilling 3-qubit GHZ states using a rate 0.118 family of lifted product QLDPC codes and obtain a threshold of 10.7% under depolarizing noise. Our results apply to larger size GHZ states as well, where we extend our technical result about a measurement property of 3-qubit GHZ states to construct a scalable GHZ purification protocol.

Bell pairs can be very large (and uncertain). While their hand-designed purification schemes have been extended by algorithmic procedures recently [9], [10], the approach still suffers from this inefficiency arising from the heralded nature of the protocol.
Our goal in this paper is to investigate a principled and systematic procedure to purify (or distill) GHZ states using quantum error correcting codes (QECCs). If one can use the same QLDPC codes that DQC will employ for FTQC ("compute code") to also store logical GHZ ancillary states, then these can potentially be directly interacted with the compute code for performing fault-tolerant (e.g., Steane) syndrome extraction and measurement-based methods for logical operations. Thus, it is very pertinent to develop a scalable GHZ purification protocol using these optimal QLDPC codes ("purification code"). While DQC is a key motivation, such a protocol serves a much wider purpose, since entanglement generation, distribution and purification form the cornerstone of quantum networking. For efficient and scalable quantum networks, one must necessarily deploy quantum repeaters whose primary function is to help entangle different subsets of parties in the network. In the long-run, third generation quantum repeaters will use quantum error correction for entanglement purification [11]. Such repeater nodes, and other nodes of the network that are not quantum computers, will still need to possess a fault-tolerant quantum memory to generate and store (shares of) high-fidelity entangled states. Therefore, if the compute nodes will deploy QLDPC codes, then QLDPC purification codes could potentially unify the functioning of different parts of the network and enable seamless integration.
Entanglement purification is a well-studied problem in quantum information, where one typically starts with n copies of a noisy Bell pair, or a general mixed state, and distills k Bell pairs of higher fidelity [12]. Several teams of researchers have worked on this problem, and the contributions range from fundamental limits [12]- [17] to simple and practical protocols [9], [13], [18], [19]. Of course, if one can distill Bell pairs, then these can be "fused" in sequence to entangle multiple parties, but direct distillation of an entangled resource between all parties can be more efficient [20]. Some purification schemes involve two-way communications between the involved parties while others only need one-way communication. We focus on oneway schemes in this paper. The connection between one-way entanglement purification protocols (1-EPPs) and QECCs was established by Bennett et al. in 1996 [13]. They showed that any QECC can be converted into a 1-EPP (and vice-versa). This framework enables systematic n-to-k protocols where the rate and average output fidelity are directly a function of the QECC rate and decoding performance, respectively. Since the recently constructed QLDPC codes have asymptotically constant rate and linear distance scaling with code size [1]- [6], our work paves the way for high-rate high-fidelity entanglement distillation.

A. Purifying Bell Pairs with QLDPC Codes
In 2007, Wilde et al. [18] showed that any classical convolutional code can be used to distill Bell pairs via their entanglement assisted 1-EPP scheme. In the development of this scheme, they mention a potentially different method to use a QECC for performing 1-EPP [ Protocol I for GHZ purification using stabilizer codes. Alice generates n copies of the ideal 3-qubit GHZ state and marks one qubit of each triple as 'A', another as 'B', and the third as 'C'. She measures the stabilizers of the code on the n qubits marked 'A' and classically communicates the results through a noiseless channel to Bob. She also sends all the remaining 2n qubits to him. First, Bob uses the "GHZ-map" to measure stabilizers of the 2n-qubit code induced by Alice's code and uses the results to correct errors on qubits 'B' and 'C'. Second, if the code is not CSS, then he performs a suitable diagonal Clifford U C on qubits 'C'. Thirdly, he measures the stabilizers of Alice's n-qubit code on qubits 'B' in order to impose the same code on qubits 'B' and 'C'. Lastly, he communicates the results to Charlie over a noiseless channel and also sends qubits 'C' to him. Charlie measures the stabilizers of Alice's code to detect and correct errors on qubits 'C'. Finally, all three parties invert the encoding to convert the logical GHZ states to physical GHZ states. The scheme is suited for a linear network of three parties, but notice the asymmetrically larger burden on Bob, which makes the protocol less scalable to larger GHZ states.
syndrome, shares her code stabilizers and syndrome with Bob through a noiseless classical channel, and sends qubits 'B' to Bob over a noisy Pauli channel. Using the transpose property, Bob measures the appropriate code stabilizers on qubits 'B', and corrects channel errors by combining his syndrome with Alice's syndrome. Finally, Alice and Bob decode their respective qubits, i.e., invert the encoding unitary (see Section I-B), to convert the k logical Bell pairs into k physical Bell pairs. Since the code corrects some errors, on average the output Bell pairs are of higher fidelity than the initial n noisy ones. As our first contribution, we elucidate this protocol for general stabilizer codes [21], [22] through the lens of the stabilizer formalism [23], using the 5-qubit perfect code [21], [24] as an example. This approach clarifies many details of the protocol, especially from an error correction standpoint, and helps adapt it to different scenarios. For the performance of the protocol, note that any error on Alice's qubits can be mapped into an equivalent error on Bob's qubits using the transpose property, in effect increasing the error rate on Bob's qubits. Therefore, since only Bob corrects errors in this protocol, the failure rate of the protocol is the same as the logical error rate (LER) of the code on the depolarizing channel, with an effective channel error rate that accounts for errors on Alice's qubits as well as Bob's qubits (as long as they do not amount to a Bell state stabilizer together). For all simulations in this work, we consider a rate 0.118 family of lifted product (LP118) QLDPC codes decoded using the sequential schedule of the min-sum algorithm (MSA) based iterative decoder with normalization factor 0.8 and maximum number of iterations set to 100 [25], [26]. The LER of this code-decoder pair is shown in Fig. 1, where we see that the threshold is about 10.6-10.7%. Also, even with n = 544, the LER is ≈ 10 −6 at depolarizing rate 10 −2 . Again, note that these curves can be interpreted as the performance of Bell pair purification when only Bob's qubits are affected by errors.

B. New Protocols to Purify GHZ States with QLDPC Codes
Given these insights, we proceed to investigate the purification of GHZ states. As in the Wilde et al. protocol, we consider only local operations and one-way classical communications (LOCC), and assume that these are noiseless. The key technical insight necessary to construct the protocol is the GHZ-equivalent of the transpose property of Bell pairs. Given n copies of the GHZ state, whose three subsystems are marked 'A', 'B' and 'C', we find that applying a matrix on qubits 'A' is equivalent to applying a "stretched" version of the matrix on qubits 'B' and 'C' together (see Lemma 3). We call this mapping to the stretched version of the matrix the GHZ-map, and prove that it is an algebra homomorphism [27], i.e., linear, multiplicative, and hence projector-preserving. Recollect from the Bell pair purification setting that we are interested in measuring stabilizers on qubits 'A' and understanding their effect on the remaining qubits. Using the properties of the GHZ-map, we show that it suffices to consider only the simple case of a single stabilizer. With this great simplification, we prove that imposing a given [[n, k, d]] stabilizer code on qubits 'A' simultaneously imposes a certain [[2n, k, d ]] stabilizer code jointly on qubits 'B' and 'C'. By performing diagonal Clifford operations on qubits 'C', which commutes with any operations on the other qubits, one can vary the distance d of the induced 'BC' code. Then, we use this core technical result to devise a natural protocol that purifies GHZ states using any stabilizer code ("Protocol I", see Fig. 2).
We perform simulations on the [ [5,1,3]] perfect code and compare the protocol failure rate to the LER of the code on the depolarizing channel, both using a maximum-likelihood decoder. In terms of error exponents, we show that it is always better for Bob to perform a local diagonal Clifford operation on Charlie's qubits, rather than Alice doing the same. We support the empirical observation with an analytical argument on the induced BC code and Charlie's code. Finally, we finish by showing that the average output k-qubit density matrix of the protocol is diagonal in the GHZ-basis, and its fidelity is directly dictated by the protocol's failure rate. While the scheme is suited for a linear network of three parties, it is obvious that there is an asymmetrically larger burden on Bob, which makes the protocol less scalable to larger GHZ states. Nevertheless, we think that this protocol still has pedagogical value in understanding the implications of the new insight on GHZ states.
Motivated by this drawback, we devise an improved protocol that avoids the additional 2n-qubit measurements of Protocol I. The new scheme ("Protocol II") is depicted in Fig. 3 for CSS codes. The protocol can be extended to general stabilizer codes through additional diagonal Clifford operations as in Protocol I but, for simplicity, we focus on CSS (QLDPC) codes here. It will also be interesting to investigate if there are any potential gains from employing non-CSS stabilizer codes in entanglement purification, because CSS codes are known to be optimal for certain aspects of fault-tolerant quantum computing [28]. When Alice measures stabilizers on qubits 'A', the new GHZ property still implies that there is a 2n-qubit code automatically induced on qubits 'B' and 'C' together. In order to split that code into individual codes on qubits 'B' and 'C', Alice performs a second round of the same (n-qubit) stabilizer measurements but this time on qubits 'B'. This enables Bob and Charlie to measure the same stabilizers on their respective qubits and correct errors induced by the channel on qubits 'B' and 'C', respectively. The flow of the protocol is naturally applicable when Alice is connected to both Bob and Charlie but those two parties are not connected directly. But we emphasize that the protocol is scalable and we summarize its extension to larger GHZ states with larger number of parties connected by any network topology; the key requirement is that the qubits of a recipient over a network edge have already been measured and projected to the code subspace before those qubits are sent over the edge.
In Fig. 4 we report simulation results for Protocol II on 3-qubit GHZ states using the same LP118 code family and MSA decoder as in Fig. 1. All data points except the first one on each curve (for depolarizing rate 0.09) were computed by collecting close to 10 4 logical errors. We observe that the threshold is very close to the single decoder case in Fig. 1, which is reassuring since the GHZ protocol needs both Bob and Charlie to run decoders. However, the logical error rates are significantly higher than those in Fig. 1. This is likely due to the fact that both decoders must succeed for the protocol to not fail. Note that there can be situations where an error on Alice's qubit cancels the errors on Bob's and Charlie's due to the new GHZ property. But it is unclear whether these have a significant effect on the protocol performance. We plan to study this carefully in future work because it is undesirable for failure rates to increase as we scale the protocol to larger number of parties.

C. Decoding QLDPC Codes under Realistic Noise Models
While these results are relevant to the "code capacity" error model, where there are only qubit errors and all operations are assumed noiseless, in a separate work a subset of the authors considered decoding this family of QLDPC codes under a "phenomenological" noise model, i.e., with an additional (classical) error model on the syndromes [26]. In that setting, motivated by practical situations, the syndromes extracted from a measurement circuit are assumed to have an additional random Gaussian noise, thereby yielding "soft" syndromes. It was shown then that the MSA decoder can be modified appropriately such that the decoding performance is almost as good as the above ideal syndrome scenario. Therefore, by reinterpreting that work in the context of entanglement purification, we highlight that the protocol can be applied to more realistic settings as well.
Since we are constructing a new GHZ purification protocol based on this new insight about GHZ states, we have considered this simple model of noiseless LOCC and noisy qubit communications. We emphasize here that, to the best of our knowledge, this is the first protocol to use quantum error correction for purifying GHZ states, and we also report simulation results of state-of-the-art QLDPC codes with an efficient iterative decoder. While the problem of noisy local operations is important and has received attention [8], [9], [15], [29], we leave these investigations on our protocol for future work.

D. Purification-Inspired Algorithm to Generate Logical Pauli Operators
In the process, inspired by stabilizer measurements on Bell/GHZ states, we have developed a new algorithm to generate logical Pauli operators for any stabilizer code (see Algorithm 3 and its explanation in Appendix IV-B). The core idea is to first simulate the generation on n Bell/GHZ states by creating a table of their 2n stabilizers. It turns out that we only need the ZZIand XXX-type stabilizers for the GHZ case, which is why we ignore the n IZZ-type stabilizers. Then we simulate the measurement of each code stabilizer on qubits 'A' using the stabilizer formalism. At the end of this process, it can be shown that the non-code-stabilizer rows in the table must be a combination of logical Pauli operators on multiple subsystems. Finally, we carefully identify the logical Pauli operators on qubits 'A' and return those as the desired operators on the given code. Fig. 3: Protocol II for GHZ purification using CSS codes. The protocol can be extended to general stabilizer codes through additional diagonal Clifford operations as in Protocol I. Alice generates n copies of the ideal 3-qubit GHZ state and marks one qubit of each triple as 'A', another as 'B', and the third as 'C'. She measures the stabilizers of the QLDPC code on qubits 'A' and classically communicates the results through a noiseless channel to Charlie. She also uses the results to appropriately measure stabilizers of (a potentially equivalent) QLDPC code on qubits 'B' and communicates these results to both Bob and Charlie, again through noiseless classical channels. She sends qubits 'B' to Bob and qubits 'C' to Charlie. Finally, both Bob and Charlie make stabilizer measurements, correct errors, and then all three parties invert the encoding to convert the logical GHZ states to physical GHZ states. Note that Bob and Charlie can perform their operations asynchronously.  III. NOTATION AND BACKGROUND The Pauli group on n qubits is denoted by P n . We denote Pauli matrices I, X, Y, Z and their tensor products using the notation E(a, b), where a, b ∈ {0, 1} n denote respectively the Xand Z-components of the n-qubit Pauli operator. The weight of a Pauli operator is the number of qubits on which it acts nontrivially (i.e., does not apply I). For example, E([0, 1, 0, 1], [0, 0, 1, 1]) = I ⊗ X ⊗ Z ⊗ Y ≡ IXZY has weight 3 and we dropped the tensor product symbol ⊗ for brevity. Two Pauli operators E(a, b), E(c, d) either commute or anticommute, and this is dictated by the symplectic inner product in the binary vector space. If [a, b], [c, d] s := ad T + bc T = 0 (resp. 1) (mod 2), then they commute (resp. anticommute).

Charlie
A stabilizer group S is generated by commuting Pauli operators ε i E(a i , b i ), i = 1, 2, . . . , r, where ε i ∈ {±1} and −I 2 n / ∈ S. The [[n, k, d]] stabilizer code defined by S is given by Q(S) = {|ψ ∈ C 2 n : g |ψ = |ψ ∀ g ∈ S}, where k = n − r. The logical Pauli operators of the code commute with all stabilizers but do not belong to S, and their minimum weight is d. The code is completely defined by its stabilizers and logical operators, or equivalently by an encoding circuit U Enc (S). The projector onto the code subspace is given by A CSS (Calderbank-Shor-Steane) code is a special type of stabilizer code for which there exists a set of stabilizer generators such that each generator is purely X-type, i.e., of the form E(a i , 0), or purely Z-type, i.e., of the form E(0, b j ). Such a code can be described by a pair of classical binary linear codes C X and C Z , where the rows of the parity-check matrix H X (resp. H Z ) for C X (resp. C Z ) are a i ∈ {0, 1} n (resp. b j ∈ {0, 1} n ). Since E(a i , 0) and E(0, b j ) must commute, the symplectic inner product constraint leads to the condition a i b T j = 0 for all i, j or, equivalently, H X H T Z = 0. A quantum (CSS) low-density parity-check (QLDPC) code is described by a pair (C X , C Z ) of classical LDPC codes, which implies that H X and H Z are sparse, i.e., each stabilizer involves few qubits and each qubit is involved in few stabilizers. It is very challenging to construct good QLDPC codes due to the constraint H X H T Z = 0 on two sparse matrices, but recent exciting work has developed optimal QLDPC codes where k and d scale linearly with n [1]- [6]. For our simulations, we chose a specific family of lifted product QLDPC codes from [25, Table II] that have asymptotic rate k/n = 0.118. To decode these codes, we use the computationally efficient min-sum algorithm (MSA) based iterative decoder under the sequential schedule [30], [31], with a normalization factor of 0.8 and maximum number of iterations set to 100 (also see the description in [25]).
A stabilizer state corresponds to a code with dimension k = 0, and can equivalently be represented by a maximal stabilizer group, i.e., with r = n. Any Pauli measurement on the state can be simulated by a well-defined set of rules to update this stabilizer group. These rules are given by the stabilizer formalism for measurements [23], [32].
For any matrix M , the Bell state |Φ AB = This property extends to n copies of the Bell state as well. When M is a projector, which is the case when we perform stabilizer measurements on qubits 'A', i.e., M = Π S , using the fact that M 2 = M we conclude from the above property that projecting qubits 'A' automatically projects qubits 'B' as well according to M T . Therefore, imposing a code on qubits 'A' simultaneously imposes the "transpose" code on qubits 'B'.
A more detailed discussion of these background concepts can be found in Appendix I.
IV. REVISITING THE BELL PAIR DISTILLATION PROTOCOL In Ref. [18], Wilde et al. described a protocol to distill Bell pairs using an arbitrary quantum stabilizer code. We reiterate this protocol here and provide more clarity on the reasons behind its working. Then, in the next section, we will generalize this protocol to distill GHZ states, i.e., the -qubit entangled state GHZ = 1 √ 2 (|00 · · · 0 + |11 · · · 1 ).
Initially, Alice generates n copies of the Bell state |Φ + (n "ebits"), rearranges the qubits as described above, and sends Bob's set of n qubits to him over a noisy channel. It is not necessary that Alice must prepare Bell pairs locally and then transmit half the qubits to Bob. Indeed, the protocol is applicable as long as Alice and Bob share some initial (noisy) Bell pairs. Then, Alice measures the stabilizers of a quantum stabilizer code defined by S = ε i E(a i , b i ); i = 1, . . . , r on her qubits, with ε i = ±1. Let her measurement results be (−1) mi , m i ∈ {0, 1}. This projects her qubits onto the codespace fixed by the stabilizers S = (−1) mi ε i E(a i , b i ); i = 1, . . . , r . Alice applies some suitable Pauli "correction" to bring her qubits back to the code subspace Q(S) (rather than Q(S )), if that is the code she desires to use. She classically communicates the chosen stabilizers, S, the measurements {m i } r i=1 , and the Pauli correction to Bob. Although we use the term "correction", there is really no error on Alice's qubits. Instead, the terminology is used to indicate that Alice brings the qubits to her desired code space. Furthermore, even if there is some error on Alice's qubits, one can map it to an equivalent error on Bob's qubits using the Bell state matrix identity.
Note that the authors of Ref. [18] do not explicitly mention that the Pauli correction needs to be communicated, but it could be necessary in situations where Alice's and Bob's decoders are not identical or have some randomness embedded in them. For the code, though any appropriate definition of logical Pauli generators works with the protocol, we employ Algorithm 3 to obtain generators that are "compatible" with our way of analyzing the protocol (using the stabilizer formalism). This phenomenon will become more clear after the [[5, 1, 3]] code example in this section. While the algorithm simulates measurements on GHZ states to define logical Paulis, an equivalent algorithm can be constructed that only simulates Bell measurements.
Remark 1. In this protocol, whenever the syndrome of Alice is non-trivial, i.e., at least one m i equals 1, she can either perform a Pauli correction or just define her code to be Q(S ) and not perform any correction. If the protocol is defined so that she always does the latter, as depicted in Fig. 5 where there is no 'Recovery' block on Alice's qubits, then Bob can adjust his processing accordingly based on the syndrome information from Alice.
Without loss of generality, we can assume that Alice sends Bob's qubits to him only after performing her measurements and any Pauli correction. So, the channel applies a Pauli error only after Bob's qubits got projected according to . . , r . Now, Bob measures the stabilizers ε i E(a i , b i ) and applies corrections on his qubits using his syndromes as well as Alice's syndromes (and the Bell matrix identity, which in particular involves the transpose). This projects his qubits to the same codespace as Alice. Finally, Alice and Bob locally apply the inverse of the encoding unitary for their code, U Enc (S) † . If Bob's correction was successful, this converts the k logical Bell pairs into k physical Bell pairs that are, on average, of higher quality than the n noisy Bell pairs initially shared between them. This protocol is shown in Figure 5 and summarized in Algorithm 1.
While the steps of the protocol are clear, it is worth considering why the logical qubits of Alice and Bob must be k copies of the Bell pair, assuming all errors were corrected successfully. To get some intuition, let us quickly consider the example of the 3-qubit bit-flip code defined by S = ZZI, IZZ . According to (57), the projector onto Q(S) is Π S = (I8+ZZI) 2 (I8+IZZ) 2 . The encoding unitary, as described in Appendix I-B, is U Enc = CNOT 1→2 CNOT 1→3 . Since Z T = Z, Alice's measurements will project Bob's qubits onto the same code subspace as her's. For convenience, assume that Alice obtains the trivial syndrome  . . , r = n − k} on her n qubits, obtains syndrome, 5 (b) sends the remaining n qubits to Bob over a noisy quantum channel, 6 (c) sends the stabilizers and syndrome (which together define Q(S)) to Bob over a noiseless classical channel. . . , r = n − k} on his n qubits, 10 (b) combines the syndrome information from Alice as well as his measurements and interprets using Section I-D, 11 (c) performs necessary Pauli corrections on his qubits to bring them to the code space of Q(S). 12 13 If the channel error was correctable, pairs of logical qubits of Alice's and Bob's codes form k Bell states 14 // If channel error was NOT correctable, some pair of logical qubits form a Bell state with an unknown Pauli error 15 Alice and Bob respectively apply the inverse of the encoding unitary for their code on their n qubits 16 // The encoding unitary is determined by the logical Pauli operators for Q(S) obtained from Algorithm 3 (+1, +1) and that the channel does not introduce any error. Then, according to (64), the resulting (unnormalized) state after Alice's measurements is (Π S ⊗ Π S ) Φ + 3 . Consider the action of (I 8 + ZZI) on a computational basis state |x , x = [x 1 , x 2 , x 3 ] ∈ F 3 2 : Hence, after the action of (Π S ⊗ Π S ) and inversion of the encoding unitary by Alice and Bob, we obtain Thus, the output is a single Bell pair and ancillary qubits on Alice and Bob. In Appendix II, we show this phenomenon for arbitrary CSS codes by generalizing the state vector approach used above.

A. Bell Pair Distillation using the 5-Qubit Code
In the remainder of this section, with the [ [5,1,3]] code [13], [24] as an example, we use the stabilizer formalism to show that the above phenomenon is true for any stabilizer code. Recall that this code is defined by S = XZZXI, IXZZX, XIXZZ, ZXIXZ .
As described in Appendix I-B, the corresponding binary stabilizer matrix is given by Initially, Alice starts with 5 copies of the standard Bell state |Φ + ⊗5 , and marks one qubit of each copy as Bob's. She does not yet send Bob's qubits to him. The stabilizer group for this joint state of 5 "ebits" (or "EPR pairs") is where e i ∈ F 5 2 is the standard basis vector with a 1 in position i and zeros elsewhere, 0 ∈ F 5 2 is the all-zeros vector, and the Xand Zcomponents in the E(a, b) notation have been split into Alice's qubits and Bob's qubits. Observe that this is a maximal stabilizer group on 10 qubits and hence, there are no non-trivial logical operators associated with this group, i.e., the normalizer of E 5 in P 10 is itself.
It will be convenient to adopt a tabular format for these generators, where the first column of each row gives the sign of the generator, the next two columns give the X-components of Alice and Bob in that generator, the subsequent two columns give the Z-components of Alice and Bob in that generator, and the last column gives the Pauli representation of that generator for clarity. Hence, the above generators are written as follows.
Given this "initialization", let us track these 10 stabilizers through each step of the protocol, as shown in Table I.
(1) Alice measures the first stabilizer generator X A1 Z A2 Z A3 X A4 , and assume that the measurement result is ε 1 ∈ {±1}. We apply the stabilizer formalism for measurements from Section I-C to update E 5 . Since there are several elements of E 5 that anticommute with this generator, we choose to remove 1 ) and replace all other  1,3]] code. Any '0' that is not part of a string represents 00000, and e i ∈ F 5 2 is the standard basis vector with a 1 in the i-th position and zeros elsewhere. Code stabilizers are typeset in boldface. An additional left arrow indicates which row is being replaced with a code stabilizer, i.e., the first row that anticommutes with the stabilizer. Other updated rows are highlighted in gray. Classical communications: A → B.
Step Sign anticommuting elements by their product with Z A4 Z B4 . Let this updated group in Step (1) of Table I be denoted as E 5 . For visual clarity, code stabilizer rows are boldfaced and binary vectors are written out in full. Now, we observe that if Bob measures the same generator on his qubits, i.e., X B1 Z B2 Z B3 X B4 , then it is trivial because it commutes with all elements in E (1) 5 and hence is already contained in E (1) 5 . This is a manifestation of the Bell state matrix identity discussed in Section I-D. Indeed, Bob's generator can be obtained by multiplying Table I. (2) Alice measures the second stabilizer generator X A2 Z A3 Z A4 X A5 , and assume that the measurement result is ε 2 ∈ {±1}.
Then, the new joint stabilizer group, E 5 , is given in Step (2) of Table I. This stabilizer generator anticommutes with the third row of the top block and the second and fifth rows of the bottom block. We have replaced Z A5 Z B5 (fifth row of the bottom block) with this generator and multiplied the other anticommuting elements with Z A5 Z B5 . It can be verified that the second stabilizer generator of Bob is already in E (3) Alice measures the third stabilizer generator X A1 X A3 Z A4 Z A5 , and assume that the measurement result is ε 3 ∈ {±1}.
Then, the new joint stabilizer group, E 5 , is given in Step (3) of Table I. Once again, it can be verified that the third stabilizer generator of Bob is already in E (4) Alice measures the final stabilizer generator Z A1 X A2 X A4 Z A5 , and assume that the measurement result is ε 4 ∈ {±1}.
Then, the new joint stabilizer group, E 5 , is given in Step (4) of Table I. As before, it can be verified that the final stabilizer generator of Bob is already in E (4) 5 . This completes all measurements of Alice, and she now sends Bob's qubits over the channel. To understand the working of the protocol in the ideal scenario, assume that no errors occur.
Since we know that all stabilizer generators of Bob are in E (4) 5 , we conveniently perform the following replacements: Recollect that for the [[5, 1, 3]] code, the logical Pauli operators are . If we used Algorithm 3, we would obtain the same Z and X = −Y 1 Z 3 Z 4 . Then, by grouping Alice's code stabilizers and Bob's code stabilizers, the group E (4) 5 can be rewritten as Using some manipulations, we recognize that the two operators on the second line in E Thus, E can be interpreted as having 8 stabilizer generators (Alice and Bob combined) and a pair of logical X A X B and logical Z A Z B operators, which implies that the pair of logical qubits shared between Alice and Bob forms a Bell pair. This can be converted into a physical Bell pair by performing the inverse of the encoding unitary on both Alice's and Bob's qubits locally. Note that this encoding unitary must be compatible with the above definition of the logical Paulis for the [[5, 1, 3]] code, i.e., when the physical X and Z on the input (logical) qubit to the encoder is conjugated by the chosen encoding unitary, the result must be the above logical Paulis X and Z, respectively, potentially multiplied by some stabilizer element.
Remark 2. In this example, we have assumed that Bob's qubits do not suffer any error, so that we can clearly show the existence of the correct logical Bell stabilizers. If, however, the channel introduced an error, then Alice and Bob can jointly deduce the error by measuring the signs of all generators of E (1) 5 and applying the necessary Pauli correction. Since there are no non-trivial logical Pauli operators, any syndrome-matched correction can differ from the true error only by a stabilizer, so any error is correctable by the joint action of Alice and Bob. But, since we prohibit non-local measurements between Alice and Bob, our error correction capability is limited to that of the code (on Bob's side). If the channel introduces a correctable Pauli error for the chosen code and Bob's decoder, then the protocol will output k perfect Bell pairs. However, if the Pauli error is miscorrected by Bob's decoder, then there will be a logical error on the code, and hence at least one of the k output Bell pairs will suffer from an unknown Pauli error.
We can arrive at the above conclusion without knowing the specific logical operators for the code. After Alice measures all her stabilizer generators, we know that Bob's stabilizer generators will also be present in the group, simply based on the Bell state matrix identity from Section I-D. For this example, the transpose in that identity did not make a difference, but for other codes this can only introduce an additional minus sign since Y T = −Y . For an [[n, k, d]] code, we now have a 2n-qubit stabilizer group E (n−k) n where 2(n − k) generating elements are Alice's and Bob's stabilizer generators. We are left with 2n − 2(n − k) = 2k elements in the generators, each of which must jointly involve Alice's and Bob's qubits. These commute with each other and with the 2(n − k) stabilizer generators of Alice and Bob, and are independent, so we can rename them as the logical X Aj X Bj and logical Z Aj Z Bj for j = 1, 2, . . . , k. Thus, by definition, the k pairs of logical qubits form k logical Bell pairs. Alice and Bob can produce physical Bell pairs by simultaneously inverting the (same) encoding unitary for the code locally. This is the key idea behind the working of the Bell pair distillation protocol employed by Wilde et al. in [18].

V. DISTILLATION OF GREENBERGER-HORNE-ZEILINGER (GHZ) STATES
In this section, we extend the above Bell pair distillation protocol to distill GHZ states, GHZ = (|00···0 +|11···1 ) . For clarity, we will specifically discuss the standard case of = 3, but the results and analysis extend to larger as well. Let n GHZ states be shared between Alice, Bob, and Charlie. We rearrange all the qubits to keep Alice's, Bob's and Charlie's qubits together respectively. Hence, this joint state can be expressed as Since the GHZ state has stabilizers Thus, we have identified the GHZ version of the basic properties of Bell states that was needed in the Bell pair distillation protocol. However, the critical part of the Wilde et al. protocol was the transpose trick that formed the Bell matrix identity in Appendix I-D. When applied to stabilizer codes, this implied that each stabilizer generator εE(a, b) of Alice is transformed into the generator εE(a, b) T = ε(−1) ab T E(a, b) (using Lemma 9(a)) for Bob. Naturally, we need to determine the equivalent phenomenon for GHZ states before we can proceed to constructing a distillation protocol.

A. GHZ State Matrix Identity
In the following lemma, we generalize the Bell state matrix identity in Appendix I-D to the GHZ state.
M xy |x y| ∈ C 2 n ×2 n be any matrix acting on Alice's qubits. Then, Proof: Similar to the Bell case, we calculate This completes the proof and establishes the identity. The above property generalizes naturally to larger -qubit GHZ states, GHZ = (|00···0 +|11···1 ) Lemma 4. Let M = x,y∈F n 2 M xy |x y| ∈ C 2 n ×2 n be any matrix acting on qubits 'A'. Then, As our next result, we prove some properties of the GHZ-map defined in the above lemma.
Lemma 5. The GHZ-map M ∈ C 2 n ×2 n → M ∈ C 2 2n ×2 2n in Lemma 3 is an algebra homomorphism [27]: (c) Projector-preserving: If M is a projector, then M is also a projector.
Proof: We prove these properties via the definition of the mapping. (a) Since M xy = x| M |y = α x| A |y + β x| B |y = αA xy + βB xy , the property follows.
(c) This simply follows from the multiplicative property via the special case A = B = M . This completes the proof and establishes the said properties of the GHZ-map.
We are interested in performing stabilizer measurements at Alice and deducing the effect on Bob's and Charlie's qubits. The above properties greatly simplify the analysis, given that the code projector for a stabilizer code (57) is a product of sums. Due to the multiplicativity of the GHZ-map M → M , we only have to analyze the case where Alice's code has a single stabilizer generator εE(a, b), i.e., her code projector is simply M = I N +εE(a,b) 2 , where N = 2 n . Now, using linearity, we just need to determine I N and E(a, b). Then, due to Lemma 9(a), we have ) . Theorem 6. Given n copies of the GHZ state shared between Alice, Bob and Charlie, measuring E(a, b) on Alice's n qubits and obtaining the result ε ∈ {±1} is equivalent to measuring the following with results +1 on the qubits of Bob and Charlie: where Z Bi (resp. Z Ci ) refers to Z on i-th qubit of Bob (resp. Charlie), and e i has a 1 in the i-th position and zeros elsewhere.
Proof: Using the discussion before the statement of the theorem, we will calculate I N and E(a, b) to establish the result. Recollect that |0 0| n=1 = I+Z 2 and hence |0 0| ⊗n = 1 2 n v∈F n 2 E(0, v). Then, using Lemma 9, we have where e i ∈ F n 2 is the standard basis vector with 1 in the i-th position and zeros elsewhere. Note that E([0, 0], [e i , e i ]) BC = Z Bi ⊗ Z Ci is the GHZ stabilizer I A Z B Z C on the i-th triple of qubits between A, B and C (15). Next, we proceed to calculate E(a, b) using a similar approach.
Thus, when Alice's measurement applies the projector M = I N +εE(a,b) 2 , Bob's and Charlie's qubits experience the projector Since the second term, I N , only corresponds to already existing stabilizers Z Bi ⊗ Z Ci for n copies of the GHZ state, the only new measurement corresponds to the Pauli operator ε(−1) ab T E(a, b) ⊗ E(a, 0). More interestingly, if we consider a Y -measurement for Alice, where a = b = 1, then E(a, b) T ⊗ E(a, 0) = Y T ⊗ X = −Y ⊗ X. Thus, assuming the measurement result is +1, the new BC stabilizers are −Y ⊗ X, Z ⊗ Z . It can be verified that the post-measurement state for this case will be (|0 +ı|1 ) , which is stabilized by the above stabilizer.
Naturally, this insight can be generalized to larger GHZ states as well.
Theorem 7. Given n copies of the -qubit GHZ state with subsystems A 1 , A 2 , . . . , A , measuring E(a, b) on the n qubits of subsystem A 1 and obtaining the result ε ∈ {±1} is equivalent to measuring the following with results +1 on the qubits of the remaining ( − 1) subsystems: x * y denotes the element-wise product of two vectors, Z At,i refers to Z on i-th qubit of subsystem A t , and e i is the standard basis vector with a 1 in the i-th position and zeros elsewhere.
, splitting the entries of b into ( − 1) disjoint groups.
As we desired, the above result shows how a Pauli measurement on one subsystem, A 1 , of (multiple copies of) the GHZ state affects the remaining subsystems. All the GHZ stabilizers involving subsystems A 2 , A 3 , . . . , A are retained. Hence, the post-measurement state is "GHZ-like" on these ( − 1) subsystems but with an additional globally entangling stabilizer. This is akin to the globally entangling all-X stabilizer for the standard GHZ state, but it depends on the Pauli operator being measured on A 1 . Note that, since the Pauli measurement randomly projects onto a subspace, the induced stabilizers given by the theorem do not uniquely determine the post-measurement state on the ( − 1) subsystems. The degrees of freedom for the state will be quantified shortly in a more general setting. One might argue that this theorem can be obtained by directly applying the stabilizer formalism to S GHZ . However, some thought clarifies that arriving at the conclusions rigorously takes at least an equal amount of effort.
In the context of measuring a set of (n − k) stabilizer generators of a code (on qubits A 1 ), the above result confirms that this induces a joint stabilizer code on the remaining ( − 1) subsystems. There are n( − 1) qubits on these subsystems and each code stabilizer generator contributes a stabilizer generator for this induced code. Besides, as stated in the theorem, there are n( − 2) GHZ stabilizers on all pairs of adjacent subsystems, {A j A j+1 ; j = 2, . . . , − 1}, independent of the code stabilizers being measured. Hence, the induced code has (n − k) + n( − 2) stabilizer generators, which means it is an [[n( − 1), k]] code and the post-measurement state has k logical degrees of freedom. The minimum distance of the induced code will depend on the minimum distance of the A 1 -code as well as the new GHZ stabilizers and the choice of {b i }.

B. Protocol I
We now have all the tools to investigate a natural stabilizer code based GHZ distillation protocol that attempts to generalize the Bell pair distillation protocol discussed in Section IV. The block diagram of this protocol was shown earlier in Fig. 2 and the protocol is summarized as an algorithm in Algorithm 2. Let us consider the 3-qubit code with stabilizers S = Y Y I, IY Y to understand the subtleties in the steps of the protocol. First, similar to the Bell pair scenario, we have the following stabilizer group for 3 copies of the GHZ state: Next, like the example for the Bell pair distillation protocol, we can evolve these stabilizers through the proposed steps of the protocol to understand its working. In Appendix III, we use such a tabular approach to elucidate the steps of this protocol. This serves as an instructive example to understand how the GHZ property influences the construction of a purification protocol for GHZ states. In particular, since the property implies that the Z-component of any non-purely-X-type stabilizer is lost in the induced code on qubits 'B' and 'C', we discuss how one can perform diagonal Clifford operations to ensure that all three subsystems are imposed the same code. The placement of these operations is critical and we detail its effects by simulating the protocol performance for the 5-qubit code.
In this protocol, Alice starts by preparing n GHZ states and measuring the n-qubit stabilizers of her code on qubits 'A'. Then, using Theorem 6, Bob proceeds by measuring the 2n-qubit stabilizers of the code induced on qubits 'B' and 'C' by Alice's choice of code on qubits 'A'. Subsequently, he also measures the same n-qubit stabilizers as Alice but on qubits 'B', so that there is a code induced just on qubits 'C' and Charlie can use that code to correct errors from the channel. If we imagine the three parties being on a linear network topology, then this protocol seems reasonable since each party retains his/her qubits and passes on all remaining qubits to the next hop in the chain. However, there is an asymmetry in the operations since Bob needs to perform two rounds of measurements and one involves twice the number of qubits. Furthermore, the protocol is (potentially) not scalable to larger number of parties with varied network topologies.

C. Distillation-Inspired Algorithm to Generate Logical Pauli Operators
While constructing and analyzing the protocol using the tabular approach, we realized that the evolution of the table under stabilizer measurements automatically reveals the logical Pauli operators of the code in an explicit manner in certain rows. Indeed, each stabilizer measurement replaces one row and alters several others that anticommute with it using the rules of the stabilizer formalism for measurements (Section I-C). After all stabilizers are measured on qubits 'A', one realizes that the non-replaced (but altered) rows in the top section of the table, i.e., the ZZI-type rows, are of the form Z where Z i denotes the logical Z operator on the i-th logical qubit of the code. Therefore, one can easily read off these logical operators (up to some subtleties that can be taken care of). A similar approach is applied to the bottom section of the table, i.e., the XXX-type rows, to obtain the logical X operators of the code. The details of the algorithm are discussed in Appendix IV-B and the algorithm itself is summarized in Algorithm 3.
⊗ E(0, e j ) C ; i = 1, 2, . . . , r = n − k , j = 1, 2, . . . , n} on the received 2n qubits, 10 (b) performs necessary Pauli corrections on all qubits to bring them to the code space of the joint BC code, 11 (c) measures the stabilizer generators {E(a i , b i ) ; i = 1, 2, . . . , r = n − k} on the n qubits of subsystem B and obtains syndrome {ε B i }; for purely Z-type stabilizers E(0, b i ) the sign is ε A i , so we set ε B i := +1 for them, 12 (d) sends the stabilizers, syndrome and logical Pauli operators to Charlie over a perfect classical channel, 13 (e) performs appropriate (see Appendix IV-A) local diagonal Clifford on qubits C, 14 (f) sends qubits C to Charlie over a noisy Pauli channel. 15 16 Charlie 17 (a) uses Q(S), Theorem 6, and Bob's syndrome to determine the signs ε

his stabilizers, and then measures the generators {ε
. . , r = n − k} on his n qubits, 18 (b) performs the necessary Pauli corrections on all qubits to bring them to the code space of his code. 19 20 // If the channel error was correctable, triples of logical qubits of Alice's, Bob's and Charlie's codes form k GHZ states 21 // If channel error was NOT correctable, some triple of logical qubits form a GHZ state with an unknown Pauli error 22 Alice, Bob, and Charlie respectively apply the inverse of the encoding unitary for their code on their n qubits 23 // The encoding unitary is determined by the logical Pauli operators obtained from Algorithm 3

D. Output Fidelity of GHZ Distillation Protocol
During the protocol, if error correction at Bob and/or Charlie miscorrects and introduces a logical error, then the final effect is a change in the signs of some of the logical GHZ stabilizers. This in turn means that after the decoding step, some of the k triples will be the standard GHZ state corrupted by an unknown Pauli operation. Hence, the output of the protocol is correct

Algorithm 3: Algorithm to generate logical Paulis of a stabilizer code through GHZ measurements (see Appendix IV-B)
Input : An [[n, k, d]] stabilizer code defined by its stabilizer generators {ε i E(a i , b i ) ; i = 1, 2, . . . , r = n − k} Output: The logical X generators {ν j E(c j , d j ) ; j = 1, . . . , k} and logical Z generators {E(0, f j ) ; j = 1, . . . , k} 1 Initialization: Form a r × (2n + 1) binary parity-check matrix H for the code, whose rows are [a i , b i , ε i ]. Preprocess the matrix so that its first 2n columns take the form H 1:2n = 0 H Z H1 H2 , where H Z is a r Z × n matrix of full rank, and H 1 is a r X × n matrix of full rank (r X + r Z = r = n − k). The rows of H Z provide the generators for all purely Z-type stabilizers of the code. While performing row operations on H 1:2n , care must be taken to adhere to Pauli multiplication arithmetic (Lemma 9(b)). 2 3 Simulate the creation of n copies of the GHZ state as follows. Create a 2n × (6n + 1) GHZ stabilizer matrix S GHZ whose first n rows take the form [0, 0, 0, e i , e i , 0, +1] and the second n rows take the form [e i , e i , e i , 0, 0, 0, +1], where i = 1, 2, . . . , n. This matrix is almost the same as Step (0) in Table II, but we have omitted the middle section. Simulate the measurement of H (p) , the p-th row of H, on subsystem A of the GHZ states, using Section I-C: 7 Replace the first anticommuting row of S GHZ with H (p) and multiply subsequent anticommuting rows by H (p) , using Lemma 9(b) 8 end 9 10 for q in the set of non-replaced rows of S GHZ with (row) index at most n do 11 if S Form a k × n matrix F whose rows are f j

34
Define the new f j 's as the rows of T −1 F 35 else 36 Retain the definitions of logical Z generators E(0, f j ) and logical X generators ν j E(c j , d j ) 37 end 38 39 return {Z j = E(0, f j ) , X j = ν j E(c j , d j ) ; j = 1, 2, . . . , k} with probability (1 − P f ), and produces at least one Pauli corrupted GHZ state with probability P f , using the notation in Fig. 6. To make this precise, denote by |GHZ 0 , |GHZ 1 , . . . , |GHZ 7 the eight possible variants of the GHZ state under Pauli operations, i.e., each variant has the stabilizer group α 1 ZZI, α 2 IZZ, α 3 XXX with α 1 , α 2 , α 3 ∈ {±1}. Then, assuming all variants are equally likely conditioned on a failure event, the density matrix representing the output of the protocol is where Similar to the case of triorthogonal codes in magic state distillation [33], it is likely useful to consider the reduced density matrix for one of the k output triples, and relate its fidelity (with respect to |GHZ GHZ|) to properties of the code and decoder. In Ref. [33], the authors adopted exactly such a strategy for distillation of T -states, under a purely Z-error model and relying on post-selection where non-trivial syndromes are discarded. In recent work [34], it has been shown that performing error correction rather than just detection (and post-selection) leads to better performance of triorthogonal codes. For our scenario of GHZ distillation, it is an interesting problem to construct codes and decoders for this protocol where we can relate the output fidelity to code properties and arrive at analytical scaling arguments with increasing code size. This would be useful for comparing with fundamental limits of entanglement distillation [17] and assessing the optimality of this protocol.

E. Protocol II
. . , r = n − k} on the n qubits 'C'; for purely Z-type stabilizers E(0, b i ) the sign is ε A i , 16 (b) uses the syndrome to run a decoder that estimates a Pauli error 17 (c) applies the estimate to qubits 'C' as the recovery operation 18 19 // If the channel error was correctable, triples of logical qubits of Alice's, Bob's and Charlie's codes form k GHZ states 20 // If channel error was NOT correctable, some triple of logical qubits form a GHZ state with an unknown Pauli error 21 Alice, Bob, and Charlie respectively apply the inverse of the encoding unitary for their code on their n qubits 22 // The encoding unitary is determined by the logical Pauli operators obtained from Algorithm 3 To address the drawbacks of Protocol I, the protocol can be modified so that Alice starts by measuring qubits 'A' and qubits 'B' separately. Though this does not circumvent the issue of performing twice the number of measurements at one of the nodes, this avoids the need of 2n-qubit measurements. Since the GHZ property implies the inducement of a 2n-qubit code on qubits 'B' and 'C', it appears that this extra round of n-qubit measurements on qubits 'B' is inevitable. So, even now, when Alice measures 0) C } jointly on qubits 'B' and 'C'. But, when she measures the same stabilizers {ε i E(a i , b i )} on qubits 'B', one can multiply with the corresponding 2n-qubit stabilizer to see that the joint stabilizers can be broken into purely 'B' and purely 'C' stabilizers. Therefore, once Alice performs the two rounds of measurements, she can send qubits 'B' to Bob and qubits 'C' to Charlie, along with the necessary classical information. As individual codes have been induced separately on qubits 'B' and qubits 'C', Bob and Charlie can still perform local n-qubit measurements to fix errors during qubit transmission. Finally, this scheme suits other network topologies such as when Alice is connected to both Bob and Charlie but those parties do not have a direct connection between them. Note that for CSS codes, for any stabilizer generator ε i E(a i , b i ), whenever a i = 0 we have b i = 0. Hence, the induced code from Theorem 6 is automatically CSS and we do not need any diagonal Clifford operation mentioned earlier in Protocol I. Since Protocol II relies on the same intuitions from Theorem 6, we do not elaborate further. We also note that there can be further variations based on other practical considerations.
This simplified protocol was shown earlier in Fig. 3. We simulated the protocol by following a tabular approach, as in Appendix III for Protocol I, using a state-of-the-art family of lifted product QLDPC codes with asymptotic rate 0.118 and an efficient iterative decoder based on the min-sum algorithm (MSA) with normalization factor 0.8. The results were shown in Fig. 4, where we can see that the threshold under depolarizing noise is about 10.7%. Comparing the results to Fig. 1, it is apparent that the threshold matches that of the underlying logical error rate of the code on this channel (i.e., no distillation but standard quantum error correction simulation). This is important because it shows that even when both Bob and Charlie run decoders to correct errors on their respective qubits, the overall threshold is unchanged from the single channel setting. On the other hand, the comparison also shows that the protocol failure rate is significantly worse for each channel parameter compared to Fig. 1. This could be the effect of requiring both decoders to succeed, but it is a cause for concern when we extend the protocol to GHZ states with > 3. Indeed, we do not want the protocol failure rates to progressively get worse, albeit with the same threshold. Therefore, we will study this phenomenon more carefully in future work and identify the best way to scale this protocol for larger . For completeness, we summarize the protocol for arbitrary .

F. Protocol II for Arbitrary
Initially, A 1 generates n ideal copies of the -qubit GHZ state, names the qubits of each copy A 1 through A , chooses some [[n, k]] code Q(S) defined by a stabilizer S, and measures the (n − k) generators of S on qubits A 1 . Then, A 1 applies Theorem 6 to determine the induced code Q S ( −1) on the remaining subsystems. Let us consider = 4 for simplicity. For tracking the protocol, we initially create a table whose rows are the binary representations of the generators of S ⊗n GHZ . Group the n Z Ai Z Bi I Ci I Di generators in the first part of the table, the n I Ai Z Bi Z Ci I Di generators in the second part, the n I Ai I Bi Z Ci Z Di in the third part, and finally the n X Ai X Bi X Ci X Di in the fourth part. If there is a purely Z-type generator, E(0, b) A , for S, then it will commute with the first three parts and only affect the last part based on the stabilizer formalism. Moreover, by an appropriate linear combination of the rows of the first part, one can produce the element E(0, b) A ⊗ E(0, b) B , which when multiplied by the new code stabilizer produces the stabilizer E(0, b) B on purely subsystem 'B'. By a similar trick in the second part and subsequently in the third part, one can produce single-subsystem stabilizers E(0, b) C and E(0, b) D as well. Hence, it suffices to only consider non-purely-Z-type stabilizers E(a, b) A , a = 0.
Such stabilizers transform into the multiple-subsystem stabilizers described by Theorem 6. Now, qubits of 'B', 'C', and 'D' need to be transmitted over a noisy channel to the respective nodes, based on the network topology. For those nodes to be able to correct errors, a code needs to be imposed purely on each subsystem before transmission of the respective qubits. Let (node) A be connected to (node) B. Then, based on the choice of b 1 , b 2 , b 3 in Theorem 6, A measures code stabilizers E(a, b 1 ) B on qubits 'B'. With some thought, one sees that these stabilizers only affect the second part of the table. Now, since ±E(a, b 1 ) B ⊗ E(a, b 2 ) C ⊗ E(a, b 3 ) D is already a stabilizer, by multiplying with E(a, b 1 ) B we obtain a code on 'B' and a residual code jointly on 'C' and 'D'. The qubits of 'B' can be transmitted to node B (along with necessary classical sign information of stabilizers), which can perform error correction.
If A is not connected to C and D, then A has to send those qubits to B. Thus, it appears that A has to perform stabilizer measurements as above not only on 'B' but on 'C' and 'D' as well. However, this can be relegated to subsequent nodes to reduce the burden on A. Let A also send qubits 'C' and 'D' to node B along with qubits 'B'. There is some joint Pauli error on 'B', 'C', and 'D', and the error correction of B only fixes the error part on 'B'. If B measures code stabilizers on 'C', then the preexisting Pauli error can be transformed into an effective Pauli error after the code was imposed on 'C'. This enables node C to correct this error as well as any error encountered while B sends qubits 'C'. A similar statement holds for D as well. Thus, the protocol can be stated as follows: for every edge connected to a node, the node performs stabilizer measurements on the respective subsystem to impose a code on the qubits of the recipient on that edge. The correctness of the protocol relies on carefully tracking signs of stabilizers based on such measurements at each node. Once all qubits are distributed, each node uses the logical Paulis of their respective codes to determine and invert the encoding unitary. This converts the k logical GHZ states into k perfect physical GHZ states, provided all error corrections were successful. The average output density matrix and average output fidelity still take the form discussed in Section V-D.

VI. CONCLUSION AND FUTURE WORK
In this work, we began by describing the Bell pair distillation protocol introduced in Ref. [18], and used the stabilizer formalism to understand its working. We identified that the Bell state matrix identity (Appendix I-D) plays a critical role in that protocol. As our first result, we proved the equivalent matrix identity for GHZ states, where we introduced the GHZ-map and showed that it is an algebra homomorphism. Using the GHZ-map, we proved our main result (Theorem 6) that describes the effect of Alice's stabilizer measurements (on qubits 'A') on qubits 'B' and qubits 'C'. Then, we constructed a natural GHZ distillation protocol whose steps were guided by the aforementioned main result. We demonstrated that the placement of a certain local Clifford on qubits 'C' in the protocol has an immense effect on the performance of the protocol. We described the relation between the probability of failure of the protocol and the output fidelity of the GHZ states. As part of our protocol, we also developed a new algorithm to generate logical Pauli operators for an arbitrary stabilizer code. To circumvent some drawbacks of the protocol, we described an alternate protocol and produced performance results using state-of-the-art QLDPC codes and an efficient iterative decoder. Finally, we discussed the scalability of the protocol for larger GHZ states involving more than 3 parties and arbitrary network topologies.
In future work, we plan to study the scaling of the logical error rate with the increase in number of parties. Since a key motivation for this work was distributed quantum computing (DQC), we will investigate a complete architecture for a distributed implementation of the recently proposed optimal families of QLDPC codes. As part of the architecture, we envisage that the QLDPC-based GHZ purification scheme proposed in this paper will play a critical role in supplying logical GHZ states encoded in the same QLDPC codes that are used for DQC. We will study the implications for fault-tolerance of such an architecture.
is the symplectic inner product between [a, b] and [c, d] in F 2n 2 , as indicated by the subscript 's'. Hence, the map γ : (P n , ·) → (F 2n 2 , ⊕) defined by E(a, b) → [a, b] is a homomorphism with kernel {ı κ I N , κ ∈ Z 4 }. For details about extending the definition of E(a, b) to Z 4 -valued arguments, see [36].

B. Stabilizer Codes and Encoding Unitaries
A stabilizer group S is a commutative subgroup of P n that does not contain −I N . If the group has r ≤ n independent generators ε i E(a i , b i ), where ε i ∈ {±1}, then S = ε i E(a i , b i ); i = 1, . . . , r has size |S| = 2 r . Since the generators are Hermitian and unitary, they have eigenvalues ±1. Recollect that commuting matrices can be diagonalized simultaneously. The stabilizer code defined by S is the common +1 eigenspace of all generators, i.e., it is the 2 k -dimensional subspace, k = n − r, fixed by all elements of S: (56) Using the homomorphism γ, we can write a r × (2n + 1) generator matrix G S for the stabilizer group: Since S must be a commutative group, the symplectic inner product between any pair of rows must be zero. Hence, the subspace of binary mappings of all elements of S, denoted γ(S), is given by the rowspace of G S .
A CSS (Calderbank-Shor-Steane) code is a special type of stabilizer code for which there exists a set of generators where either b i = 0 or a i = 0 in each generator, i.e., the generators are purely X-type and purely Z-type operators. Clearly, for such a code, G S has a block diagonal form where we can express the X-type (resp. Z-type) operators as the rowspace of a matrix [H X , 0] (resp. [0, H Z ]), and 0 represents the all-zeros matrix (of appropriate size). In this case, the commutativity condition for stabilizers is equivalent to the condition H X H T Z = 0. Therefore, H X and H Z can be thought of as generating two classical linear codes C X and C Z .
The projector onto the +1 eigenspace of a Pauli matrix E(a, b) is I N +E(a,b)

2
. Therefore, since Q(S) is the simultaneous +1 eigenspace of r commuting matrices ε i E(a i , b i ), the projector onto the code subspace Q(S) is While the stabilizer group S defines the code space, an encoding unitary U Enc (S) fully specifies the mapping from logical k-qubit states to physical n-qubit code states in Q(S). The n input qubits to U Enc (S) can be split into k logical qubits, whose joint state is arbitrary, and r = n − k ancillary qubits, each of which is initialized in some specific state such as |0 . If ancillas are initialized in the |0 state, then the stabilizer group for these n input qubits is generated by {Z i ; i = k + 1, k + 2, . . . , n}, since Z |0 = |0 . If we conjugate each of these r generators by the encoding unitary U Enc (S), then we will obtain r generators U Enc (S) Z i U Enc (S) † of S. Similarly, if we conjugate the X i and Z i operations on the k logical qubits -which can be used to express arbitrary operations on them since Pauli operators form a basis -by U Enc (S), then we will obtain the generators of logical X and Z operators compatible with the chosen U Enc (S). Therefore, an alternative method to specify U Enc (S) is to specify the generators of S as well as the generators of logical X and Z operators. Since we are requiring U Enc (S) to map Paulis to Paulis, it is always Clifford [23]. Note that U Enc (S) is still not unique since we are not specifying how X i on the ancillas must be mapped, but we do not care about these additional mappings.
There are at least two algorithms provided in the literature for generating the logical Pauli operators of stabilizer codes. One is by Gottesman [21], [37], where the idea is to construct the normalizer of the stabilizer group inside the Pauli group, and then perform suitable row operations on the generators of the normalizer. The other is by Wilde [38], where he performs a symplectic Gram-Schmidt orthogonalization procedure to arrive at the generators of logical X and logical Z operators. In this work, as part of our GHZ distillation protocol, we provide a new algorithm to generate logical X and Z operators for any stabilizer code (see Algorithm 3). The output of the algorithm is compatible with the way logical Paulis must be defined for our analysis of the protocol. Additionally, the logical Z operators from our algorithm are always guaranteed to be purely Z-type operators for any stabilizer code. If the code is CSS, then the logical X operators are always purely X-type.

C. Stabilizer Formalism
When a stabilizer group on n qubits has n independent generators, Q(S) is a 1-dimensional subspace that corresponds to a unique quantum state |ψ(S) (up to an irrelevant global phase), commonly referred to as a stabilizer state [21]. The actions of unitary operations and measurements on |ψ(S) can be tracked by updating these n generators accordingly [23], [32]. For any element g of S, and an arbitrary unitary operation U on |ψ(S) , we observe that U |ψ(S) = U · g · |ψ(S) = (U gU † ) · U |ψ(S) , so the stabilizer element g has evolved into the element g = U gU † after the action of U . Of course, only if U is a Clifford operation we have that g is also a Pauli matrix. Thus, in this case the evolution of the state can be tracked efficiently by simply transforming G S (and the associated signs) through binary operations (see "CHP" algorithm [32]).
The stabilizer formalism also provides a method to systematically update the stabilizers under Pauli measurements of the state |ψ(S) . Assume that we have n generators for the stabilizer group, namely ε i E(a i , b i ), i = 1, . . . , n, and that we are measuring the Pauli operator µE(u, v) to obtain the measurement (−1) m , m ∈ {0, 1}. Then, we have the following cases.
1) If [u, v], [a i , b i ] s = 0 for all i, then either E(u, v) or −E(u, v) already belongs to S, so there is nothing to update.

D. Bell State Matrix Identity
Let n standard Bell pairs be shared between Alice and Bob. We rearrange the 2n qubits to keep Alice's qubits together and Bob's qubits together. We can write the joint state as Let M = x,y∈F n 2 M xy |x y| ∈ C 2 n ×2 n be any matrix acting on Alice's qubits. Then, it has been observed that [13], [18] ( Example 3. As in Example 2, consider the standard Bell pair and measure Y ⊗ I. If the measurement result is +1, then the projector P Y = I2+Y 2 gets applied to the first qubit. Then, according to the above identity, this is equivalent to applying on the second qubit. This exactly agrees with the post-measurement state (|0 +ı|1 ) If Alice measures the generators of a stabilizer group S = ε i E(a i , b i ); i = 1, . . . , r , and obtains results (−1) mi , m i ∈ {0, 1}, then M = Π S is the projector onto the subspace Q(S ) of the stabilizer code defined by S = (−1) mi ε i E(a i , b i ); i = 1, . . . , r . According to the above identity, this is equivalent to projecting Bob's qubits onto the stabilizer code defined by where we have applied Lemma 9(a). Note that, in such cases where M is a projector, we can write so that the action of Alice can be interpreted as both Alice and Bob projecting their own qubits simultaneously.

APPENDIX II LOGICAL BELL PAIRS FOR ARBITRARY CSS CODES
In this appendix, we show that when n raw Bell pairs are projected onto the subspace of a CSS code through stabilizer measurements, the induced logical state is that of k Bell pairs. We take a meet-in-the-middle approach where we first consider k Bell pairs and show how their encoded state looks like, and then we project n Bell pairs to prove that the resulting state is the same as the aforesaid encoded state.
Let C 1 , C 2 be two binary linear codes such that C 2 ⊂ C 1 . For the [[n, k]] CSS code defined by these codes, C 2 produces the X-stabilizers, C 1 produces the logical X operators, C ⊥ 1 produces the Z-stabilizers, and C ⊥ 2 produces the logical Z operators. Let G C1/C2 denote a generator matrix for the quotient group C 1 /C 2 that represents the "pure" logical X operators that do not have any X-stabilizer component. In other words, the rows of G C1/C2 give the generators of logical X operators for the CSS code. Let U Enc denote an encoding unitary for the code. Then, the encoded state of k Bell pairs is [37], [39] For the other direction, we start with |Φ + n AB and then apply the projector Π CSS for the code on Alice's qubits. By the Bell state matrix identity (Section I-D), this means that we are effectively simultaneously applying Π T CSS = Π CSS on Bob's qubits as well. Here, the transpose has no effect because the stabilizer generators for CSS codes are purely X-type or purely Z-type, and only such operators appear in the expression for the code projector (57). Let G 2 and G ⊥ 1 represent generator matrices for the codes C 2 and C ⊥ 1 , respectively. Then, we have For any z ∈ F n 2 , since E(0, v) |z = (−1) zv T |z , we have (I N + E(0, v)) |z = 2 |z if zv T = 0 and (I N + E(0, v)) |z = 0 otherwise. This implies that Π Z |z = |z or 0 depending on whether z ∈ C 1 or not, respectively. Similarly, it is easy to check that Π X |z = 1 |C2| y∈C2 |z ⊕ y . Putting these together, we observe that This state must be normalized by the square root of the probability that we get the all +1 syndrome, which corresponds to the subspace of the considered CSS code. It can be checked that all syndromes are equally likely, so the probability is 1/2 n−k . Dividing (73) by 1/ √ 2 n−k , we arrive at exactly the same state in (66). This establishes that when CSS stabilizer measurements are performed on n Bell pairs, the resulting code state corresponds to k logical Bell pairs.

APPENDIX III PROTOCOL I WITH THE 3-QUBIT CODE
The steps of the protocol, with this particular [[3, 1, 1]] code as an example, are shown in Table II. Again, we use the stabilizer formalism for measurements from Section I-C. We will explain each step below and discuss the potential subtleties that can arise. It could be useful to imagine the three parties as being three nodes A -B -C on a linear network chain.
(0) Alice locally prepares 3 copies of the perfect GHZ state and groups her qubits together for further processing. She keeps aside the grouped qubits of Bob's and Charlie's but does not send those to them yet. She also writes down the parity check matrix for the 9 qubits, based on only GHZ stabilizers, along with signs, as shown in Step (0) of Table II. (1) Alice measures the stabilizer Y A1 Y A2 I A3 = E([(e 1 + e 2 ) A , 0 B , 0 C ], [(e 1 + e 2 ) A , 0 B , 0 C ]) and the group G 3 gets updated as shown in Step (1) of Table II, assuming that the measurement result is ε A 1 ∈ {±1}. Based on the stabilizer formalism (Section I-C), the measured stabilizer replaces the first row (as indicated by the left arrow) and the second row is multiplied with the previous first row. For visual clarity, code stabilizer rows are boldfaced and binary vectors are written out in full. Furthermore, as per Theorem 6, this measurement of E(e 1 + e 2 , e 1 + e 2 ) by Alice should imply that ε A 1 E(e 1 + e 2 , e 1 + e 2 ) T B ⊗ E(e 1 + e 2 , 0) C = ε A 1 E(e 1 + e 2 , e 1 + e 2 ) B ⊗ E(e 1 + e 2 , 0) C automatically belongs to the (new) stabilizer group. Indeed, this element can be produced by multiplying the elements  1,1]] code defined by S = Y Y I, IY Y . Each '0' below represents 000, and e i ∈ F 3 2 is the standard basis vector with a 1 in the i-th position and zeros elsewhere. Code stabilizers are typeset in boldface. An additional left arrow indicates which row is being replaced with a code stabilizer, i.e., the first row that anticommutes with the stabilizer. Other updated rows are highlighted in gray. Classical communications: A → B, B → C.
Step Sign (±1) Alice creates 3 copies of the ideal GHZ state and writes down their 9 stabilizer generators; she replaces XXX with −Y Y X for convenience Alice measures code stabilizers on her qubits, updates first 3 rows, and rewrites −Y Y X entries to introduce joint BC stabilizers by Theorem 6 000 000 011 000 000 Alice applies P † = 1 0 0 −ı to all 3 qubits of C system to change X's to Y 's, then sends B's and C's qubits to Bob over noisy Pauli channel 000 000 011 000 000 110 110 000 110 110 Bob measures the 5 joint BC stabilizers, corrects errors (signs), measures B stabilizers, and sends C's qubits to Charlie over noisy Pauli channel 000 000 011 000 000 011 000 000 011 000 1 000 000 110 000 000 110 2 000 000 011 000 000 011 Charlie measures his code stabilizers, corrects errors, inverts the encoding unitary; Alice and Bob invert their encoding after sending qubits Lemma 9(b). This is exactly how the seventh row gets updated.
(2) Alice measures the second stabilizer I A1 Y A2 Y A3 = E([(e 2 +e 3 ) A , 0 B , 0 C ], [(e 2 +e 3 ) A , 0 B , 0 C ]) and the group gets updated as shown in Step (2) of Table II. The procedure is very similar to that in Step (1).
Since Alice has measured all her stabilizer generators, and the stabilizer formalism preserves the commutativity of the elements in the group, the third row in the first block of 3 rows must necessarily commute with Alice's stabilizers. Thus, the Alice component of the third row must form a logical operator for Alice's code, and we define it to be the logical Z operator, i.e., Z A = ZZZ = E(0, e 1 + e 2 + e 3 ). We will see shortly that Bob's qubits get the same code (possibly with sign changes for the stabilizers), so this third row can be written as the logical GHZ stabilizer Z A Z B I C .
This phenomenon also generalizes to any [[n, k, d]] stabilizer code, with some caveats when the code has some purely Z-type stabilizers, and we determine the logical Z operators either after Alice's set of measurements or apriori using some linear algebraic arguments (see Appendix IV-B for details). Note that it is convenient to choose the logical Z operators such that they respect the GHZ structure of our analysis, e.g., Z = IIY will not be compatible here.
(3) If we consider the parity-check matrix after Step (2), we see that rows 4 through 8 are the stabilizers promised by Theorem 6 that act only on B and C systems. However, due to the same result, the C parts of rows 7 and 8 only have X-s instead of Y -s. So, to change them back to Y -s, Alice applies the inverse of the Phase (i.e., √ Z) gate to all 3 qubits of the C system. She specifically applies the inverse, rather than √ Z itself, to get rid of the −1 sign for the last row.
For a general stabilizer code, the appropriate diagonal Clifford must be chosen as discussed in Appendix IV-A. This operation converts the BC stabilizers ε which ensures that Charlie gets the same code (up to signs of stabilizers) as Alice and Bob. Since the Clifford is guaranteed to be diagonal, it leaves purely Z-type stabilizers unchanged. Later, in Section III-A, we show that it is better for Bob to perform this Clifford on Charlie's qubits, rather than Alice.
Though we have used the −Y Y X GHZ stabilizer here for convenience, for a general code we can simply continue to use XXX. After Alice has measured her stabilizer generators, this last block of 3 rows could have changed but they still commute with the generators. Since the middle block never gets affected by Alice's measurements, we can guarantee using Theorem 6 that two of the last 3 rows must be the joint BC stabilizers induced by Alice's two generators. Hence, the remaining row's Alice component must form a logical operator for Alice's code, and will be distinct from the previously defined logical Z operator. We define this to be the logical X operator, i.e., X A = IIY = Y 3 = E(e 3 , e 3 ). As we will see shortly, both Bob and Charlie get the same code, so this last row of the third block can be written as the logical GHZ stabilizer X A X B X C . The generalization to arbitrary stabilizer codes is discussed in Appendix IV-B.
Then, she sends Bob both his qubits as well as Charlie's qubits over a noisy Pauli channel, which introduces the signs η, ν i , µ i ∈ {±1}, i = 1, 2, 3. She also classically communicates the code stabilizers, her syndromes {ε A 1 , ε A 2 }, and the logical Z and X operators to him.
(4) Now, based on Alice's classical communication, Bob applies Theorem 6 to obtain the stabilizer generators B ⊗ E(e 2 + e 3 , e 2 + e 3 ) C for the induced joint code on his as well as Charlie's qubits. As per Theorem 6, he also includes the IZZ-type GHZ stabilizers, i.e., i , e C i ]) ; i = 1, 2, 3}, as stabilizer generators for the [ [6,1]] joint code on BC systems. He measures these 5 stabilizers to deduce and correct the error introduced by the channel on the 6 qubits sent by Alice. Assuming perfect error correction, the signs will be back to the ones in Alice's final parity-check matrix.
When Bob sends Charlie's qubits to him, the channel might introduce errors on those 3 qubits. To deduce and correct these errors, there must have been a code induced on Charlie's qubits even before the transmission. Hence, after correcting errors on the 6 qubits of BC systems, Bob measures the same stabilizers as Alice's code but on his qubits. This produces rows 4 and 5 in Step (4) of Table II, and row 6 gets updated as per the stabilizer formalism. When looking at rows 7 and 8 of Step (3), it is evident that one can correspondingly multiply them with these new rows 4 and 5 to produce the same code just on Charlie's qubits.
This phenomenon extends to general stabilizer codes as well, where the joint stabilizers ε XZZXI, IXZZX, XIXZZ, ZXIXZ . The result is shown as the green curve marked 'purely X code for C' in Fig. 6. When compared to the standard QEC performance of this code on the depolarizing channel, we see that the exponent is worse. This is because, by the arguments above, all non-purely Z-type stabilizers of Q(S) get converted into a purely X-type stabilizer for Charlie's code. To mitigate this, we make Alice perform a local diagonal Clifford operation on qubits C to transform E(a i , 0) C into E(a i , b i ) C later. The performance of this is shown as the solid blue curve marked 'same code for all', which is equally worse. This time the reason is that the BC code Q(S ) has stabilizers of the form E(a i , b i ) B ⊗ E(a i , b i ) C , which means that the code cannot distinguish the i-th qubit of B and the i-th qubit of C. Finally, we make Bob perform the aforesaid diagonal Clifford on qubits C, and this produces the solid dark red curve marked 'local Clifford by B'. Clearly, the exponent now approaches the "fundamental limit" set by standard QEC on the depolarizing channel.

APPENDIX IV IMPLEMENTATION DETAILS OF PROTOCOL I
In this appendix we provide additional details regarding the implementation of Protocol I, which is also relevant for Protocol II. Specifically, we discuss how to identify the appropriate diagonal Clifford to be applied on qubits C, and explain the reasoning behind Algorithm 3 that generates logical Pauli operators for any stabilizer code. The discussion related to Algorithm 3 will also clarify some aspects of the distillation protocol.

A. Diagonal Clifford on qubits C
A diagonal Clifford unitary on n qubits can be described using an n × n binary symmetric matrix R as [40], [41] (74) Its action on a Pauli matrix E(a, b) is given by [36], [41] U where b * aR is the entrywise product of the two binary vectors. It is well-known that any diagonal Clifford operator can be formed using the phase gate, P , and the controlled-Z gate, CZ, defined as Set n = 2. Then, the R matrices for P 1 = (P ⊗ I), P 2 = (I ⊗ P ), and CZ are respectively R P1 = 1 0 0 0 , R P2 = 0 0 0 1 , and R CZ = 0 1 1 0 .
This generalizes naturally to more than 2 qubits. The diagonal entries of R describe which qubits get acted upon by P , and the pairwise off-diagonal entries describe which pairs of qubits get acted upon by CZ [41]. In our protocol, from Theorem 6 we observed that Alice's stabilizers of the form E(a i , b i ), with a i = 0, become the stabilizer E(a i , 0) for Charlie. This also means that any logical X operator of the form E(c j , d j ) would have transformed into E(c j , 0) (e.g., see last row of Step (2) in Table II and compare it to the last row of Step (3)). Therefore, the purpose of the diagonal Clifford on qubits C is to convert the stabilizers E(a i , 0) back into E(a i , b i ) and the logical X operators E(c j , 0) back into E(c j , d j ). Given the above insight into diagonal Clifford operators, we want to find a binary symmetric matrix R such that U R E(a i , 0)U † R = E(a i , b i ) for all i = 1, 2, . . . , r X (using notation in Algorithm 3) and U R E(c j , 0)U † R = E(c j , d j ) for all j = 1, 2, . . . , k. Thus, we need a feasible solution R for {a i R = b i , i = 1, 2, . . . , r X } and {c j R = d j , j = 1, 2, . . . , k}.
We solve this system of linear equations on a binary symmetric matrix as follows. First, let A be the matrix whose rows are {a i } and {c j }, and let B be the matrix whose rows are {b i } and {d j }. Then, we have the system AR = B. We recall the vectorization property of matrices, which implies that vec(QU V ) = (V T ⊗ Q)vec(U ).
Here, vectorization of a matrix is the operation of reading the matrix entries columnwise, top to bottom, and forming a vector (e.g., this is done through the command U(:) in MATLAB). Setting Q = A, U = R, V = I, we get (I ⊗ A)vec(R) = vec(B), which is a standard linear algebra problem for the unknown vector vec(R). However, we desire a binary symmetric matrix R. We impose this constraint as (I − W )vec(R) = 0, where 0 denotes the all-zeros vector of length n 2 , and W is the permutation matrix which transforms vec(Q) into vec(Q T ) for any matrix Q. In summary, we obtain the desired R (or equivalently the diagonal Since R is symmetric, it has n(n + 1)/2 degrees of freedom, which accounts for the second constraint. The matrix A has r X + k < n rows and the Kronecker product with I results in n(r X + k) constraints on n(n + 1)/2 variables. It remains to be shown if there is always a feasible solution for any valid A and B. Note that [A, B] represents a matrix whose rows are stabilizers and logical X operators. This means any pair of rows must be orthogonal with respect to the symplectic inner product, which implies that AB T + BA T = 0. Thus, a given A and B is valid if and only if AB T is symmetric.

B. Logical Paulis from GHZ Measurements
The procedure in Algorithm 3 to determine logical X and Z generators of a stabilizer code is inspired by the stabilizer measurements on Bell or GHZ states, viewed through the lens of the stabilizer formalism for measurements (Section I-C). Though the algorithm could have been constructed just using measurements on Bell states, we preferred GHZ states because there can be an additional non-trivial sign for the logical X operators due to an odd number of subsystems in the GHZ state. Of course, logical operators obtained using GHZ states will also apply to the Bell protocol since a negative sign on an even number of subsystems (A and B in Bell states) leads to an overall positive sign for X A X B and Z A Z B .
We have a code Q(S) defined by its stabilizer group S = ε i E(a i , b i ) ; i = 1, 2, . . . , r = n − k . Define the r × (2n + 1) stabilizer (or parity-check) matrix H whose i-th row is [a i , b i , ε i ]. First, we bring the first 2n columns of the stabilizer (or parity-check) matrix of the code to the following standard form: Here, the r Z rows of H Z form all generators for the purely Z-type stabilizers of the code. The bottom part of the matrix is such that the r X × n matrix H 1 has full rank (r X + r Z = r = n − k). While performing row operations on the initial parity-check matrix H , one has to account for the Pauli multiplication rule in Lemma 9(b), and not simply perform binary sums of (the first 2n columns of the) rows, i.e., the last column of H must be updated to reflect changes in signs.
Next, we simulate the creation of n GHZ states by creating a 2n × (6n + 1) GHZ stabilizer matrix S GHZ , whose first n rows are [0, 0, 0, e i , e i , 0, +1] and the second n rows are [e i , e i , e i , 0, 0, 0, +1]. This matrix is the same as Step (0) of Table II, but we have omitted the middle section since the measurements on subsystem A trivially commute with entries I A Z B Z C of this section. Now, we use the stabilizer formalism for measurements (Section I-C) to simulate measurements of the rows of H on subsystem A of S GHZ . Clearly, the stabilizers from [0, H Z ] commute with the first n rows, so these will only replace r Z rows in the bottom half of S GHZ . The stabilizers from [H 1 , H 2 ] will necessarily anticommute with at least one of the first n rows of S GHZ , and these r X rows get replaced. This can be established by counting the dimension of purely Z-type operators with which each row of [H 1 , H 2 ] can commute, one after the other. Crucially, the stabilizer formalism guarantees that all rows of the evolved S GHZ remain linearly independent and always commute.
The (n − r X ) non-replaced rows within the first n rows can be divided into two types. Before we simulate any stabilizer measurements, the first n rows have standard basis vectors e i for the Z-parts of A and B. These can be rewritten such that we have r Z rows of the form [0, 0, 0, z, z, 0, +1], where z corresponds to rows of H Z , all of which are linearly independent by assumption of the standard form. Since these correspond to code stabilizers (on A as well as B), the measurement of rows of [H 1 , H 2 ] will not replace these. After the measurements, when r Z rows in the bottom half have been replaced by [0, 0, 0, z, 0, 0, ε z ], we can multiply with the corresponding rows of the top half, i.e., [0, 0, 0, z, z, 0, +1], to produce purely Z-type stabilizers on subsystem B, which later define Bob's code. These Z-operators on B in the top half form the first type of r Z rows. The remaining (n − r X ) − r Z = k rows form the second type, and they have to form logical Z Aj Z Bj , for j = 1, 2, . . . , k, since they commute with all code stabilizers and the columns of subsystem C remain zero. Thus, the Z-component of subsystem A of these k rows produce the logical Z generators of the code, and they always have sign +1.
A similar argument applies to the bottom half of the evolved S GHZ matrix. The stabilizer measurements from H Z replace r Z rows out of the n rows. The remaining (n − r Z ) rows can again be divided into two types. The first type of rows give operators that can be rewritten as the BC stabilizers guaranteed by Theorem 6. Specifically, these can be identified by the fact that their A-parts will be linearly dependent on the A-parts of the other rows of the evolved S GHZ matrix. Indeed, this is how one can cancel the A-parts of these r X rows to produce the r X BC stabilizers corresponding to [H 1 , H 2 ]. The remaining (n − r Z ) − r X = k rows of the bottom half form the second type, and they have to form logical X Aj X Bj , for j = 1, 2, . . . , k, since they commute with all code stabilizers and are linearly independent from all other rows. The A-parts of these rows are used to define the logical X operators of the code. The primes on these logical operators indicate that they might not exactly pair up with the corresponding logical Z operators defined earlier. This is because they are only guaranteed to be logical operators independent of the logical Z operators, but not to be the appropriate pairs {X j } of the previously determined {Z j }.
Once these pseudo logical X operators are determined, we can easily find the necessary pairs for the logical Z operators. Let the logical Z operators be E(0, f i ), i = 1, 2, . . . , k, and let these pseudo logical X operators be ν j E(c j , d j ), j = 1, 2, . . . , k. Therefore, if F is the matrix whose rows are f i and [C, D] is the matrix whose rows are [c j , d j ], then we need If T = I k , then we can simply pre-multiply the equation by T −1 (mod 2) to achieve the desired result. In this case, we redefine the logical Z operators to be given by the rows of T −1 [0, F ]. This completes the reasoning behind Algorithm 3.