Approximating Decoherence Processes for the Design and Simulation of Quantum Error Correction Codes on Classical Computers

Quantum information is prone to suffer from errors caused by the so-called decoherence, which describes the loss in coherence of quantum states associated to their interactions with the surrounding environment. This decoherence phenomenon is present in every quantum information task, be it transmission, processing or even storage of quantum information. Consequently, the protection of quantum information via quantum error correction codes (QECC) is of paramount importance to construct fully operational quantum computers. Understanding environmental decoherence processes and the way they are modeled is fundamental in order to construct effective error correction methods capable of protecting quantum information. Moreover, quantum channel models that are efficiently implementable and manageable on classical computers are required in order to design and simulate such error correction schemes. In this article, we present a survey of decoherence models, reviewing the manner in which these models can be approximated into quantum Pauli channel models, which can be efficiently implemented on classical computers. We also explain how certain families of quantum error correction codes can be entirely simulated in the classical domain, without the explicit need of a quantum computer. A quantum error correction code for the approximated channel is also a correctable code for the original channel, and its performance can be obtained by Monte Carlo simulations on a classical computer.


I. INTRODUCTION
Since Richard Feynman's original and ground-breaking proposal in [1] of constructing computers that follow the laws of quantum mechanics to simulate physical systems that obey said laws, the scientific community has gone to extraordinary lengths in order to build an operational quantum computer. Following the introduction of these novel ideas about constructing quantum machines, research has shown that application of quantum physics theory is not only useful for simulation of complex quantum mechanical systems such as The associate editor coordinating the review of this manuscript and approving it for publication was Abdullah Iliyasu . macromolecules for drug discovery [2]- [5], but also to efficiently solve tasks which are computationally unmanageable in a reasonable amount of time for classical computers. The most prominent of these tasks are the factorization of prime numbers and the discrete logarithm problem [6], Byzantine agreement [7], or searching an unstructured database or an unordered list [8], [9]. Furthermore, quantum computing is a powerful asset for secure communications. For instance, Quantum Key Distribution (QKD) protocols enable two parties to create a shared random secret key only known to them. The key remains secure given that these protocols allow the aforementioned parties to detect if a malicious entity is trying to gain knowledge of it, which enables them VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ to adapt their message exchange so that the eavesdropper extracts no information. The best-known cryptographic QKD protocols are the BB84 protocol proposed by Bennett and Brassard in [10] and the E91 protocol by Ekert presented in [11]. In light of the astonishing potential of quantum frameworks, the construction of devices capable of exploiting the benefits offered by the paradigm of quantum computing represents a step forward in the advance of technology, regardless of these machines being in the form of fully operational quantum computers [12], quantum processors as accelerators of classical computers [13], or specific devices that perform QKD [14]. Unfortunately, quantum information is vulnerable to errors that arise and corrupt quantum states while they are being processed, which oftentimes leads to incorrect algorithm outcomes. The frailty of quantum information is caused by the phenomenon known as quantum decoherence, which refers to the destruction of the superposition of quantum states from the interaction that they have with the environment [15]. This quantum noise arises during every task related to the quantum computing paradigm: information storage, processing or communication. Hence, it is necessary to invoke Quantum Error Correction Codes (QECC) in order to have qubits with sufficiently long coherence times 1 for practical applications. Quantum information is so sensitive to decoherence that many think that quantum computation is unfeasible without the aid of quantum error correction tools.
The earliest formulation of QECCs appeared in 1995 [15] when Shor proposed a 9-qubit code, which was later aptly named after him. This code is capable of correcting errors of weight one and it does not saturate the Quantum Hamming Bound (QHB) [16], which implies that the same task could be achieved with codes of shorter length. Nevertheless, it stands as the first proposal of an error correcting code for the quantum computing paradigm. Since then, research efforts have focused on deriving QECC schemes that approach the quantum capacity limits [17] at a reasonable complexity cost. Encoding and decoding gate depths play an important role in the paradigm of correcting quantum errors due to the fact that quantum gates introduce additional errors. Additionally, the runtimes of the decoding algorithms should not exceed the coherence times of the qubits that are being processed, else these states would suffer from new decoherence effects while being corrected. A major breakthrough in the development of QECCs came in Gottesman's Ph.D. thesis [19], where he proposed the theory of Quantum Stabilizer Codes (QSC), which is a very useful framework that facilitates the construction of QECC families from classical binary and quaternary error correction codes. This formalism has led to the development of several promising QECC families such as Quantum Reed-Muller codes [20], Quantum Low Density Parity Check (QLDPC) codes [21], [22], Quantum Low Density Generator Matrix (QLDGM) codes [23], [24], Quantum Convolutional Codes (QCC) [25], Quantum Turbo Codes (QTC) [26]- [29], and Quantum Topological Codes [30], [31].
QECC design requires the assumption of an error model in the form of a quantum channel that accurately represents the decoherence processes that affect quantum information. Based on these error models, appropriate strategies to combat the effects of decoherence can be derived. In order for the designed QECCs to be applicable in realistic quantum devices, the quantum channels should capture the essential characteristics of the physical processes that make qubits lose their coherence. The depolarizing model is a widespread quantum error model used to evaluate the error correcting abilities of QECC families [15], [16], [19]- [22], [24]- [31]. This decoherence model is especially useful due to the fact that it makes the system fulfill the Gottesman-Knill theorem, 2 and thus it can be efficiently simulated on a classical computer [32], [33]. At the time of writing, the availability of quantum computers for researchers is limited and the accessible machines operate on a reduced number of qubits. Therefore, classical resources remain an invaluable tool for the design of advanced QECCs that will be used beyond the Noisy Intermediate-Scale Quantum (NISQ) era. The NISQ era [34] is a term coined by Preskill that makes reference to the time when quantum computers will be able to perform tasks that classical computers are incapable of, but will still be too small (in qubit number) to provide fault-tolerant implementations of quantum algorithms.
The motivation for this article is to provide a clear description of how the physical processes that lead to the corruption of quantum information can be approximated so that they are tractable using classical resources. To do so, we survey the techniques existing in the literature to approximate general quantum channels (resulting from decoherence phenomena) by the widely used Pauli channels. Consequently, we review the decoherence processes that affect quantum information compromising the integrity of quantum algorithms, and their mathematical description, in the form of quantum channels. We then study the use of the technique called twirling [35]- [40] in order to obtain approximations of those quantum channels, resulting in the asymmetric Pauli channel and the depolarizing channel. We also justify that any error correction methods designed for the twirled approximated channels will also be valid for the more realistic original channels. Furthermore, we provide an extensive discussion regarding the Pauli channel and its symmetric instance, commonly known as the depolarizing channel, as these approximations are widely used for QECC constructions with classical resources. Additionally, even though the memoryless version is generally considered, we discuss the emergence and impact of memory effects 2 The depolarizing channel can be implemented by applying random Pauli gates to the encoded quantum states. The Gottesman-Knill theorem states that Pauli gates are quantum computations that can be efficiently simulated on classical computers. on those channels [41], [42]. Finally, we survey the way in which twirled approximations are implemented on a classical computer, and discuss how they can be employed to simulate the performance of QECCs, as done in [24], [28], [29].
The main contribution of this article is two-fold. On the one hand, we clearly explain how it is possible to construct and evaluate certain families of quantum error correction codes without the explicit need of a quantum computer, since the error models and quantum codes can be entirely simulated in the classical domain. On the other hand, we describe how the performance of the aforementioned quantum error correcting schemes under general decoherence models can be assessed on classical computers. The key is that decoherence models can be approximated as Pauli channels. which can be implemented with classical resources, and that the performance of quantum error correcting codes over the approximated channels is equivalent to that obtained for the original decoherence model.
The rest of the paper is organized as follows: section II describes the typical decoherence models that are used for QECC design; section III approximates such models into channels that can be implemented efficiently on classical computers via twirling techniques; section IV discusses the use of memory in the twirled models in order to use more realistic approximations when correlations between different quantum information units exist; section V describes the basic theory of stabilizer coding and the way the so called Pauli-to-binary isomorphism is used to efficiently simulate such QECCs under the corruption of the presented twirled channels. Finally, section VI concludes the paper.

II. DECOHERENCE AND THE DAMPING CHANNELS
Environmental decoherence is described as the undesired interaction 3 of a qubit with the environment, resulting in the perturbation of the coherent superposition of the basis states of such a quantum information unit. Decoherence mechanisms arise from several physical processes. Thus, these physical interactions are described mathematically based on different representations making use of quantum information theory.
The physical processes that generate this quantum noise depend on the technology that is used to construct the qubits of the quantum computer. Some of the most promising technologies currently in use are transmon superconducting qubits [43], trapped ion qubits [44] and quantum dot qubits [45] among others. Decoherence also takes place during the transmission of quantum states through fiber optics or via laser beams. However, despite the diversity of the quantum mechanical processes that generate errors in the quantum information that is being processed, they can all be mathematically represented by the so-called combined amplitude and phase damping channel. 3 In quantum mechanical terms this interaction is more specifically called entanglement.

A. AMPLITUDE DAMPING CHANNEL
Quantum systems suffer from energy dissipation, which is the name given to the effects that quantum mechanical states suffer as a result of spontaneous energy losses. This can happen when an atom in an excited state emits a photon and returns to its ground state, when a spin system at high temperature approaches equilibrium with its environment, when a photon in an interferometer is subjected to scattering and attenuation; or when a photon is absorbed during its transmission through an optical fiber. An example of why this energy loss occurs can be appreciated by observing the energy profile of a transmon qubit, which is shown in figure 1. Presently, transmon superconducting qubits are an extremely popular technology for building experimental quantum computers (IBM or Rigetti computing use this technology for their prototypes). In the energy level chart presented for said qubits, it is easy to see that the states defining the qubits (or qudits if more than two energy levels are taken) are given by the possible quantized energy states of such a superconducting transmon. Therefore, the |0 state corresponds to the ground state of the system, and the |1 state is related to an excited state. Due to the fact that an isolated system tends to drop into its lowest energy state, decoherence for this type of qubit takes place when energy is lost to the environment as a consequence of the excited state approaching the ground state. |0 state corresponds to the ground state of the system while the |1 state is an excited state of such a transmon. Source: [46].
Every energy dissipation process has its own unique features, but in the quantum information framework their behaviour can be well described mathematically using the so-called amplitude damping channel. Consider the two-level transmon qubit described before, whose ground state is denoted by |0 T and its excited state by |1 T . Additionally, consider the mathematical abstraction that the environment is also a two-level system with a vacuum state |0 E and an excited state |1 E , that is always initialized in the vacuum state. Then the amplitude damping channel describes the evolution of the composite system as [47] where γ refers to the damping probability or the probability that the system loses energy to the environment when it is in its excited state. The expressions in (1) mathematically describe, in a general and unified way, what was previously discussed regarding the processes of energy dissipation that occur for the various methods used for the construction of qubits. Qubits are superpositions of states, so the effect that the amplitude damping channel induces in a general qubit 4 |ψ = α|0 + β|1 is described by the following transformation: Following the rationale of (1), what (2) is describing is the fact that state |ψ decoheres to state with probability 1−γβ 2 leaving the environment unchanged; or it decoheres to |0 with probability γβ 2 releasing a photon (or energy quanta) to the environment, and so leaving it in state |1 E .
It is important to remark that in the discussion above |ψ has been assumed to be an isolated qubit. However, this is not true in general, since in a real system more than one qubit is present, and they will most assuredly interact with each other. This interaction is reflected in the form of entanglement, which implies that the decoherence processes affect the set of qubits in a combined fashion, and not in the isolated manner described before. Regardless, the assumption that each of the qubits interacts with the environment independently is still reasonable [47], simplifying the description of the amplitude damping channel to that of (2).
Up to this point, the description has been based on pure quantum states. However, a quantum system can be represented in a more general fashion as a statistical ensemble of pure states. In the latter representation, quantum states should be described by means of the so-called density matrices. It is important to mention that the formulation of quantum mechanics for density matrices is equivalent to the one of pure states [32], but that each is more convenient for specific tasks of quantum information theory. In the context of the density matrix formulation, a quantum channel N is a completely-positive, trace-preserving (CPTP) linear map between the spaces of operators. Such CPTP maps can be 4 Note that the subscript T is no longer used. The rationale now refers to every qubit technology, and not just transmon superconducting qubits. Its prior use was intended as a particular instance for explanatory purposes.
represented by the so-called Choi-Kraus decomposition or operator-sum representation, which means that the input density matrix ρ is mapped by quantum channel N as where the E k matrices are operators on the state space of the principal system and receive the name of Kraus operators or error operators of the channel. As quantum channels are trace-preserving operators, Kraus operators must fulfill k E † k E k = I, with I being the identity matrix. The minimum number of Kraus operators is called the Kraus rank of the quantum channel N . Channels with Kraus rank 1 are called pure.
For this description of quantum channels, the amplitude damping channel N AD describing the decoherence processes that cause energy loss in qubits acts on an arbitrary quantum state with density matrix ρ as where the error operators of the channel are It is important to look at how the Kraus operators are related to the description that has been given about the amplitude damping channel in (2). If we apply each of the error operators to the pure state |ψ , we will see that each of the possible outcomes that were present in (2) are obtained. Applying E 0 , which, up to normalization, is the state that occurs with probability |E 0 |ψ | 2 = (1 − γβ 2 ) in (2) if the environment is left unchanged. If E 1 is applied, which, up to normalization, is the state that occurs with probability |E 1 |ψ | 2 = γβ 2 in (2) if an energy quanta is lost to the environment. In consequence, we can see the equivalence relationship between the Kraus operators of the CPTP map of the amplitude damping channel and the description of this channel that was given for pure states.
To conclude with the description of the amplitude damping channel as a mathematical model that represents the decoherence effects on qubits related to energy dissipation processes, we should relate the damping probability with some physical parameter that represents the constructed realistic qubits. The temporal evolution of the damping probability is given by [32] where T 1 is the qubit relaxation time. The definition of T 1 is system specific, i.e it depends on the technology that is being used, but it can be measured for the employed qubits [48], and so it is related to the physical system via the mathematical model of the amplitude damping channel. Table 1 shows the relaxation time of the quantum machines

B. DEPHASING CHANNEL
Another set of physical mechanisms that affect quantum information are those encompassed by the terms dephasing or phase damping. These processes are uniquely quantum mechanical and describe the loss of quantum information without loss of energy. This kind of environmental decoherence can arise when a photon scatters randomly as it travels through a waveguide or when electronic states are perturbated due to the action of stray electrical charges. What happens during these events is that the system evolves for an amount of time which is not known with precision, and so partial information about its quantum phase is lost. A simple model to mathematically describe this particular set of processes can be constructed using phase kicks on the qubits [32]. Phase kicks are rotations applied to the qubit |ψ with a random angle θ. Assuming that the random variable θ follows a Gaussian distribution with mean 0 and variance 2λ, then the dephasing channel or phase damping channel N PD acts on the qubit with density matrix ρ as where the Kraus operators of the channel are and λ is interpreted as the scattering probability of a photon without loss of energy. In a similar fashion to what was done for the amplitude damping channel, we apply the error operators of the channel to an arbitrary pure state |ψ = α|0 + β|1 . The application of E 0 results in with probability |E 0 |ψ | 2 = (1 − λβ 2 ). The effect of E 1 produces which is the state |1 , up to normalization, occurring with probability λβ 2 .
To conclude with the description of the phase damping channel, we have to relate the scattering probability λ to physical parameters that can be measured for the qubits of each of the possible existing technologies, as was previously done for the damping probability. The temporal evolution of the scattering probability that defines the dephasing channel [35] is where T 1 is the qubit relaxation time as defined for the amplitude damping channel and T 2 is called the qubit dephasing time. T 2 is often referred to as the 'spin-spin' relaxation time, and it is considered one of the subtlest processes in the quantum information paradigm, which has led, unsurprisingly, to extensive research on the topic [32].

C. COMBINED AMPLITUDE AND PHASE DAMPING CHANNEL
If the previously described amplitude and phase damping channels are combined, a fairly complete mathematical model of qubit decoherence is obtained, as several physical processes that corrupt quantum information are encompassed by this abstraction. We will designate this amalgamation of channels as the combined amplitude and phase damping channel N APD . The effect of this combined channel is captured by the next operator-sum representation with Kraus rank 3 [35]: with the following error operators: where I, X, Y, Z are the Pauli matrices. We use P n to refer to the set of n-fold tensor products of the Pauli matrices P n = {I, X, Y, Z} ⊗n [55]. The damping, γ , and the scattering, λ, probabilities, as well as their corresponding time evolutions, are defined as in earlier sections of this work.
If we now consider the combination of (8), (13), (14) and (15), the time evolution of an arbitrary quantum state with density matrix ρ can be obtained as where ρ ij corresponds to the element of the matrix ρ in row i and column j. This last expression shows that the qubits are likely to decohere if the operation time (either information transmission, processing or storage) t is of the same order of magnitude as either the relaxation time (t ≈ T 1 ) or the dephasing time (t ≈ T 2 ). This means that the reliability time of a qubit is defined by T 1 and T 2 , and so the run times of the algorithms that will be executed on such devices will be limited to the aforementioned reliability time. 5

III. TWIRL APPROXIMATIONS OF QUANTUM CHANNELS
The combined amplitude and phase damping channel introduced in the previous section, embodies an accurate mathematical depiction of several decoherence processes that corrupt the quantum information that is processed in a quantum device, transmitted between machines or stored in a quantum memory. Thus, the application of this model to create QECCs that will guarantee the stability of quantum information in the presence of decoherence seems reasonable. Unfortunately, the fact that the dimensions of the Hilbert spaces of N -qubit composite systems scale with a factor of 2 N makes matters increasingly more complex, since this implies that the classical simulation of quantum algorithms or error correction mechanisms may ultimately turn into an intractable problem. For example, a simple surface code of distance d = 5 needs 25 physical qubits, which is a system that has a Hilbert space with 33 millions of dimensions [35]. As a result, the error dynamics of QECC schemes cannot be efficiently modeled on classical computers 6 by using (16), meaning that an approximation that is manageable employing conventional methods must be found.
In this section we present how approximated channels that can be efficiently simulated on classical computers are obtained by using a quantum information theory technique called twirling [35]- [40]. Additionally, we justify why QECCs designed for the twirled channels will also be valid for the original channel, and why those approximated channels are indeed tractable as classical problems. This means that the quantum computing community has access to the necessary design tools working with algorithms that rely only on the classical machines that are available in this day and age.

A. THE GOTTESMANN-KNILL THEOREM AND THE PAULI CHANNEL
The decoherence model presented in section II in the form of the combined amplitude and phase damping channel cannot be efficiently implemented on a classical computer when multiqubit systems are considered. The quandary then becomes which quantum systems can be efficiently simulated using classical methods, if any at all.
Such a question is answered by the Gottesman-Knill theorem for stabilizer circuits 7 [32], [33]. The theorem states the following [32]: 5 The motivation behind QECCs is to extend the reliability times of qubits so that time demanding algorithms can be effectively executed in the quantum devices. 6 We assume that these QECC systems will be designed on classical computers due to the lack of practical quantum machines at the time of writing. 7 We will use the term ''stabilizer formalism'' to refer to the quantum circuits that are constructed according to the restrictions that are described in the Gottesman-Knill theorem.

Theorem 1 (Gottesman-Knill Theorem): Suppose a quantum computation is performed using only the following elements: state preparations in the computational basis, Hadamard gates, phase gates, controlled-NOT gates, Pauli gates, and measurements of observables in the Pauli group
(which includes measurement in the computational basis as a special case), together with the possibility of classical control conditioned on the outcomes of such measurements. Such computation can be efficiently simulated on a classical computer.
This theorem shows that certain quantum computations that involve very complex and highly entangled states are actually tractable problems for classical computers. The stabilizer formalism does not describe all possible quantum computations, but it does so for a sizeable quantity of these operations. For instance, notice that all the QECC families [20]- [31] that are constructed within the Quantum Stabilizer Code (QSC) framework [19] will fulfill the conditions of the theorem (including the entanglement-assisted versions), and so we will be able to implement them in regular computers.
As a consequence of Theorem 1, it is desirable to operate with a quantum channel whose dynamics are described by the stabilizer formalism, so that the error model of the system can be efficiently included in the classical simulation. In quantum information theory, the Pauli channel, N P , is the CPTP map that transforms the quantum state with density matrix ρ as (17) where the constants {p x , p y , p z } are interpreted as the probabilities of each specific operator impinging on the state ρ. The Kraus operators of this channel are √ p A A, with A ∈ P 1 .
The operation of each of the error operators on a qubit |ψ is: • Operator I leaves the state unchanged • Operator Z causes a phase-flip • Operator X causes a bit-flip • Operator Y causes both bit-flip and phase-flip 8 Following (17), we see that N P maps qubits |ψ onto a linear combination of the original state (I), the phase-flipped 8 Note that Y = iXZ. state (Z), the bit-flipped state (X) and the phase-andbit-flipped state (Y), where the sum is weighted by the probabilities p A , A ∈ P 1 . The symmetric instance of this channel 9 is a widely used model for decoherence in quantum information theory [20]- [31] and it is known as the depolarizing channel N D . The expression (17) then simplifies to and p receives the name of depolarizing probability. This symmetric instance of the Pauli channel represents the worst case scenario for the family of channels, as all the possible errors are equally likely to happen [47]. The generalization of the Pauli channel from (17) to an nqubit system is mathematically described by where p A is the probability for each of the n-fold Pauli operators to occur.
Observing (17), (22) and (23), it is clear that the dynamics of these channels are described by only using Pauli operators. By the Gottesman-Knill theorem (theorem 1), these quantum channel instances can be efficiently simulated on classical computers, and consequently, it is worthwhile to relate the decoherence processes described in section II, and mathematically modeled via the combined amplitude and phase damping channel, to the transformations described by the Pauli channel. It is clear that the Pauli channel is incapable of capturing the exact manner in which the combined amplitude and damping channel corrupts input states, as there is no combination of p x , p y and p z that makes N P (ρ) = N APD (ρ). However, we can use techniques of quantum information theory to engender approximations of N APD that are in the form of a Pauli channel.

B. APPROXIMATING QUANTUM CHANNELS WITH TWIRLING
In the prior section, we discussed the fact that quantum operations require an exponential number of parameters to completely describe them. Naturally, this leads to the conclusion that the simulation of such quantum constructs using classical resources cannot be efficiently performed. However, tools aiming at overcoming the issue of exponential parameter growth, 10 by successfully extracting useful partial information of the quantum dynamics, have been developed in the literature. This information extraction can be performed using a technique called twirling. Twirling is an extensively used method in quantum information theory to study the average effect of general quantum noise models via their mapping to more symmetric versions of themselves [35]- [40]. Generally, twirling a quantum channel N comes down to averaging the composition U † • N • U for unitary operators U(ρ) = UρU † that are randomly chosen according to some probability measure µ(U) [36], [39]. The resulting channel receives the name of twirled channel. A particularly important scenario, and the one we will consider in this article, is the case where the distribution over the unitaries µ(U) is discrete. For such a case, the dynamics of the twirled channel are given byN where µ(U l ) is a probability distribution over U l . Twirling a channel over some unitary is useful to analyze the average effects that said quantum channel produces. However, the serviceability of the information that the twirled channel provides is not clear, nor is the concept of how an approximated more symmetric channel constructed via twirling can be employed to design QECC schemes capable of fighting realistic decoherence. This conundrum is answered by the following lemma [36], which introduces a fundamental way in which the approximations of decoherence models can be realised to construct QECCs with the twirling method presented in this article.
Lemma 1: Any correctable code for the twirled channelN is a correctable code for the original channel N up to an additional unitary correction.
The proof of lemma 1 is given in [36]. The importance of the lemma resides in the idea that designing QECCs that correct errors for an approximated channel obtained by twirling will also correct errors of the original channel up to unitary correction. This allows us to design error correction codes for the twirled channel, which, due to the lemma, will still be valid for the original channel. The analysis conducted for the approximated channels will not be precise due to the fact that parts of the evolutions of the density matrices will be lost when operating with the twirled channels instead of the original ones. Nevertheless, work in the literature has shown that the estimations obtained using the twirling methods presented in this article are accurate [56].

1) PAULI TWIRL APPROXIMATION (PTA)
An extensively used type of twirl is the reputed Pauli twirl. This twirl consists in averaging the original channel N with the elements of the n-fold set of Pauli operators P n weighted in an equiprobable manner. Using (25), the Pauli twirled channel will then be given bȳ In other words, the off-diagonal elements of χ are eliminated [36], and the resulting channel has the form of an n-qubit Pauli channel as in (23). The coefficients χ P l , which are the probabilities of the errors of the Pauli channel, are calculated as in [37] where the coefficients α (k) P l come from writing the error operators of the Choi-Kraus decomposition in (3) in the Pauli basis as and are calculated as For this reason, applying the Pauli twirl to a general quantum channel N will result in a twirled channel with the structure of the Pauli channel. The probabilities for each of the elements of the n-qubit Pauli group can be calculated using (27) and (29). Additionally, we know from the Gottesman-Knill theorem that quantum CPTP maps with such a structure can be simulated efficiently on classical computers. Lemma 1 implies that we can design QECCs using this twirled approximation in order to construct QECCs for the original channel N , and that such error correction methods will maintain their error correcting capabilities over the original channel. As a result, we will be able to design QECC families for realistic quantum computers by designing and simulating them on classical computers based on approximating the channels by Pauli twirling.
We now proceed with the application of the Pauli twirl to the model of decoherence presented in section II in the form of the combined amplitude and phase damping channel N APD . For this purpose, we first write the transformation of the density matrix induced by the channel presented in (16), as a function of the Pauli matrices YρX, (30) where (15) is used to express the Kraus operators in the Pauli basis. As stated in (26), twirling the combined amplitude and 172630 VOLUME 8, 2020 phase damping channel by the Pauli operators P 1 will remove the off-diagonal elements of (30), resulting in the twirled channel ZρZ. (31) It is obvious that the channel shown in (31) exhibits the form of a Pauli channel (17). The probabilities for each of the possible Pauli errors in (17) are where (8) and (13) have been used in order to obtain the time dependencies of the damping and scattering probabilities [35], [40]. This asymmetric channel takes the form of the depolarizing channel (22) when the relaxation time and the dephasing time are the same T 1 = T 2 . Consequently, twirling a general quantum channel by the set of n-fold Pauli operators P n results in an approximated channel with the same form as a Pauli channel that will be efficiently implementable on a classical machine thanks to the Gottesman-Knill theorem (theorem 1). We will use the term Pauli Twirl Approximation (PTA) [35] to refer to the family of channels that are obtained using this symmetrization process. The process is simple, as it is based on rewriting the Kraus operators of the original channel in the Pauli basis and then eliminating the off-diagonal terms of the density matrix evolution equation.
With the goal of obtaining the PTA channels, it is interesting to define a widely used parameter when asymmetric Pauli channels are considered [40], [51], the asymmetry coefficient α. This parameter is introduced because after Pauli twirling the combined amplitude and phase damping channel, the resulting probabilities for the bit-flip (X) and phase-andbit-flip (Y) turn out to be equal, and so the degrees of freedom of the resulting PTA are reduced to two. Hence, it is possible to re-express the channel (17) using the parameter α, which is defined as the ratio [40] We can then rewrite the expression for the Pauli channel (17) as where p = p x +p y +p z . Note that α in (33) is a time-dependent parameter. However, if the coherent time duration t of the task is assumed to be negligible compared to the relaxation time T 1 , it can be approximated as Based on this approximation, the data from tables 1 and 2 can be used to determine the degree of asymmetry of existing quantum devices under the PTA. Table 3 shows these values, which indicate that the PTAs obtained for each of the existing technologies manifest drastically varying asymmetry levels. Note the vast difference between the channels, ranging from the depolarizing channel (α ≈ 1) to channels presenting very strong asymmetry levels (α ≈ 10 6 ).

2) CLIFFORD TWIRL APPROXIMATION (CTA)
Another twirling operation used for the purpose of channel approximation is known as Clifford twirling. For this twirl, the general channel N is averaged with the elements of the n-fold Clifford group C ⊗n 1 [57] weighted equiprobably, as for the Pauli twirl. Using the discrete twirling case as in (25), the Clifford twirled channel will be [37] where the operator-sum representation of N has been used to reach the second step and | · | refers to the cardinality of a set. Expression (36) can be further developed by noting that the n-fold Clifford group C ⊗n 1 is the semi-direct product of the n-fold set of Pauli operators P n and the symplectic group of dimension n, S ⊗n 1 [58], [59]. Thus, the elements C l ∈ C ⊗n 1 can be written as C l = P j S m with P j ∈ P n and S m ∈ S ⊗n 1 [37], [58]. This way, the Clifford twirl can be written as with S m ∈ S ⊗n 1 and P j ∈ P n . The cardinality of the Clifford group is |C ⊗n 1 | = |P n ||S ⊗n 1 |. VOLUME 8, 2020 If expression (37) is observed in detail, it can be concluded that Clifford twirling a quantum channel N is analogous to first Pauli twirling such a channel and then applying the symplectic twirl to the twirled channel [37], that is Consequently, the decomposition of twirls presented in (38) can be used in order to study the overall action of Clifford twirling. The action of each step can then be summarized as [37]: 1) P n -twirl: the effect of the Pauli twirl has been extensively studied in section III-B1, and so it is known that twirling the general quantum channel N by P n will lead to a twirled channel with the form of a Pauli channel as in (26). 2) S ⊗n 1 -twirl: the effect of the symplectic twirl on the PTA obtained in step 1 is the mapping of each of the non-identity Pauli operators to a uniform sum over the 3 non-identity Pauli operators. The application of the symplectic twirl to the PTA channel obtained in step one is: where the coefficients χ P j are obtained as in section III-B1. In order to analyze the effect of the symplectic twirl, we rewrite (39) as [37] N C ⊗n where the Pauli operators have been split by indexes ω which represent the weight 11 of a Pauli operator P j , ν ω counts the number of distinct ways that ω non-identity Pauli operators can be distributed over the n factor space, and i ω = {i 1 , · · · , i ω } with i β = {1, 2, 3} denotes which of the non-identity Pauli operators occupies the β th occupied site. Now fixing ω and ν ω , the action of the symplectic twirl for each term is so one can observe that the symplectic twirl maps each of the non-identity Pauli operators to a uniform sum 11 Number of non-identity factors of the tensor product.
over the 3 non-identity Pauli operators, as stated earlier.
For this reason, the expression for the Clifford twirled channel in (40) is with coefficients ξ ω,ν ω = 1 3 ω 3 ω i ω χ ω,ν ω ,i ω . From equation, (42) it is easy to deduce that applying the Clifford twirl to a general quantum channel N will produce a channel in the form of a Pauli channel, as in the PTA case, except for the presence of an additional degree of symmetrization, as the probabilities of Pauli operators that share the same weight and spatial position of the non-identity operators will be uniformly summed over the 3 non-identity Pauli operators. Through this symmetrization process, a quantum channel that requires the definition of less parameters is obtained, as the probabilities of the error operators that share the weight and spatial distribution of non-identity factors will now have the same probabilities of occurring, which was not true, in general, for the PTA channels. As in the PTA case, the twirled channels obtained by twirling a channel with the Clifford group will fulfill the Gottesman-Knill theorem, and so it will be possible to efficiently implement them on classical computers.
We now apply the Clifford twirl to the combined amplitude and phase damping channel N APD . As done for the PTA channels, we use the description of the channel in the Pauli basis of (30) and then apply the two steps introduced earlier to such a channel: 1) P n -twirl: the application of the Pauli twirl has been addressed in section III-B1, and it leaves N APD in the form of a Pauli channel like the one given in (31). 2) S ⊗n 1 -twirl: from the previous discussion, we know that the symplectic twirl maps each of the non-identity Pauli operators to the uniform sum over the three non-identity Pauli operators. In consequence, and following (42), the Pauli channel of (31) results in is the depolarizing probability and comes from summing the p x , p y , p z from (32). It is clear that (43) has the form of the symmetric instance of the Pauli channel. In other words, it has the structure of the so-called depolarizing channel of (22). Moreover, using the temporal expressions of the probabilities of the PTA from (32), the depolarizing probability of the Clifford twirled channel is related to the physical parameters of a quantum device as where T 1 and T 2 are the relaxation and dephasing times, as presented before. Consequently, twirling a general quantum noise model N by the n-fold Clifford group C ⊗ 1 will result in an approximated channel with the form of a Pauli channel whose principal trait is that the error operators that have the same weight and non-identity operator distribution in the tensor product will have the same probability of occurring. For the one qubit case, this implies that the resulting channel will be a depolarizing channel. We will use the term Clifford Twirl Approximation (CTA) to refer to the family of channels obtained with this symmetrization method.
As indicated before, the most important aspect of the PTA and CTA channels is the fact that they are efficiently implementable on classical computers, as they fulfill the Gottesman-Knill theorem (theorem 1). This is especially interesting in the field of QECC design, given that from lemma 1 we will be able to design these codes without the need for actual quantum computers.

IV. QUANTUM MEMORYLESS CHANNELS AND CHANNELS WITH MEMORY
In the previous section, we presented a symmetrization method called twirling that can be used to approximate general quantum channels to PTA and CTA Pauli channels that can be efficiently implemented on classical computers. Generally speaking, PTA and CTA channels act collectively on n qubits, see (31) and (43), and they are both reduced to n-qubit Pauli channels with different degrees of symmetrization. More degrees of symmetrization can be obtained for the n-qubit twirled channels with extra twirling operations such as permutation twirls [36], [37]. These twirled channels are also efficiently implementable, but in order to obtain the probability distributions that define the specific Pauli channel that comes out as a PTA or CTA, the original channel N must be modeled. The interaction processes between the elements of an n-qubit system can be quite subtle, and, as stated in section II, it is reasonable to assume that, if certain conditions apply, each of the qubits of the system will interact with the environment in an independent and similar way. Additionally, simplified memory effects can be employed to express the interactions that each of the qubits has with its counterparts. Consequently, we will be able to use the expressions obtained for the PTA and CTA of the combined amplitude and phase damping channel for 1 qubit in order to simulate decoherence over complex n-qubit quantum states.
Most of the research related to quantum channels considers a scenario in which the corruption of the input quantum states occur independently and identically [20]- [31], [40], [51]. If two different quantum states are transmitted through the same channel, the noise applied to the first transmitted state is assumed to be independent of the noise applied to the second state. This can be seen as the channel having no ''memory'' of previous events, which is why such a configuration is known as a memoryless quantum channel. The memoryless assumption simplifies the noise induced input-output mapping of quantum states, and provides an accurate portrayal of particular physical events. For instance, communication schemes whose signalling rate is low to allow the environment to reset, or scenarios in which a magnetic field is applied to reset the memory of the channel, can be accurately modeled using the memoryless configuration [41].
However, memoryless quantum channels cannot be used to model all quantum communication systems. In optical fibers (a typical medium to transmit quantum information), sufficiently high signalling rates cause the environment to respond to each successive transmission in a way that is correlated to previous ones [60], [61]. In quantum information processors, quantum bits can be so closely spaced that cross-talk may occur and the channel noise might be correlated [62], [63]. These two scenarios exemplify situations in which applying the popular memoryless noise configuration would not provide a realistic description of the corresponding physical events. Instead, these systems require a quantum channel model that integrates memory effects. We introduce such a model in the sequel, where we begin by discussing memoryless quantum channels and proceed by extending them to the more general model of quantum channels exhibiting memory.

A. MEMORYLESS QUANTUM CHANNELS
In [64], the authors describe a specific scenario in which applying a memoryless quantum channel is appropriate. They assume that the sequence of signals is transmitted over the physical medium in an ordered fashion and at a constant speed and that the environment undergoes a dissipative process which resets it to a stable configuration on a timescale τ , known as the channel relaxation time 12 of the medium. Then, if the rate at which each successive signal is transmitted is much lower than the inverse of the channel relaxation time, i.e R < 1 τ , the physical event can be represented by a memoryless quantum channel.
Generally speaking, a memoryless channel affecting a system of n qubits is comprised by the tensor product of the channels affecting each of the individual qubits of the system. Mathematically, this is expressed as where N j refers to the noise channels that individually affect each of the qubits. Furthermore, it is a common assumption when dealing with this type of channels that the interaction that each qubit has with the environment is identical.

Thus, the memoryless expression (45) is reduced to
where now N is the one-qubit decoherence model that will equally affect each of the elements of the quantum system. Therefore, for the decoherence model based on the combined amplitude and phase damping channel, the memoryless channel for an n-qubit system will be expressed as where N APD is described by the error operators in (15). However, as explained in section III, we will not be able to efficiently implement this channel on a classical computer. The solution is to work with the corresponding PTA and CTA channels, which are known to fulfill the Gottesman-Knill theorem and to have the structure of n-qubit Pauli channels.
As seen in sections III-B1 and III-B2, the one-qubit channels obtained by twirling the combined amplitude and phase damping channel are the Pauli channel for the PTA case and its symmetric instance or depolarizing channel for the CTA. Accordingly, the n-qubit channels that arise in an identical and memoryless manner from these channels will have the structure of the n-qubit Pauli channels in (23) (48) with the probability distribution p A defined as Since A ∈ P n is A = n j=1 A j , A j ∈ P 1 . Thus, the total probability of the word will be the product of the individual probabilities because the channel acts independently on each of the qubits.

B. QUANTUM CHANNELS WITH MEMORY
Quantum information can sometimes be exposed to physical events that exhibit spatial or temporal correlations. The outcome of the general twirling approximations of n-qubit channels also shows features related to these phenomena. To obtain simplified channels capable of representing such interactions, twirled channels with the ability to integrate memory effects are necessary. For this purpose, we must analyze how to model correlations in PTA and CTA channels. This implies that we need to model the conditioned probability distribution P(A n |A n−1 ⊗ A n−2 ⊗ · · · ⊗ A 1 ), A j ∈ P 1 between the current 1-qubit PTA or CTA that will be conditioned by the last channel use. This way, the expression for the general n-fold Pauli channel described in (23) can be written as (50) where A = A n ⊗ A n−1 ⊗ · · · ⊗ A 1 with A j ∈ P 1 refers to each of the possible n-fold Pauli operators. Note that from (50), every PTA or CTA approximation can be seen as a channel with memory. We will now show how to model memory effects considering the 1-qubit approximations to derive memory models that can be incorporated in such efficiently implementable channels.
The most studied class of quantum channels with memory is the family of channels with Markovian correlated noise [65]- [67], which considers noise models in which quantum objects are transformed via the application of maps whose elements are randomly generated by a classical Markov process. Markov processes describe sequences of events in which the probability of each event depends only on the previous event, i.e. P(A j |A j−1 ⊗ A j−2 ⊗ · · · ⊗ A 1 ) = P(A j |A j−1 ). For this type of channel, the correlation between the single qubit Pauli operators of the set P 1 can be described by means of a 4-state Markov chain [68]- [71]. This can be seen in Figure 2, where each state corresponds to one of the single qubit Pauli operators. The transition probability from a previous error state A j−1 to the current error state A j is denoted by A possible way to capture the temporal correlation of Pauli operators over a channel with memory is to introduce the so called correlation parameter µ ∈ [0, 1], where µ = 0 indicates zero correlation and µ = 1 indicates perfect correlation. This correlation parameter was originally presented in [72], where the authors derive the first comprehensive quantum information characterization of memory effects in a continuous variable setup. Specifically, a continuous variable model of quantum memory channels that accurately portrays the transmission of quantum objects (photons encoded with information) through optical fibers characterized by finite relaxation times is proposed. The model describes each channel use as an independent bosonic mode, and along with two parameters, enables the description of multiple communication scenarios. In terms of the time delay between each successive transmission of a photon t and the relaxation time τ of the optical fiber, the model enables the representation of a memoryless configuration ( t τ ), intersymbol interference memory ( t ∼ τ ) [73], or a perfect memory configuration ( t τ ) [74]. As shown in [72] and [75], the transmissivity of the optical fiber (in the beam-splitter configuration considered for the model) plays the role of a memory parameter and it can be approximated by ≈ e − t τ . For more abstract memory modeling, we can define the correlation parameter µ and assume, as in previous works [64], [71], [75], that it can be quantified in terms of the time delay between successive channel uses and the relaxation time of the environment: µ ≈ e − t τ . A model introduced in [67] considers that two successive channel uses are related by the transition probability q A j |A j−1 defined as where δ A j−1 A j is the Kronecker delta function and p A j is the probability that Pauli operator A j ∈ P 1 is imposed on the transmitted quantum state. As a consequence, considering a quantum channel with memory defined by Markovian correlated noise, the joint probability p A is given by for each n-fold Pauli operator A ∈ P n . The resulting Markovian Pauli channel N MP will thus have the following structure where A = A n ⊗ A n−1 ⊗ · · · ⊗ A 1 with A j ∈ P 1 . In this manner, Markovian Pauli memory channels, N MP , can be created from the PTA and CTA approximations, which allow for efficient implementation on a classical computer.
General Markovian quantum channels with memory can be built by applying the memory effects to the error-sum operator representation of general quantum channels [42], [64]. These models will represent the general effects that decoherence processes with memory generate over quantum information with increased precision. An increase in the research interest regarding this topic is expected, including the development of memory models that go beyond the Markovian paradigm. Nonetheless, at the time of writing, the discussion provided in this section completely describes the current state of affairs with regard to channels that are efficiently implementable on classical computers and that have been used for QECC construction.

V. IMPLEMENTING PAULI CHANNELS FOR QECC DESIGN AND SIMULATION ON CLASSICAL COMPUTERS
Previous sections of this article have described how to model decoherence processes mathematically, and how to obtain approximated channels from these general models so that they can be appropriately simulated on classical computers. This is a significant step forward in the endeavour of the quantum information community to construct QECC families for the post-NISQ-era, before quantum machines that can take advantage of these error correction methods can be constructed. The methods discussed here employ the quantum information technique referred to as twirling to obtain channels with the structure of Pauli channels. The reason for doing so is that such channels fulfill the conditions of the Gottesman-Knill theorem, and therefore they will be implementable on classical computers. In consequence, the twirled channels we construct will be implementable in run-of-the-mill classical machines available at this moment in time.
Design and simulation of error correction methods for the Pauli channel have been exhaustively studied by the quantum information community and are well documented [20]- [31], [40], [51], [55], [71], [76]. The key to implement Pauli channels on classical machines is the Pauli-to-binary isomorphism, which maps elements of the Pauli operator set onto binary strings [47], [76]. Establishing this isomorphism between these groups enables researchers to assess the performance of designed QECC families via Monte Carlo simulations, as described in [28], [29], in a fashion reminiscent of classical coding theory.
In this section we describe the Pauli-to-binary isomorphism and the way in which this mapping can be used to simulate the performance of QECCs on classical computers. With this goal in mind, we also discuss basic stabilizer coding, syndrome measurements, and error discretization. Additionally, we provide a simple example to see how the Word Error Rate (WER) and the QuBit Error Rate (QBER) of a quantum error correction method can be obtained via Monte Carlo simulation. Both WER and QBER are used as figures of merit when gauging the error correcting capability of the QECC schemes.

A. ERROR DISCRETIZATION
Before presenting the theory of stabilizer codes, some important results about the general theory of quantum error correction must be discussed, so that the construction of the stabilizer family can be introduced in a logical and coherent manner.
We begin with the following theorem, which defines the conditions that a quantum error correcting code must fulfill to protect quantum information from a particular noise process N . We refer to this theorem as the quantum error correction conditions or the Knill-Lafflame theorem [32], [78].
Theorem 2 (Knill-Lafflame theorem): Let C ⊂ H ⊗n 2 be a quantum code defined as a subspace of the n-dimensional Hilbert space that is the state space of the n-qubit system to be protected, and let P be the projector onto C. Let N be a noise operation defined by Kraus operators {E j }. A necessary and sufficient condition for the existence of an error-correction FIGURE 3. General schematic of a stabilizer QECC. U is the unitary that performs the encoding operation, and takes |ψ ∈ H ⊗k 2 to the codespace as |ψ ∈ C(S) ⊂ H ⊗n 2 with the aid of the ancilla qubits |0 ⊗n−k . The transmitted codeword is subjected to the noise operation N resulting in the noisy codeword |ψ N ∈ H ⊗n 2 . This quantum state is processed by the syndrome extractor, so that the syndromes associated to the error operator can be obtained. The quantum state is not altered by such an operation, and the syndrome is processed by a classical algorithm to estimate a recovery operation that will be applied to the quantum state after going through the inverse encoder U † . This way we obtain the corrected state |ψ ∈ H ⊗k 2 .
operation R correcting N on C is that for some Hermitian matrix of complex numbers α. We call the Kraus operators {E j } errors. If the recovery operation R exists, we say that {E j } constitutes a correctable set of errors.
Theorem 2 defines the precise conditions that a method designed to suppress quantum errors must fulfill in order to correct a specific error operator N . For the particular approximations that permit efficient simulation of these quantum channels, the error operators are n-qubit Pauli operators. Hence, it appears that the ability these codes have to correct errors will be limited to a set of Pauli errors, which, at first glance, seems to be an important constraint. However, QECCs are actually significantly more powerful due to a theorem known as discretization of errors [32], [78].
Theorem 3 (Discretization of Errors): Suppose C is a quantum code and R is the error-correction operation to perform recovery from a noise process N with Kraus operators {E j }. Suppose that F is another noise process with Kraus operators {F k } which are linear combinations of the E j operators, i.e. F k = j ξ kj E j for some matrix ξ of complex numbers. Then the error correction operation R also corrects for the effects of the noise process F on the code C.
The primordial consequence of theorem 3 is that codes constructed for a specific noise process will also be able to correct noise events whose error operators are linear combinations of the error operators of the original noise process. This result is momentous, as it proves that the constructed codes will be able to correct an infinite set of noise processes, while only considering a finite amount of Kraus operators during their analysis and design. For the approximated channels presented in this article, the discretization of errors implies that all errors that are linear combinations of the Pauli errors that can be corrected by a quantum code will also be correctable by the code.
This result is invaluable when designing and simulating QECCs efficiently on a classical computer, since it allows us to do so by just considering discrete sets of errors represented by n-fold Pauli operators. Considering that the set of errors that have to be tested is discrete, it will be possible to represent the errors using binary notation, i.e. in a way that classical computers can manage.

B. STABILIZER CODES
The set of n-fold Pauli operators P n together with the overall factors {±1, ±i} forms a group under multiplication. This group is known as the n-fold Pauli group G n [55]. A stabilizer code is defined by an abelian subgroup S ⊂ G n . For a [[n, k, d]] unassisted 13 stabilizer code encoding k logical qubits into n physical qubits with distance d, the stabilizer set has 2 n−k distinct elements up to an overall phase, and it is generated by n − k independent generators. 14 The stabilizer code C(S) associated with the stabilizer set S is then defined as i.e. the simultaneous +1-eigenspace defined by the elements 15 of S. With H ⊗n 2 we refer to the complex Hilbert space of dimension 2 n , which is the state space of systems formed by n qubits. Figure 3 presents the general scheme for stabilizer codes that correct the noise operation N . The unitary U maps the input information word |ψ ∈ H ⊗k 2 to the codespace |ψ ∈ C(S) ⊂ H ⊗n 2 with the aid of the ancilla qubits |0 ⊗n−k . The existence of such a unitary that takes the arbitrary input states 13 In this article we discuss stabilizer coding without entanglementassistance. However, the methods presented here are equally valid for the entanglement-assisted formulation.
14 Stabilizer codes are represented by n − k generators, given that the rest are combinations of them. This provides a compact representation of the code. 15 Note that the simultaneous eigenspace is generated by the n − k independent generators, and so the code is defined by them. |ψ |0 ⊗n−k to the codespace is guaranteed [55]. The encoded state then experiences the action of a quantum channel N , which is described as one of the approximated twirl channels with the structure of a Pauli channel that have been introduced earlier in this work. The noisy quantum state, |ψ N , must then be corrected by a recovery operation. The fundamentals of quantum mechanics establish that measuring a quantum state forces its superposition state to collapse, which causes the loss of the quantum information contained in the original state. Therefore, the theory of quantum error correction codes must circumvent this issue. It achieves this by measuring the syndrome of the errors in an indirect way, avoiding the destruction of the quantum state and allowing us to garner information about the error that will be used to estimate the best recovery operation R.
The error syndromes is defined as a binary vector of length n − k that captures the commutation relationship between the generators of the stabilizer set S and the error E ∈ G n that perturbed the encoded quantum state as a consequence of the channel action. The error operators are assumed to be elements of the Pauli group, since the twirled approximated channels have the structure of a Pauli channel. From error discretization, we know that the designed codes will also correct errors that are linear combinations of the errors which they were originally intended to correct. 16 It is common knowledge that any two elements of the n-fold Pauli group either commute or anticommute. For this reason, the error operator E will either commute or anticommute with each of the generators S j , j ∈ {1, · · · , n − k}, of the stabilizer set, and the components s j of the error syndrome capture their relationship as To construct a circuit that can extract the syndrome for each of the generators S j , consider the received noisy quantum state |ψ N = E|ψ . Then the vector |ψ N is an eigenstate of each of the generators of the stabilizer set S associated to the ±1 eigenvalues, i.e.
where the particular value out of the two possible eigenvalues depends on the commutation relationship between the channel error and the stabilizer generator. As a result, in order to determine the syndrome of the error that takes place in a specific channel instance, the eigenvalue of S j must be measured [55]. Figure 4 presents a circuit built to perform such measurements.
Once the error syndrome is obtained, this classical information is used to decode or estimate the recovery operation that will correct the corrupted quantum information state. Syndrome decoders of quantum stabilizer codes depend on the specific QECC family construction [20]- [31], [40], [51], [55], [71], [76]. In general, optimal decoding of quantum stabilizer codes must consider the phenomenon known as degeneracy [77], a feature which has no equivalence in classical coding. In degenerate quantum codes, there are sets of different error operators that have the same effect in the transmitted codeword. Although it should be possible to exploit this property in the decoding process (and also for the design of better codes), the design of decoders that can efficiently perform the so-called Degenerate Quantum Maximum Likelihood Decoding (DQMLD) remains unsolved for stabilizer codes in general. 17 Current decoding of these codes is approached in terms of Quantum Maximum Likelihood Decoding (QMLD), which ignores error degeneracy and undertakes the decoding task as in the field of classical linear block codes. Building effective degeneracy-exploiting decoders is out of the scope of this article, so it will not be discussed further.
The last steps in the error correction operation are the application of the inverse of the encoding operation (unitary U † ), followed by the use of the recovery operation R, which depends on the results of the syndrome decoder.

C. PAULI-TO-BINARY ISOMORPHISM
In this section we present the Pauli-to-binary isomorphism as a way to express the operation of QECCs as binary vectors, which are generally referred to as symplectic strings. We start by establishing the map for the set of Pauli matrices P 1 . Notice that Pauli matrices are Hermitian unitary matrices with eigenvalues equal to ±1.   As mentioned earlier, we are interested in expressing the Pauli matrices in binary notation. Consider the binary words of length two (F 2 ) 2 = {00, 01, 10, 11}, with the usual modulo-2 addition as the field operation. Table 6 portrays the relationships between the elements of this Abelian group under modulo-2 addition. In the following, the notation u = (z|x), with z, x ∈ F 2 will be used to refer to the elements of (F 2 ) 2 . The group is also a vector space over the field F 2 , over which the bilinear form called symplectic form or symplectic product can be defined [76]. The symplectic product is the application : where u = (z|x) and v = (z |x ). The results of computing the symplectic product for the elements of (F 2 ) 2 are presented in table 7. We proceed by defining the map ϒ : (F 2 ) 2 → P 1 , which can be seen in table 8. Upon closer inspection of this mapping, it is easy to observe that it is specifically chosen so that the elements ϒ (z|x) and Z z X x are equal up to a phase factor, i.e. so that they are equal under the equivalence class [·], as defined before. As a result, and after closely examining tables 4, 5, 6 and 7, we reach the following key conclusions [76]: • After thorough observation of both tables 5 and 6, we see that the map [ϒ] : (F 2 ) 2 → [G 1 ] is an isomorphism as • Analyzing tables 4 and 7, we observe that the commutation relationships between the Pauli matrices are captured by the symplectic product for binary vectors of length two as Note that by making use of this map we are able to represent Pauli operators as binary vectors of length two, maintaining the ability to capture the commutation relationships via the symplectic product.
The next step is to extend the isomorphism to the n-fold set of Pauli operators in order to be able to represent systems that have n qubits. The elements of the n-fold set of Pauli operators P n are tensor products of individual Pauli matrices, and so the equivalence class [·] for these elements of the Pauli group will be defined in the same way, i.e. the set [G n ] = {[A] : A ∈ P n }. The operation of this equivalence class will be the product defined as , and this set will also be an Abelian group under multiplication, as it was the case for the 1 qubit scenario.
Having defined the n-fold effective Pauli group, we now consider the Abelian group under modulo-2 summation and the vector space composed by the binary vectors of length 2n, (F 2 ) 2n . As before, the elements of the group will be written asū = (z|x), wherez = z 1 · · · z n ∈ (F 2 ) n andx = x 1 · · · x n ∈ (F 2 ) n . Now the symplectic product for these binary strings of length 2n can be defined as an application : whereū = (z|x),v = (z |x ), u j = (z j |x j ) ∈ (F 2 ) 2 , v j = (z j |x j ) ∈ (F 2 ) 2 and T denotes the transpose. Notice that we are representing the strings as row vectors. Consequently, the symplectic inner product of the 2n length binary string will be the Boolean sum of the symplectic product of the n (F 2 ) 2 elements that form such a string. The map ϒ is now performed individually for each of the elements of the n-fold tensor product that form each of the elements of P n , namely, the map ϒ : (F 2 ) 2n → P n is taken as ϒū = ϒ u 1 ⊗ϒ u 2 ⊗· · ·⊗ϒ u n . It is trivial to see that this map has been selected so that ϒ (z|x) is equal to Z¯zXx = Z z 1 X x 1 ⊗ Z z 2 X x 2 ⊗ · · · ⊗ Z z n X x n up to an overall phase.
Considering this discussion, we can obtain similar key conclusions for the n-fold map as we did earlier for the Pauli matrices [76]: • The commutation relationships of the n-fold Pauli matrices are captured by the symplectic product as Thanks to these results, we will be able to represent the nfold Pauli matrices, which act as error operators of the Pauli channels, as binary vectors of length 2n. Stabilizer codes can be described via binary matrices using this map, as the generators of the stabilizer set S that define the code are nfold Pauli matrices. The binary representation of a stabilizer code receives the name of parity check matrix (PCM) as in classical coding. The PCM of a stabilizer code is then where each of the rows of the matrix is obtained by applying the map ϒ to each of the stabilizer generators. Moreover, as the symplectic inner product captures the commutation relationships that the n-fold Pauli matrices share between themselves, this binary representation can be used in order to extract the error syndrome that Pauli errors have with respect to the specific stabilizer code represented by the PCM H . The syndromes(ē) of a specific error E whose binary representation is given byē, is calculated as where the symplectic product between the matrix and the error vector is performed row-wise, that is to say, each of the components of the syndrome is obtained by computing the symplectic product of the associated row of the PCM and the error operator. As a direct result of this, we will be able to evaluate the performance of the target stabilizer code by conducting Monte Carlo simulations [28], [29]. An error correction round of the system can be broken down into the following steps: 1) Generate a binary error pattern,ē, of length 2n, following the probability distribution of the specific Pauli channel under consideration. The map ϒ is used to obtain the probability distributions of the binary vectors. 2) Calculate the syndrome associated to such an errors(ē) by using the symplectic product as in (65).
3) Run the specific syndrome decoder that has been implemented for the particular stabilizer code being used as a QECC. The syndrome decoder depends on the stabilizer code under consideration [20]- [31], [40], [51], [55], [71], [76], but it is always a purely classical algorithm. The quantum operation implemented for final recovery depends on the result of the syndrome decoder. 4) Check if the estimated recovery operation in this error correction round has been successful in its attempt to revert the channel error. This depends on whether error degeneracy is taken into account when defining the ''success'' of the decoding scheme [77]. This happens because for QECCs, the physical error that is sometimes estimated does not always match the error that actually occurred during transmission, but its associated logical error does indeed match the logical error related to the real physical error. This implies that even if the physical error has not been correctly estimated, quantum information will be successfully recovered.
Multiple error correction rounds are repeated to obtain metrics that quantify the error correcting capability of the code under consideration. The most common metrics used to evaluate the quality of a quantum error correcting code are: • Quantum Word Error Rate: it represents the probability that at least one qubit in the block is incorrectly decoded, i.e., a decoding failure is accounted for if just one of the operators of the length n error operators is incorrectly estimated.
• QuBit Error Rate: it is the probability that an individual qubit is incorrectly decoded, i.e. in each transmission round, each of the operators of the n length error operator is individually considered. Note that, as explained in point 4, these operational figures will depend on whether the physical errors or their associated logical errors are considered. Both figures of merit are used in the literature when assessing the performance of stabilizer codes [20]- [31], [40], [51], [55], [71], [76], and in both cases the WER and QBER obtained without considering degeneracy will be upper bounds for the WER and QBER that could be achieved if degeneracy were exploited. The particular number of iterations required to obtain precise estimates of these metrics is given by the theory of Monte Carlo simulations [28], [29].

D. EXAMPLE: [[5,1,3]] STABILIZER CODE
In previous sections we have explained how to approximate quantum channels to channels that are efficiently implementable on classical computers. We have also discussed how the performance of stabilizer codes can be simulated under the action of such channels by using the Pauli-to-binary isomorphism. In this section we present a simple example of an error correction round, aiming at providing a clear and Lookup table used to perform syndrome decoding of the stabilizer code considered in the example. Since the channel is assumed to be a memoryless depolarizing channel and degeneracy is ignored, the most probable error for each syndrome will be the one with the lowest weight that produces said syndrome. transparent portrayal of how classical simulation of QECCs works.
Consider the [ [5,1,3]] stabilizer code defined by the following stabilizer generators: The next step is to generate binary error patterns and then obtain their syndromes so that error estimations can be made from those syndromes. The probability distributions used to simulate the channel are the distributions for the twirled approximated channels presented in sections III and IV, where the n-fold Pauli operators are mapped to binary vectors via the map ϒ. For our particular example we consider a memoryless depolarizing channel. A decoding round will be performed for each one of the following error patterns: •ē 1 = 0 0 0 0 0 0 0 1 0 0 which maps to the physical error E 1 = I ⊗ I ⊗ X ⊗ I ⊗ I.
•ē 3 = 1 0 1 0 1 0 1 1 0 0 which maps to the physical error E 3 = Z ⊗ X ⊗ Y ⊗ I ⊗ Z. The associated error syndromess(ē j ) are calculated as shown in equation (65), which for the errors we are considering will be: •s(ē 1 ) = H ē 1 = 1 1 0 1 . •s(ē 2 ) = H ē 2 = 0 1 1 1 . •s(ē 3 ) = H ē 3 = 0 0 0 0 . Next, the decoding algorithm must be fed the error syndromes so that the best recovery operations considering the syndrome information can be estimated. Decoding stabilizer codes is substantially nuanced given that the optimal task differs from the classical decoding version due to degeneracy [77], and it is heavily conditioned by the code and the channel model under consideration. Nonetheless, QMLD is still a valid decoding algorithm even if it is not the optimal one. For this example, we will ignore degeneracy for the estimation of the recovery operation and so the syndrome decoder will be given by the lookup table shown in table 9.
Notice that the recovery operations in the lookup table  (table 9) are Pauli operators. This is because n-fold Pauli operators are Hermitian and unitary matrices, implying that A 2 = I, A ∈ P n . Consequently, the decoder estimates the QMLD most probable physical operator associated to the measured syndrome, and it applies such a Pauli operator to the noisy quantum information |ψ N so that the action of the channel is suppressed. 18 The recovery operation for the error words considered for our example are obtained from table 9: •s(ē 1 ) = 1 1 0 1 → R = I ⊗2 ⊗ X ⊗ I ⊗2 . The resulting state after recovery will be (I ⊗2 ⊗ X ⊗ I ⊗2 )|ψ N = (I ⊗2 ⊗ X ⊗ I ⊗2 )E 1 |ψ = |ψ . As a result, the code will succeed in correcting the channel error E 1 .
•s(ē 2 ) = 0 1 1 1 → R = I ⊗ Z ⊗ I ⊗3 . The resulting state after recovery will be (I⊗Z⊗I ⊗3 )|ψ N = (I ⊗ Z ⊗ I ⊗3 )E 2 |ψ = (I ⊗ Z ⊗ I ⊗3 )(Y ⊗ X ⊗ I ⊗2 ⊗ Z)|ψ = (Y ⊗ Y ⊗ I ⊗ I ⊗ Z)|ψ . As can be seen, the error correction operation is unsuccessful and so this channel error produces a decoding failure. Following recovery, this iteration results in an error event for the WER. However, it will be reflected as 4 error events and one success event for the QBER, since the 4 th qubit will not suffer from an error (the operator acting on it is the identity operator).
•s(ē 3 ) = 0 0 0 0 → R = I ⊗5 . The resulting state after recovery will be (I ⊗5 )|ψ N = (I ⊗5 )E 3 |ψ = (Z ⊗ X ⊗ Y ⊗ I ⊗ Z)|ψ . At first glance, it looks like the decoder has failed to recover the correct quantum information state, as the estimated error and the one that actually occurred are different, which results in a non-identity operator acting on the encoded state. This would entail the addition of an error event for the WER and 4 error events for the QBER. However, a closer look at the error operator that acts on the quantum state that lays in the codespace shows that it is actually the generator S 2 of the stabilizer set that defines the code.
As the codespace is defined as the +1-simultaneous eigenspace of the stabilizer generators, the statement (Z ⊗ X ⊗ Y ⊗ I ⊗ Z)|ψ = |ψ will be true, and so the resulting state will be error-free. This last error belongs 18 Note that Pauli matrix product is the modulo 2 sum for the binary representation under the map ϒ.
to the class of errors known as degenerate errors. Hence, if degeneracy is being considered for the computation of the performance metrics, as it should for optimal decoding, this last event will count towards a successful decoding round.
By randomizing the generated error pattern according to the probability distributions of the approximated twirl channels presented earlier, we will be able to simulate the performance of stabilizer codes for channel models that provide a realistic representation of decoherence processes. Lastly, as mentioned previously, to derive accurate estimates of the performance metrics, the rules of Monte Carlo simulations must be followed [28], [29].

VI. CONCLUSION
We have described the most basic physical phenomena that cause the corruption of quantum information by considering the basic interactions that define decoherence. The way in which these effects are mathematically modeled via CPTP operators called quantum channels has also been presented. Then, we discussed the issue of efficient implementation of general quantum channels on classical computers by approximating these channels to versions that fullfil the Gottesman-Knill theorem. To obtain these approximations, we use Pauli and Clifford twirling techniques, which approximate channels into different instances of Pauli channels with varying degrees of asymmetry. Once the twirled approximations have been introduced, we discussed the memoryless instances of quantum channels, as well as the inclusion of memory effects to provide a more realistic depiction of decoherence. Specifically, we introduced the conditions that the transmission rate must satisfy so that the memoryless instance of the quantum channel provides an appropriate representation of reality. Then, a Markovian-based modelling of channels that do not fulfill the memoryless conditions is presented. This modelling strategy relates the correlation parameter with the transmission time spacing that is present between different qubits. Using Markovian memory to model memory effects is the de facto approach for quantum channels in the literature, but more advanced models are expected in future research.
After reviewing decoherence and its modelling via approximated quantum channels, we provide an overview of stabilizer codes, introducing the way in which they are simulated and assessed on classical computers. We begin by presenting the error correction conditions and the discretizing of errors for general QECCs. Afterwards, we briefly review the encoding, the syndrome measurement and the decoding as well as the recovery operations in stabilizer codes. Then, we introduce the Pauli-to-binary isomorphism and the way in which this mapping can be employed in order to simulate the WER and QBER of a specific stabilizer code. Finally, we provide an example of a decoding round for different error words in a specific stabilizer code, aiming at clarifying the earlier theoretical description.
In summary, we have reviewed how to efficiently test stabilizer codes on classical computers for realistic approximated decoherence models obtained using twirling techniques. This allows the design and testing of complex error correction models that will be key in post-NISQ-era quantum computers and devices.