We are given: “In quantum error correction, you encode quantum states into logical states made of many qubits in order to improve their resilience to errors. In quantum error detection, you do the same but can only detect the presence of errors and not correct them. In this problem, we will consider a single [[4,2,2]] quantum error detection code, which encodes two logical qubits into four physical qubits, and investigate how robust logical quantum operations in this code are to quantum errors.

Our convention is that the four physical qubits in the [[4,2,2]] code are labelled 0,1,2,3. The two logical qubits are labelled A and B. The stabilizers are \(XXXX\) and \(ZZZZ\), where \(X\) and \(Z\) are Pauli matrices. The logical \(X\) and \(Z\) operators on the two qubits are \(X_A = XIXI\), \(X_B=XXII\), \(Z_A = ZZII\), \(Z_B = ZIZI\), up to multiplication by stabilizers.

We will consider different state preparation circuits consisting of controlled not \(CNOT_{ij}\) gates, where \(CNOT_{ij}\) has control qubit \(i\) and target qubit \(j\). As a simple model of quantum errors in hardware, we will suppose that each \(CNOT_{ij}\) gate in the circuit has a two qubit depolarizing error channel following it that produces one of the 15 non-identity two-qubit Paulis with equal probability \(p/15\). The probability \(p\) indicates the probability of an error in a single two-qubit gate. We will assess the logical infidelity of certain state preparation protocols as a function of the physical infidelity \(p\).

Checkpoint 1

Suppose that we wish to prepare a logical two-qubit GHZ state \((|00\rangle_{AB}+|11\rangle_{AB})/\sqrt{2}\) in the [[4,2,2]] code. To do so, we use the following state preparation circuit:

Note that this equation is written in matrix multiplication order, while the quantum operations in the circuit occur in the reverse order (from right-to-left in the above equation). \(H\) is a single-qubit Hadamard gate.

What is the physical state fidelity of the final physical 4-qubit state at the end of the circuit as a function of the two-qubit gate error rate \(p\)?”

Interpretation: We have a circuit that prepares a logical GHZ state. The circuit is given as: (CNOT_{03}) (H_0) (CNOT_{21}) (H_2). But careful: The equation is written in matrix multiplication order, but the quantum operations occur in reverse order. That means the actual circuit order from initial state to final state is: first apply H_2, then CNOT_{21}, then H_0, then CNOT_{03}. So the circuit steps: Step 1: Apply H on qubit 2. Step 2: Apply CNOT with control qubit 2 and target qubit 1. Step 3: Apply H on qubit 0. Step 4: Apply CNOT with control qubit 0 and target qubit 3.

We assume initial state? Typically, we start from |0000>? But careful: The problem says: “prepare a logical two-qubit GHZ state (|00>{AB}+|11>{AB})/√2”. So the intended final state is that logical state. So we want to compute the fidelity of the actual state produced by the circuit when each CNOT gate is followed by a depolarizing error channel. But note: There are two CNOT gates: CNOT_{21} and CNOT_{03}. And there are Hadamard gates: H_2 and H_0. The problem says: “each CNOT_{ij} gate in the circuit has a two qubit depolarizing error channel following it that produces one of the 15 non-identity two-qubit Paulis with equal probability p/15.” So the errors occur on the two-qubit gates. The Hadamard gates are single-qubit gates. They are assumed perfect? Probably yes. So only errors on CNOT gates.

So we have two CNOT gates. Each CNOT gate is followed by an error channel that with probability (1-p) does nothing, and with probability p, applies a random two-qubit Pauli error from the set of 15 non-identity Paulis on the two qubits involved. And these errors are applied after the CNOT gate. So the circuit becomes: Start: |0000>. Apply H_2: becomes: |0>0, |0>_1, (|0>+|1>)/√2 on qubit 2, |0>_3. So state: 1/√2 (|0000> + |0010>) but careful: qubits: 0,1,2,3. So state: |0>_0|0>_1 (|0>_2+|1>_2)/√2 |0>_3 = 1/√2 (|0000> + |0010>). Then apply CNOT{21}. This gate: control qubit 2, target qubit 1. So it flips qubit 1 if qubit 2 is 1. So after CNOT_{21}: state becomes: 1/√2 (|0000> + |0110>) because: if state is |0010>, then qubit2=1 so flip qubit1: from 0 to 1, so becomes |0110>. But wait: Actually, careful: The qubits order: 0,1,2,3. So state: |q0,q1,q2,q3>. So initial: |0,0,0,0> and |0,0,1,0>. After CNOT_{21}: For |0,0,0,0>: qubit2=0 so no flip: remains |0,0,0,0>. For |0,0,1,0>: qubit2=1 so flip qubit1: becomes |0,1,1,0>. So state: 1/√2 (|0000> + |0110>).

Then after CNOT_{21}, we have an error channel. But careful: The error channel acts on the two qubits involved: qubits 2 and 1. So it applies a Pauli error from the set of 15 non-identity Paulis on qubits (2,1) with probability p/15 each, and with probability 1-p, identity. So after this step, the state becomes a mixture. But we want the fidelity of the final physical state. The fidelity is defined as the overlap between the ideal state (which is the state without errors) and the actual state. But the actual state is a mixed state due to errors. So we want F = <ψ_ideal| ρ_actual |ψ_ideal>. And we want it as a function of p.

Then next: Apply H_0. That is a perfect gate. So after H_0: Ideal: H_0 on qubit 0. So state becomes: 1/√2 ( (H|0>)_0 |0>_1 |0>_2 |0>_3 + (H|0>)_0 |1>_1 |1>_2 |0>_3 ) = 1/√2 ( (|0>+|1>)/√2 |000> + (|0>+|1>)/√2 |110> ) = 1/2 ( |0000> + |1000> + |0110> + |1110> )? Wait careful: Actually: First term: |0000> becomes: H on qubit0: (|0>+|1>)/√2 on qubit0, so becomes: 1/√2 (|0000> + |1000>). Second term: |0110> becomes: H on qubit0: (|0>+|1>)/√2 on qubit0, so becomes: 1/√2 (|0110> + |1110>). So overall: 1/√2 * 1/√2 (|0000>+|1000>+|0110>+|1110>) = 1/2 (|0000>+|1000>+|0110>+|1110>).

But wait: That is the state after H_0 if no errors occurred. But we have errors from previous step. But we will compute the overall effect later.

Then after H_0, we apply CNOT_{03}. This gate: control qubit 0, target qubit 3. So ideal: State: 1/2 (|0000>+|1000>+|0110>+|1110>). After CNOT_{03}: For terms with qubit0=0: no flip on qubit3: so |0000> remains |0000>, |0110> remains |0110>. For terms with qubit0=1: flip qubit3: so |1000> becomes |1001>, |1110> becomes |1111>. So ideal final state: 1/2 (|0000> + |1001> + |0110> + |1111>). But is that the logical GHZ state? The logical GHZ state is (|00>{AB}+|11>{AB})/√2. And we have encoding: The logical operators: X_A = XIXI, X_B=XXII, Z_A=ZZII, Z_B=ZIZI. So what are the logical basis states? For the [[4,2,2]] code, the code space is 4-dimensional. The logical states are typically: |00>_L, |01>_L, |10>_L, |11>_L. But here we want GHZ: (|00>+|11>)/√2. Our final state: 1/2 (|0000> + |1001> + |0110> + |1111>). Check: Is that in the code space? The stabilizers are XXXX and ZZZZ. Check: XXXX applied to this state: XXXX|0000> = |1111> so that gives coefficient? Actually, the state is a superposition. For it to be in the code space, it must be an eigenvector of XXXX with eigenvalue +1. Check: XXXX on state: XXXX (|0000>+|1001>+|0110>+|1111>) = |1111> + |0110> + |1001> + |0000> = same state. So it is +1 eigenvector. ZZZZ: ZZZZ on state: ZZZZ|0000> = |0000> (since all zeros, phase +1) ZZZZ|1001> = (-1)^(1+0+0+1?) Actually careful: Z on each qubit: qubit0: if 1 then -1, so |1001>: qubit0=1 -> factor -1, qubit1=0 -> +1, qubit2=0 -> +1, qubit3=1 -> -1, so product = (-1)(-1)=+1? Actually: (-1)^(1+1)= (-1)^2=1. So |1001> gives +1. |0110>: qubit0=0, qubit1=1 -> -1, qubit2=1 -> -1, qubit3=0 -> +1, so product = (-1)(-1)=+1. |1111>: all ones: product = (-1)^4=1. So indeed it is +1 eigenvector of ZZZZ. So it is in the code space. Now, what are the logical states? We have logical operators: X_A = XIXI. Apply X_A to this state: X_A: flip qubit0 and qubit2? Actually: X_A = X on qubit0? Actually: X_A = XIXI. That means: on qubit0: X, qubit1: I, qubit2: X, qubit3: I. So X_A |0000> = |1010> which is not in our state. But wait, the logical operators are defined up to stabilizers. So we might represent logical states differently. Actually, the state we got: 1/2(|0000>+|1001>+|0110>+|1111>) might be the logical GHZ state. Check: The logical GHZ state is (|00>_L + |11>_L)/√2. So if we identify |00>_L and |11>_L appropriately. For instance, maybe |00>_L = (|0000>+|1111>)/√2 and |11>_L = (|1001>+|0110>)/√2? Then our state becomes: 1/2(|0000>+|1111>+|1001>+|0110>) = 1/√2 ( (|0000>+|1111>)/√2 + (|1001>+|0110>)/√2 ) = 1/√2 (|00>_L + |11>_L). So that works. So indeed the circuit prepares the logical GHZ state.

Now, we want the physical state fidelity. That is: F = <ψ_ideal| ρ_actual |ψ_ideal>, where ρ_actual is the state after the circuit including errors on CNOT gates.

We have two CNOT gates: CNOT_{21} and CNOT_{03}. Each CNOT gate is followed by an error channel. The error channel: After the CNOT gate, with probability 1-p, no error; with probability p, an error from the set of 15 two-qubit Paulis on the two qubits involved. And these errors are applied after the gate. So the overall evolution: We start with |0000>. Then apply H_2 (perfect). Then apply CNOT_{21} followed by error channel E1 on qubits (2,1). Then apply H_0 (perfect). Then apply CNOT_{03} followed by error channel E2 on qubits (0,3).

So the actual state is a mixture over error patterns. We want the fidelity relative to the ideal state (which is the state if no errors occur).

Let’s denote: Step 0: |ψ0> = |0000>. Step 1: Apply H_2: |ψ1> = (|0000> + |0010>)/√2. Step 2: Apply CNOT_{21} ideally: |ψ2> = (|0000> + |0110>)/√2. But then error channel on qubits (2,1): It applies a Pauli operator from the set {I, and 15 errors} with probabilities: P(I) = 1-p, and for each non-identity Pauli E (on two qubits), probability p/15. So after step 2, the state becomes: ρ_2 = (1-p) |ψ2><ψ2| + (p/15) * sum_{E in Paulis on qubits (2,1), E ≠ I} E |ψ2><ψ2| E†. But careful: The error is applied after the CNOT gate. But note: The CNOT gate itself is unitary. But the error channel acts on the state after the CNOT. So indeed, if we denote U = CNOT_{21}, then the operation is: U, then error channel on qubits (2,1). So the state after step 2 is: ρ = sum_{E} p_E (E U) |ψ1><ψ1| (U† E†) where p_I = 1-p, and p_E = p/15 for E ≠ I. But we can compute the effect on the fidelity later. Alternatively, we can compute the overall transformation as: The final ideal state is: |ψ_ideal> = U_total |0000>, where U_total = CNOT_{03} * H_0 * CNOT_{21} * H_2. And with errors: Actually, careful: The errors occur after each CNOT. So the actual operation is: Actual = (E2 * CNOT_{03}) * H_0 * (E1 * CNOT_{21}) * H_2, where E1 is a random Pauli on qubits (2,1) and E2 is a random Pauli on qubits (0,3). And these are independent. So the final state is: ρ = average_{E1, E2} (U_E1,E2) |0000><0000| (U_E1,E2)†, where U_E1,E2 = (E2 * CNOT_{03}) * H_0 * (E1 * CNOT_{21}) * H_2. And we want F = <ψ_ideal| ρ |ψ_ideal> = average_{E1,E2} |<ψ_ideal| U_E1,E2 |0000>|^2.

So we need to compute for each error pattern (E1 on qubits (2,1) and E2 on qubits (0,3)), the overlap squared between the ideal state and the state produced by that error pattern. Then average over the probabilities: For E1: probability: if E1 = I, probability = 1-p; if E1 is any of the 15 non-identity, probability = p/15 each. Similarly for E2.

So F = (1-p)^2 * (overlap for both no errors) + (1-p)(p/15) (sum over E2 errors? Actually careful: There are two independent error channels. So overall: F = Σ_{E1 in Paulis on {1,2}} Σ_{E2 in Paulis on {0,3}} (Prob(E1)* (Prob(E2)) * |<ψ_ideal| U(E1,E2) |0000>|^2, where U(E1,E2) = CNOT_{03} * E2? Actually careful: The order: U(E1,E2) = (CNOT_{03} * (error E2 on (0,3))? But wait: The error E2 is applied after CNOT_{03}. But then it is: U(E1,E2) = (E2 * CNOT_{03}) * H_0 * (E1 * CNOT_{21}) * H_2. But note that CNOT_{03} and H_0 commute? Not necessarily. But we can compute step by step.

Maybe it’s easier to compute the effect of errors on the state step by step. But since the errors are Pauli errors, and the ideal state is a specific state, we can compute how each Pauli error affects the overlap.

Alternate approach: The fidelity F is the probability that no error occurs that takes the state out of the code space? But careful: Even if an error occurs, it might still leave the state in the ideal state if the error is a stabilizer? But wait: The ideal state is a specific state. If an error E is such that E|ψ_ideal> is proportional to |ψ_ideal>, then the overlap is 1. But if it is orthogonal, then overlap 0. But also there might be cases where the overlap is something else? But since the errors are Pauli operators, and the ideal state is a specific state in the code space, then for a Pauli error E (which acts on the physical qubits), either E|ψ_ideal> is either exactly |ψ_ideal> (if E is a stabilizer or logical operator that leaves that state invariant?) Actually, careful: The ideal state is not necessarily an eigenstate of all logical operators. For example, the GHZ state: (|00>_L+|11>_L)/√2. That state is not an eigenstate of logical X on one qubit? Actually, it might be. But anyway, we can compute the overlap explicitly.

Maybe we can compute the transformation step by step and then compute the effect of errors.

Let’s denote the ideal transformation without errors: U_ideal = CNOT_{03} * H_0 * CNOT_{21} * H_2. We want to compute for a given error pattern (E1 on qubits (2,1) and E2 on qubits (0,3)), what is U(E1,E2) |0000>? And then compute its overlap with |ψ_ideal> = U_ideal|0000>.

Maybe we can compute the state after each step in terms of basis states. But since errors are Pauli, they just flip bits and phases. Perhaps we can compute the effect on the final state.

Alternatively, note that the circuit is relatively small. We can compute the state after each step explicitly.

Let’s denote the qubits as: q0, q1, q2, q3. Initial: |0000>.

Step 1: H_2: becomes: 1/√2 (|0000> + |0010>).

Step 2: CNOT_{21} (control 2, target 1): |0000> -> |0000> |0010> -> |0110> So state: 1/√2 (|0000> + |0110>).

Now, error E1 on qubits (2,1). E1 is a Pauli operator on qubits 2 and 1. It can be written as: E1 = σ_a ⊗ σ_b, where σ_a, σ_b are from {I, X, Y, Z} but not both I? Actually, includes I? But error channel: with probability 1-p, we take I. And with probability p/15 each, we take any of the 15 non-identity Paulis. So we consider E1 as any Pauli on two qubits. But we’ll compute the effect on the state. The state after step 2 (before error) is: |ψ2> = 1/√2 (|0000> + |0110>). Now, apply a Pauli operator on qubits 1 and 2. Let’s denote the Pauli as: P = P_1 ⊗ P_2, acting on qubit1 and qubit2. But careful: The order: qubits: 0,1,2,3. So error on qubits (2,1) means: it acts on qubit2 and qubit1. But note: The order in the tensor product: usually we list qubits in order: 0,1,2,3. So if we have an operator on qubits 1 and 2, it is I_0 ⊗ (operator on qubit1) ⊗ (operator on qubit2) ⊗ I_3. But careful: The error E1 is applied after CNOT_{21}. So the state becomes: E1 |ψ2>. So let’s compute: For each Pauli on two qubits, what is E1 |ψ2>? But note: The state |ψ2> is a superposition of two basis states: |0,0,0,0> and |0,1,1,0>. So we can compute the effect of a Pauli on qubits 1 and 2 on these states. Let’s denote the basis states as |q0, q1, q2, q3>. For the two states: State A: |0,0,0,0> State B: |0,1,1,0>

A Pauli on qubits 1 and 2: It will map these to some states with possible phases. For example, if E1 = I ⊗ I, then it’s identity. If E1 = I ⊗ X (on qubit2 only? But careful: The Pauli on two qubits: They are products of Paulis on each qubit. So we list all 16 Paulis? But we only need the overlap with the ideal final state later. Perhaps we can compute the final state after the entire circuit with an error E1. That might be easier: The circuit after step 2 error is: Then apply H_0 (on qubit0) and then CNOT_{03} with error E2. So the final state (for a given error pattern) is: |ψ_final(E1,E2)> = (CNOT_{03} * E2) * H_0 * (E1 * CNOT_{21} * H_2) |0000>? But careful: The order: Actually, the operations are: H_2, then CNOT_{21}, then E1, then H_0, then CNOT_{03}, then E2. So if we denote U1 = CNOT_{21} * H_2, then state after step2 (before error) is U1|0000>. Then after E1, it’s E1 U1|0000>. Then H_0 gives H_0 E1 U1|0000>. Then CNOT_{03} gives CNOT_{03} H_0 E1 U1|0000>. Then E2 gives E2 CNOT_{03} H_0 E1 U1|0000>. So final state: |ψ_final> = E2 * CNOT_{03} * H_0 * E1 * CNOT_{21} * H_2 |0000>. And the ideal is: |ψ_ideal> = CNOT_{03} * H_0 * CNOT_{21} * H_2 |0000>. So if we define V = CNOT_{03} * H_0 * CNOT_{21} * H_2, then ideal = V|0000>. And with errors: |ψ_final> = E2 * CNOT_{03} * H_0 * E1 * (CNOT_{21} * H_2) |0000>. But note that CNOT_{21} * H_2 is not necessarily commuting with E1? But we can try to move E1 through H_0 and CNOT_{03}? Because E1 acts on qubits 1 and 2. And H_0 acts on qubit0, so it commutes with E1. And CNOT_{03} acts on qubits 0 and 3, so it commutes with E1 as well because they act on disjoint sets? Actually, careful: E1 acts on qubits 1 and 2. CNOT_{03} acts on qubits 0 and 3. So they commute. So we have: |ψ_final> = E2 * CNOT_{03} * H_0 * E1 * CNOT_{21} * H_2 |0000> = E2 * E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2 |0000>? Check: E1 commutes with H_0 and CNOT_{03} because they act on different qubits. But careful: E1 commutes with CNOT_{21}? Actually, E1 is applied after CNOT_{21} so it’s not that we can move E1 before CNOT_{21} because they don’t necessarily commute? But wait: The order is: … * H_0 * E1 * CNOT_{21} * H_2. But E1 and CNOT_{21} act on the same qubits? Actually, E1 acts on qubits 1 and 2. CNOT_{21} acts on qubits 2 and 1. So they do not commute in general. So we cannot simply move E1 past CNOT_{21} because that would change the order. But we can move E1 past H_0 and CNOT_{03} because those act on qubits 0 and 3, so they commute with E1. So we get: |ψ_final> = E2 * CNOT_{03} * H_0 * E1 * CNOT_{21} * H_2 |0000> = E2 * E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2 |0000>? That would be if E1 commutes with CNOT_{21}? But they don’t necessarily commute. Actually, check: The operations: We have: H_0 * E1 * CNOT_{21}. Since H_0 commutes with E1 (different qubits), we can write H_0 * E1 = E1 * H_0. So then it becomes: E2 * CNOT_{03} * E1 * H_0 * CNOT_{21} * H_2. And then, does E1 commute with CNOT_{03}? Yes, because CNOT_{03} acts on qubits 0 and 3, and E1 on qubits 1 and 2, so they commute. So then we get: E2 * E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2 = E2 * E1 * V. So actually, it seems that E1 does commute with CNOT_{21}? Wait, careful: The order is: … * H_0 * E1 * CNOT_{21} * H_2. We want to move E1 to the left. But E1 and CNOT_{21} do not necessarily commute because they act on the same qubits. However, note that H_0 is to the left of E1. But H_0 commutes with E1, so we can swap H_0 and E1: That gives: … * E1 * H_0 * CNOT_{21} * H_2. So then we have E1 and then H_0 * CNOT_{21} * H_2. But then can we move E1 past CNOT_{03}? CNOT_{03} is to the left of that whole block? Actually, the full expression is: E2 * CNOT_{03} * (H_0 * E1 * CNOT_{21} * H_2). And we already swapped H_0 and E1: = E2 * CNOT_{03} * (E1 * H_0 * CNOT_{21} * H_2) = E2 * CNOT_{03} * E1 * (H_0 * CNOT_{21} * H_2). And since CNOT_{03} and E1 commute (disjoint sets), we get: E2 * E1 * CNOT_{03} * (H_0 * CNOT_{21} * H_2) = E2 * E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2 = E2 * E1 * V. So indeed, it seems that we can move E1 all the way to the left? But wait: Check: The original order: … * H_0 * E1 * CNOT_{21} * H_2. Is it valid to swap E1 and CNOT_{21}? We did not swap E1 and CNOT_{21}. We only swapped H_0 and E1. That is valid because H_0 acts on qubit0 and E1 on qubits 1,2. Then we have CNOT_{03} * E1. Those commute because CNOT_{03} acts on qubits 0,3 and E1 on 1,2. So then we get E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2. But then we have E1 and then CNOT_{03} * H_0 * CNOT_{21} * H_2. But note: CNOT_{03} * H_0 * CNOT_{21} * H_2 is V. But wait, is it? V = CNOT_{03} * H_0 * CNOT_{21} * H_2. So then we get: |ψ_final> = E2 * E1 * V |0000>? That would imply that the error E1 effectively acts on the final state? But careful: That manipulation would be valid if E1 commuted with CNOT_{21}? But we never swapped E1 and CNOT_{21}. Let’s do it step by step carefully:

We have: U_total with errors = E2 * CNOT_{03} * H_0 * E1 * CNOT_{21} * H_2. Write this as: U = E2 * CNOT_{03} * H_0 * E1 * CNOT_{21} * H_2. Consider the product: H_0 * E1. Since H_0 acts on qubit0 and E1 on qubits 1,2, they commute. So H_0 * E1 = E1 * H_0. So then: U = E2 * CNOT_{03} * (E1 * H_0) * CNOT_{21} * H_2 = E2 * CNOT_{03} * E1 * H_0 * CNOT_{21} * H_2. Now, CNOT_{03} and E1: They act on disjoint sets (CNOT_{03} on qubits 0,3; E1 on qubits 1,2), so they commute. So CNOT_{03} * E1 = E1 * CNOT_{03}. So then: U = E2 * (E1 * CNOT_{03}) * H_0 * CNOT_{21} * H_2 = E2 * E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2 = E2 * E1 * V. So indeed, U = E2 * E1 * V. That is surprising because it suggests that the error E1, which is applied after CNOT_{21}, actually becomes like an error on the final state? But wait, is that valid? Check with a simple example: Suppose we have a circuit: CNOT followed by an error E on the same qubits. That is: E * CNOT. And if we had instead applied the error before the CNOT, it would be CNOT * E’ where E’ is something else. They are not necessarily the same. But here, we are not moving E1 before CNOT_{21}; we are moving it past H_0 and CNOT_{03}. And indeed, H_0 and CNOT_{03} act on different qubits, so they commute with E1. So it is valid: U = E2 * CNOT_{03} * H_0 * E1 * CNOT_{21} * H_2. = E2 * CNOT_{03} * (H_0 * E1) * CNOT_{21} * H_2 = E2 * CNOT_{03} * (E1 * H_0) * CNOT_{21} * H_2 = E2 * (CNOT_{03} * E1) * H_0 * CNOT_{21} * H_2 = E2 * (E1 * CNOT_{03}) * H_0 * CNOT_{21} * H_2 = E2 * E1 * CNOT_{03} * H_0 * CNOT_{21} * H_2 = E2 * E1 * V. So indeed, the effect of the error E1 (on qubits 1,2) applied after CNOT_{21} is equivalent to applying E1 at the very end? But careful: That would be if the operations after E1 commute with E1. And they do because H_0 acts on qubit0 and CNOT_{03} acts on qubits 0,3. So yes, that is valid. So then the final state with errors is: |ψ_final> = (E2 * E1) V |0000> = (E2 * E1) |ψ_ideal>. But wait: Is that always true? Check with a specific error: Say E1 = X on qubit1 (so Pauli X on qubit1, and I on qubit2). Then originally: After CNOT_{21}, state is 1/√2 (|0000> + |0110>). Then E1 (X on qubit1) gives: 1/√2 (|0100> + |0010>). Then H_0: H on qubit0: becomes: 1/√2 ( (|0>+|1>)/√2 on qubit0 for first term? Actually, careful: State after E1: |0100> and |0010>. Apply H_0: |0100> -> (|0100> + |1100>)/√2. |0010> -> (|0010> + |1010>)/√2. So state = 1/√2 * 1/√2 (|0100>+|1100>+|0010>+|1010>) = 1/2 (|0100>+|1100>+|0010>+|1010>). Then CNOT_{03}: |0100>: qubit0=0 -> no flip: remains |0100> |1100>: qubit0=1 -> flip qubit3: becomes |1101> |0010>: qubit0=0 -> remains |0010> |1010>: qubit0=1 -> becomes |1011> So state = 1/2 (|0100>+|1101>+|0010>+|1011>). Now, if we instead apply E1 at the end on the ideal state? The ideal state is: |ψ_ideal> = 1/2 (|0000>+|1001>+|0110>+|1111>). Applying E1 (X on qubit1) to that: X on qubit1: |0000> -> |0100> |1001> -> |1101> |0110> -> |0010> (because qubit1: 1->0) |1111> -> |1011> So indeed we get the same state. So it works.

Similarly, E2 is applied after CNOT_{03} but since there is nothing after, it’s just E2 at the end. So indeed, the overall effect is that the final state is: (E2 * E1) |ψ_ideal>, where E1 is the error on qubits (2,1) and E2 is the error on qubits (0,3). But wait: Is that always true? Check: What if E1 acts on qubit2? For example, E1 = X on qubit2. Then after CNOT_{21}, state: 1/√2 (|0000>+|0110>). Apply E1: X on qubit2: becomes: 1/√2 (|0010>+|0100>). Then H_0: |0010> -> (|0010>+|1010>)/√2 |0100> -> (|0100>+|1100>)/√2 So state = 1/√2 * 1/√2 (|0010>+|1010>+|0100>+|1100>) = 1/2 (|0010>+|1010>+|0100>+|1100>) Then CNOT_{03}: |0010>: qubit0=0 -> remains |0010> |1010>: qubit0=1 -> becomes |1011> |0100>: qubit0=0 -> remains |0100> |1100>: qubit0=1 -> becomes |1101> So state = 1/2 (|0010>+|1011>+|0100>+|1101>). Now, if we apply E1 at the end to ideal state: E1 = X on qubit2. Ideal: 1/2 (|0000>+|1001>+|0110>+|1111>) Apply X on qubit2: |0000> -> |0010> |1001> -> |1011> |0110> -> |0100> (qubit2: 1->0) |1111> -> |1101> So indeed same. So it works.

Thus, the final state with errors is: |ψ_final(E1,E2)> = (E2 E1) |ψ_ideal>, where E1 is a Pauli on qubits {1,2} and E2 is a Pauli on qubits {0,3}. And these errors occur with probabilities: For E1: P(I) = 1-p, and for each of the 15 non-identity, p/15. Similarly for E2.

Then the fidelity is: F = average_{E1,E2} |<ψ_ideal| (E2 E1) |ψ_ideal>|^2. But careful: The errors E1 and E2 are independent. So: F = Σ_{E1} Σ_{E2} P(E1) P(E2) |<ψ_ideal| E2 E1 |ψ_ideal>|^2. And note that E2 and E1 act on disjoint sets, so they commute. So E2 E1 is a Pauli operator on all 4 qubits. But it’s not arbitrary: It is of the form: (Pauli on qubits 0,3) times (Pauli on qubits 1,2). So there are 16*16 = 256 such operators. But many will give the same overlap.

So we need to compute: For a given Pauli operator P (on 4 qubits) that factors as P = P_{03} ⊗ P_{12} (with the order: qubit0,3 and qubit1,2 separately), what is |<ψ_ideal| P |ψ_ideal>|^2? And then average over the distribution of P. But careful: The distribution: E1 is on qubits 1,2: probability: For each Pauli on two qubits, probability: if it is identity, then probability = 1-p; if non-identity, then probability = p/15. Similarly, E2 on qubits 0,3: same. And then the overall Pauli is P = E2 * E1. But note: This product is not independent over all 256 because the probabilities factor: P(P) = (probability for E2) * (probability for E1). So we can compute: F = (Σ_{E2 on {0,3}} P(E2) ) * (Σ_{E1 on {1,2}} P(E1) |<ψ_ideal| E2 E1 |ψ_ideal>|^2 )? But careful: The overlap |<ψ_ideal| E2 E1 |ψ_ideal>|^2 depends on both E2 and E1. So we need to sum over both.

Maybe we can compute: For a given Pauli operator P on all 4 qubits that is of the form: P = P_0 ⊗ P_1 ⊗ P_2 ⊗ P_3, but with the constraint that P_0 and P_3 come from E2 (so they are correlated? Actually, E2 is a Pauli on qubits 0 and 3. That means it is of the form: σ_a ⊗ σ_b on qubits 0 and 3. But note: There is no constraint that the Pauli on qubit0 and qubit3 are the same type? They can be independent. So E2 can be, for example, X on qubit0 and Z on qubit3. So indeed, E2 is any element of the Pauli group on two qubits. Similarly, E1 is any element on qubits 1 and 2. So then the overall Pauli is: P = (E2) ⊗ (E1) but careful: The order: qubits: 0,1,2,3. So if E2 = σ_{a} on qubit0 and σ_{b} on qubit3, and E1 = σ_{c} on qubit1 and σ_{d} on qubit2, then P = σ_a ⊗ σ_c ⊗ σ_d ⊗ σ_b. So it’s a Pauli on 4 qubits. So we want to compute the overlap squared: |<ψ_ideal| P |ψ_ideal>|^2.

Now, what is |ψ_ideal>? |ψ_ideal> = 1/2 (|0000> + |1001> + |0110> + |1111>). Let’s write it explicitly: |ψ_ideal> = 1/2 ( |0,0,0,0> + |1,0,0,1> + |0,1,1,0> + |1,1,1,1> ).

We want to compute <ψ_ideal| P |ψ_ideal> for a Pauli P. Since P is a Pauli, it will map basis states to other basis states up to a phase. So we can compute: <ψ_ideal| P |ψ_ideal> = 1/4 sum_{x,y in {the four basis states}} . But maybe we can determine conditions for P to leave |ψ_ideal> invariant (up to a phase). Because then the overlap squared is 1. Otherwise, it might be 0? But careful: Could it be something like 1/2? Let’s check: Suppose P = I on all qubits, then overlap = 1. What if P is a stabilizer? The stabilizers are XXXX and ZZZZ. Check: XXXX|ψ_ideal> = ? We computed: XXXX|ψ_ideal> = |ψ_ideal> so that gives overlap = 1. What about ZZZZ? ZZZZ|ψ_ideal> = |ψ_ideal> as computed, so overlap = 1. What about something like: P = X on qubit0 only? That is not of the form we consider? But wait, our errors are only on the pairs: E2 on (0,3) and E1 on (1,2). So P will always have the property that the Pauli on qubit0 and qubit3 are from the same two-qubit Pauli? Actually, careful: E2 is a Pauli on qubits 0 and 3. So if E2 = X ⊗ I, then that means on qubit0: X, on qubit3: I. And E1 could be I ⊗ I. Then overall P = X on qubit0, I on others. So such an operator is possible. So let’s compute: P = X_0 (i.e., Pauli X on qubit0 only). Then P|ψ_ideal> = 1/2 (X|0000> + X|1001> + X|0110> + X|1111>) = 1/2 (|1000> + |0001> + |1110> + |0111>). The overlap with |ψ_ideal> is: <ψ_ideal| = 1/2 (<0000|+<1001|+<0110|+<1111|). Dot product: <0000|1000> = 0, <0000|0001> = 0, <0000|1110> = 0, <0000|0111> = 0. <1001|1000> = 0, <1001|0001> = 0? Actually, <1001|0001> = 0 because first qubit: 1 vs 0, <1001|1110> = 0, <1001|0111> = 0. <0110|1000> = 0, <0110|0001> = 0, <0110|1110> = 0, <0110|0111> = 0, <1111|1000> = 0, <1111|0001> = 0, <1111|1110> = 0, <1111|0111> = 0. So overlap = 0. So |<ψ_ideal| X_0 |ψ_ideal>|^2 = 0. What about P = Z on qubit0? Then P|ψ_ideal> = 1/2 (Z|0000> + Z|1001> + Z|0110> + Z|1111>) = 1/2 (|0000> - |1001> + |0110> - |1111>). Overlap with |ψ_ideal> = 1/2 (1 -1 +1 -1)/2? Actually, compute: <ψ_ideal| P|ψ_ideal> = 1/4 ( <0000|0000> + <1001|1001>(-1) + <0110|0110> + <1111|1111>(-1) )? But careful: <ψ_ideal| = 1/2 (<0000|+<1001|+<0110|+<1111|) So <ψ_ideal| P|ψ_ideal> = 1/4 [ <0000|0000> + <1001|(-|1001>) + <0110|0110> + <1111|(-|1111>) ] = 1/4 (1 -1+1-1)=0. So that gives 0. What about P = X on qubit0 and X on qubit3? That is E2 = X ⊗ X. And E1 = I. Then P = X_0 X_3. Then P|ψ_ideal> = 1/2 (X_0X_3|0000> + X_0X_3|1001> + X_0X_3|0110> + X_0X_3|1111>) = |0000> -> X_0X_3|0000> = |1001>? Actually, careful: X on qubit0 and X on qubit3: |0000> -> |1001> |1001> -> |0000> |0110> -> |1110>? Wait: |0110>: qubit0=0 -> becomes 1, qubit3=0 -> becomes 1, so becomes |1110>? But then qubit1 and qubit2 remain: 1 and 1, so actually |1110>. |1111> -> |0111>? Actually, |1111>: qubit0=1 -> becomes 0, qubit3=1 -> becomes 0, so becomes |0111>. So P|ψ_ideal> = 1/2 (|1001> + |0000> + |1110> + |0111>). That is not exactly |ψ_ideal> because |ψ_ideal> has |0110> and |1111> instead of |1110> and |0111>. So overlap? <ψ_ideal| = 1/2 (<0000|+<1001|+<0110|+<1111|). Dot product: <0000| from P: coefficient 1/2 from second term gives 1/2 * 1/2? Actually, careful: <ψ_ideal| P|ψ_ideal> = (1/2)(1/2) sum_{i,j} but maybe compute directly: <ψ_ideal| P|ψ_ideal> = 1/4 [ <0000|P|0000> + <0000|P|1001> + … Actually, easier: The states in |ψ_ideal> are: A = |0000>, B = |1001>, C = |0110>, D = |1111>. P maps: A -> |1001> = B, B -> |0000> = A, C -> |1110> (call this E, not in set), D -> |0111> (F, not in set). So then <ψ_ideal| P|ψ_ideal> = 1/4 ( + ? Actually, careful: P|ψ_ideal> = 1/2 (B + A + E + F). Then <ψ_ideal| = 1/2 (A†+B†+C†+D†). So inner product = 1/4 (A†(B+A+E+F) + … = 1/4 (0+1+0+0 + from B: 1+0+0+0 + from C: 0+0+0+0 + from D: 0+0+0+0) = 1/4*(1+1)=1/2. So then |<ψ_ideal| P|ψ_ideal>|^2 = 1/4. So that is not 0 or 1. So we need to compute the overlap squared for each Pauli P that is of the form: P = (Pauli on qubits 0,3) ⊗ (Pauli on qubits 1,2). And then average over the probabilities.

So let’s denote: P = P_{03} ⊗ P_{12}, where P_{03} is a Pauli on qubits 0 and 3, and P_{12} is a Pauli on qubits 1 and 2. We want: O(P) = |<ψ_ideal| P |ψ_ideal>|^2. And |ψ_ideal> = 1/2 (|0000> + |1001> + |0110> + |1111>).

Let’s denote basis states as: A = |0,0,0,0> B = |1,0,0,1> C = |0,1,1,0> D = |1,1,1,1>

So |ψ_ideal> = (A+B+C+D)/2.

Now, P acts on these qubits. Write P = P_0 ⊗ P_1 ⊗ P_2 ⊗ P_3. But note that our P is not arbitrary: It factors as (P_0 ⊗ P_3) and (P_1 ⊗ P_2) independently. But that is actually all Paulis? Actually, any Pauli on 4 qubits can be written as a product of Paulis on each qubit. But here the probability distribution is such that the Pauli on qubits 0 and 3 comes from the set of 16 with a specific distribution (1-p for identity, p/15 for each non-identity) and similarly for qubits 1 and 2. So we need to compute O(P) for each such P.

Since P is a Pauli, it will map each computational basis state to another computational basis state, multiplied by a phase ±1 or ±i. But Paulis are Hermitian and square to identity, so they have eigenvalues ±1. Actually, careful: Pauli matrices are Hermitian, so they have real eigenvalues. But Y gives i? Actually, Y = iσ_y? Actually, Pauli matrices are usually defined as X, Y, Z where Y = iσ_y? Actually, in quantum computing, Pauli matrices are often taken as X, Y, Z with Y = [[0, -i], [i, 0]] so that they are Hermitian. So indeed, X, Y, Z are Hermitian. So phases are ±1. But wait: Y = iσ_y? Actually, standard Pauli matrices: σ_x, σ_y, σ_z. σ_y is [[0, -i], [i, 0]]. That is Hermitian? Check: (σ_y)† = transpose conjugate = [[0, -i], [i, 0]]^T = [[0, i], [-i, 0]] which is not the same. Actually, σ_y is not Hermitian? Wait, Hermitian means A = A†. For σ_y, we have σ_y = [[0, -i], [i, 0]]. Its conjugate transpose is: take conjugate: [[0, i], [-i, 0]], then transpose: [[0, -i], [i, 0]] which is actually the same? Let’s do carefully: σ_y = [ [0, -i], [i, 0] ]. Conjugate: [ [0, i], [-i, 0] ]. Transpose: [ [0, -i], [i, 0] ]. So indeed, σ_y is Hermitian. So all Paulis are Hermitian and have eigenvalues ±1. So the phase factors are ±1. So when P acts on a basis state, it gives ± another basis state.

So we want: <ψ_ideal| P |ψ_ideal> = 1/4 sum_{s,t in {A,B,C,D}} . But since P is a Pauli, it will map each basis state to a specific basis state (maybe not necessarily one of A,B,C,D? But note: A,B,C,D have specific patterns: A: 0,0,0,0 B: 1,0,0,1 C: 0,1,1,0 D: 1,1,1,1 So these states are characterized by: qubit0 and qubit3 are the same? Actually, A: 0 and 0; B: 1 and 1; C: 0 and 0? Actually, C: qubit0=0, qubit3=0; D: qubit0=1, qubit3=1. And qubit1 and qubit2: A: 0,0; B: 0,0; C: 1,1; D: 1,1. So the states in the support of |ψ_ideal> are those where qubit0 = qubit3 and qubit1 = qubit2. So if P maps one of these states to another state, it will be in the support if the image also satisfies that condition. So we need to determine for a given Pauli P (which is a product of Paulis on each qubit), what is the action on the set {A,B,C,D}? And then compute the matrix element.

Alternatively, we can compute: <ψ_ideal| P |ψ_ideal> = 1/4 ( + + … + ) but note that P will map each basis state to a basis state. So many terms are zero unless P maps exactly one of these to another. So it’s easier: Let f(P) = <ψ_ideal| P |ψ_ideal>. This is the expectation value of P in the state |ψ_ideal>. And since |ψ_ideal> is a pure state, f(P) is a number. And then we want |f(P)|^2.

Maybe we can compute f(P) by noting that |ψ_ideal> is an eigenstate of the stabilizers XXXX and ZZZZ. Also, it is a GHZ state in the logical basis. Perhaps we can express P in terms of logical operators? But careful: P is a physical Pauli. And the code is [[4,2,2]]. The logical operators are: X_A = XIXI, X_B = XXII, Z_A = ZZII, Z_B = ZIZI. But our state is not a logical basis state but a GHZ state. Alternatively, we can compute f(P) by direct computation. Since P factors as (P on qubits 0,3) and (P on qubits 1,2), and the state |ψ_ideal> is a product state between the pair (0,3) and (1,2)? Actually, check: |ψ_ideal> = 1/2 (|00>{03}|00>{12} + |11>{03}|00>{12}? That doesn’t work because: A: |00>{03}|00>{12}; B: |11>{03}|00>{12}; C: |00>{03}|11>{12}; D: |11>{03}|11>{12}. So indeed, |ψ_ideal> = (|00>{03} + |11>{03})/√2 ⊗ (|00>{12} + |11>{12})/√2? But careful: That would be 1/2 (|00>|00> + |00>|11> + |11>|00> + |11>|11>) but our state is only the terms where the two pairs are either both 00 or both 11. So actually, |ψ_ideal> is not a product state between the two pairs; it is entangled across them? Actually, it is: |ψ_ideal> = 1/2 (|00>{03}|00>{12} + |11>{03}|00>{12} + |00>{03}|11>{12} + |11>{03}|11>{12}). That is actually a product state? Because (|00>+|11>)/√2 on first pair and (|00>+|11>)/√2 on second pair would give: 1/2 (|00>|00> + |00>|11> + |11>|00> + |11>|11>). So indeed, |ψ_ideal> = (|00>{03}+|11>{03})/√2 ⊗ (|00>{12}+|11>{12})/√2. Wait, check: (|00>+|11>)/√2 on qubits 0,3 gives: 1/√2 (|00> + |11>). And on qubits 1,2: 1/√2 (|00> + |11>). Their tensor product is 1/2 (|0000> + |0011> + |1100> + |1111>). But our |ψ_ideal> is: 1/2 (|0000> + |1001> + |0110> + |1111>). These are different: Our state has |1001> and |0110> instead of |0011> and |1100>. So it’s not a product across the partition (0,3) and (1,2). Actually, note: In our state, qubit0 and qubit3 are correlated: they are the same. And qubit1 and qubit2 are correlated: they are the same. But if we take the pair (0,3) and (1,2), then the state is: 1/2 (|00,00> + |11,00> + |00,11> + |11,11>) but careful: |11,00> means qubits 0,3 = 11 and qubits 1,2 = 00, that is |1,0,0,1> which is B. And |00,11> means qubits 0,3 = 00 and qubits 1,2 = 11, that is |0,1,1,0> which is C. So indeed, if we label the pair (0,3) as one system and (1,2) as the other, then the state is: 1/2 (|00,00> + |11,00> + |00,11> + |11,11>). That is actually a product state? Because a product state would be (a|00>+b|11>) ⊗ (c|00>+d|11>) = ac|00,00> + ad|00,11> + bc|11,00> + bd|11,11>. For this to equal 1/2 (|00,00>+|11,00>+|00,11>+|11,11>), we need ac=1/2, ad=1/2, bc=1/2, bd=1/2. That implies a, b, c, d are all nonzero and a/b =? Actually, from ad=1/2 and bc=1/2, we get a/b = (1/(2d))/(1/(2c)) = c/d. And from ac=1/2 and bd=1/2, we get a/b = (1/(2c))/(1/(2d)) = d/c. So we need c/d = d/c, so c^2 = d^2, so d = ± c. And then ac=1/2 gives a=1/(2c), and ad=1/2 gives a=1/(2d) so that forces c=d. So then a=1/(2c), and bc=1/2 gives b=1/(2c), so a=b. Then the product state is (a|00>+a|11>) ⊗ (c|00>+c|11>) = a c (|00>+|11>)⊗(|00>+|11>) = ac (|00,00>+|00,11>+|11,00>+|11,11>). For this to be 1/2, we need ac=1/2. So indeed, if we take a = 1/√2 and c = 1/√2, then we get the state. So actually, |ψ_ideal> is a product state between the two pairs? But careful: The pairs are (0,3) and (1,2). But then the state would be: (|00>{03}+|11>{03})/√2 ⊗ (|00>{12}+|11>{12})/√2 = 1/2 (|00,00>+|00,11>+|11,00>+|11,11>). But that gives: |00,00> = |0,0,0,0> = A, good. |00,11> = |0,1,1,0> = C, good. |11,00> = |1,0,0,1> = B, good. |11,11> = |1,1,1,1> = D, good. So indeed, |ψ_ideal> = (|00>{03}+|11>{03})/√2 ⊗ (|00>{12}+|11>{12})/√2. That is a product state! But wait, then the logical GHZ state is actually a product of two Bell states? That is interesting. So we have: |ψ_ideal> = |Φ^+>{03} ⊗ |Φ^+>{12}, where |Φ^+> = (|00>+|11>)/√2. That is a very nice factorization. So then the expectation value of a Pauli operator P that factors as P = P_{03} ⊗ P_{12} is: <ψ_ideal| P |ψ_ideal> = (<Φ^+|{03} ⊗ <Φ^+|{12}|) (P_{03} ⊗ P_{12}) (|Φ^+>{03} ⊗ |Φ^+>{12}) = (<Φ^+|{03}| P{03} |Φ^+>{03}) * (<Φ^+|{12}| P_{12} |Φ^+>{12}). So then |<ψ_ideal| P |ψ_ideal>|^2 = |<Φ^+| P{03} |Φ^+>|^2 * |<Φ^+| P_{12} |Φ^+>|^2.

And we know the expectation values for two-qubit Paulis on a Bell state |Φ^+> = (|00>+|11>)/√2. For a two-qubit Pauli operator Q on two qubits, what is <Φ^+| Q |Φ^+>? Recall: |Φ^+> = (|00>+|11>)/√2. For Pauli operators: I ⊗ I: <Φ^+| I|Φ^+> = 1. X ⊗ I: <Φ^+| (X⊗I) |Φ^+> = ? (X⊗I)|Φ^+> = (|10>+|01>)/√2 = |Ψ^+> actually. And <Φ^+|Ψ^+> = 0. So = 0. Similarly, I ⊗ X: 0. Z ⊗ I: (Z⊗I)|Φ^+> = (|00> - |11>)/√2 = |Φ^-⟩, so overlap 0. I ⊗ Z: similarly 0. X ⊗ X: (X⊗X)|Φ^+> = (|11>+|00>)/√2 = |Φ^+>, so <Φ^+| = 1. Y ⊗ Y: (Y⊗Y)|Φ^+> = ? Y|0> = i|1>, Y|1> = -i|0>, so (Y⊗Y)|00> = (i)(-i)|11> = |11>, and (Y⊗Y)|11> = (-i)(i)|00> = |00>, so actually (Y⊗Y)|Φ^+> = |Φ^+>, so <Φ^+| = 1? But careful: There might be a phase? Actually, Y|0> = i|1>, Y|1> = -i|0>, so (Y⊗Y)|00> = i * (-i) |11> = 1|11>, and (Y⊗Y)|11> = (-i)(i)|00> = 1|00>. So indeed, Y⊗Y gives 1. Z ⊗ Z: (Z⊗Z)|Φ^+> = (|00>+(-1)(-1)|11>)/√2? Actually, Z⊗Z|00> = |00>, Z⊗Z|11> = (-1)(-1)|11> = |11>, so = |Φ^+>, so <Φ^+| = 1. What about cross terms like X⊗Z? (X⊗Z)|Φ^+> = (|10> + (-1)|01>)/√2 = (|10> - |01>)/√2 = -|Ψ^->, so overlap with |Φ^+> is 0. Similarly, any Pauli that is not of the form where the two Paulis are the same? Actually, more generally: For |Φ^+>, we have: <Φ^+| (σ ⊗ τ) |Φ^+> = Tr(σ^T τ)/2? Actually, there is a formula: <Φ^+| (A⊗B) |Φ^+> = Tr(A^T B)/2. So then if A and B are Paulis, then A^T is just A (since Paulis are symmetric except Y? Actually, X and Z are symmetric, Y is antisymmetric? Y^T = -Y. So careful: For Paulis, X^T = X, Z^T = Z, Y^T = -Y. So then <Φ^+| (X⊗X) |Φ^+> = Tr(X X)/2 = Tr(I)/2 = 1. For (Y⊗Y): Tr(Y^T Y)/2 = Tr((-Y)Y)/2 = -Tr(Y^2)/2 = -Tr(I)/2 = -1? But we computed earlier that (Y⊗Y)|Φ^+> seemed to give |Φ^+>? Let’s recompute (Y⊗Y) carefully: Y = [[0, -i], [i, 0]]. So Y|0> = i|1>, Y|1> = -i|0>. So (Y⊗Y)|00> = Y|0> ⊗ Y|0> = (i|1>)⊗(i|1>) = i^2 |11> = -|11>. (Y⊗Y)|11> = Y|1>⊗Y|1> = (-i|0>)⊗(-i|0>) = i^2 |00> = -|00>. So (Y⊗Y)|Φ^+> = (Y⊗Y)(|00>+|11>)/√2 = (-|11> - |00>)/√2 = -|Φ^+>. So <Φ^+| (Y⊗Y)|Φ^+> = -1. But then |<Φ^+| (Y⊗Y)|Φ^+>|^2 = 1. So for the squared overlap, the sign doesn’t matter. But careful: Our quantity is |<ψ_ideal| P |ψ_ideal>|^2. And that is the square of the expectation value. So we only care about the magnitude squared. So for a two-qubit Pauli Q on the Bell state, we have: If Q = I⊗I, then expectation = 1, so squared = 1. If Q is such that it is a product of two Paulis that are the same? Actually, from the formula: <Φ^+| (σ ⊗ τ) |Φ^+> = Tr(σ^T τ)/2. For Paulis, if σ and τ are both I, then 1. If σ = X, τ = X, then Tr(X^T X)=Tr(X^2)=Tr(I)=2, so expectation = 1. If σ = Y, τ = Y, then Y^T = -Y, so Tr((-Y)Y)= -Tr(Y^2)= -2, so expectation = -1, squared = 1. If σ = Z, τ = Z, then Tr(Z^T Z)=Tr(Z^2)=2, so expectation = 1. What if σ and τ are different? For example, X⊗Z: Tr(X^T Z)=Tr(X Z)=Tr( -iY ) = 0. So indeed, the expectation is nonzero only if the Pauli on the two qubits is either I⊗I, or X⊗X, Y⊗Y, Z⊗Z. But wait: What about Paulis like X⊗Y? Then Tr(X^T Y)=Tr(X Y)= Tr(iZ)=0. So yes. Also, what about Paulis that are not tensor products? Actually, any two-qubit Pauli is a tensor product up to a phase? Actually, all Paulis are of the form ± σ_a ⊗ σ_b, and since we care about the operator itself (and not the phase? But careful: Our errors are Pauli operators. They are exactly these: They are from the set {I,X,Y,Z} on each qubit, but with a possible overall phase? But Pauli group elements can have factors of i. But typically, we take the Pauli matrices to be Hermitian, so they are X, Y, Z. And then products are Hermitian if an even number of Y’s? Actually, X, Y, Z are Hermitian. The product of two Hermitian matrices is Hermitian if they commute. But X and Y anticommute, so X⊗Y is not Hermitian? Actually, check: (X⊗Y)† = X†⊗Y† = X⊗Y, so it is Hermitian. Because (X⊗Y)† = X†⊗Y† = X⊗Y since X and Y are Hermitian. So indeed, all products of Paulis are Hermitian. But they might have eigenvalues ±1. So it’s fine. So, for a two-qubit Pauli Q on the pair, we have: <Φ^+| Q |Φ^+> is nonzero only if Q is one of: I⊗I, X⊗X, Y⊗Y, Z⊗Z. And in those cases, the expectation is ±1. So then |<Φ^+| Q |Φ^+>|^2 = 1. For any other two-qubit Pauli, the expectation is 0. But wait: There is also the possibility: What about Q = -I⊗I? That is not considered because Pauli errors are usually taken as the matrices themselves, not with a minus sign? But actually, the Pauli group includes ±1, ±i. But in the error model, we said “one of the 15 non-identity two-qubit Paulis”. Typically, the depolarizing channel is defined as: with probability p, apply an error from the set {P} where P are all Pauli operators on two qubits except identity, and they are equally likely. And these Pauli operators are usually taken as the matrices themselves (like X, Y, Z on each qubit, and their products). But note: Could there be a minus sign? For example, -I⊗I is identity up to a phase? Actually, -I is not considered because it’s just a global phase. So we consider the Pauli operators as the set {I, X, Y, Z} on each qubit, and then the two-qubit Paulis are of the form σ_a ⊗ σ_b, where σ_a, σ_b ∈ {I, X, Y, Z}. And there are 16 such operators. And the error channel: with probability 1-p, no error (which is I⊗I). With probability p/15 each, we take any of the other 15. So that is our distribution.

So then, for a given error pattern: E2 is a two-qubit Pauli on qubits (0,3). And E1 is a two-qubit Pauli on qubits (1,2). Then the overall effect is P = E2 ⊗ E1 (but careful: The order: qubits: 0,1,2,3, so actually it is: On qubit0: from E2, on qubit3: from E2, on qubit1: from E1, on qubit2: from E1). So then: <ψ_ideal| P |ψ_ideal> = (<Φ^+|{03}| E2 |Φ^+>{03}) * (<Φ^+|{12}| E1 |Φ^+>{12}). So then |<ψ_ideal| P |ψ_ideal>|^2 = |<Φ^+| E2 |Φ^+>|^2 * |<Φ^+| E1 |Φ^+>|^2.

Now, for a two-qubit Pauli Q (on a pair), we have: If Q = I⊗I, then |<Φ^+| Q |Φ^+>|^2 = 1. If Q is such that it is of the form: X⊗X, Y⊗Y, Z⊗Z, then also |<Φ^+| Q |Φ^+>|^2 = 1. For any other Q (like I⊗X, X⊗I, etc.), |<Φ^+| Q |Φ^+>|^2 = 0.

So, for E2 (on qubits 0,3), the overlap squared is 1 if E2 is either the identity or if E2 is a Pauli that is the same on both qubits (i.e., X⊗X, Y⊗Y, Z⊗Z). How many such Paulis are there? Identity: 1, and then X⊗X, Y⊗Y, Z⊗Z: that’s 3. So total 4 Paulis out of 16 give nonzero overlap (and actually give overlap squared = 1). For any other Pauli (there are 12 of them), the overlap is 0. Similarly for E1 on qubits (1,2): exactly 4 Paulis give overlap squared = 1, and 12 give 0.

And note that these errors are independent. So then the fidelity F is: F = Σ_{E2} Σ_{E1} P(E2) P(E1) * (1 if E2 is in S and E1 is in S, else 0), where S is the set of two-qubit Paulis that are either identity or of the form σ⊗σ (with σ ∈ {X,Y,Z}). So then F = (Probability that E2 is in S) * (Probability that E1 is in S).

Now, what is the probability that a given two-qubit Pauli error (from the depolarizing channel on that gate) is in S? For a given CNOT gate, the error channel: With probability 1-p, we have no error (which is identity, and identity is in S). With probability p, an error occurs, and then it is uniformly distributed among the 15 non-identity Paulis. Among these 15, how many are in S? The ones in S besides identity are: X⊗X, Y⊗Y, Z⊗Z. So that’s 3. So the probability that the error is in S (i.e., the Pauli is such that it doesn’t ruin the overlap) is: (1-p) + (p/15)3 = 1-p + p/5 = 1 - (4/5)p. Wait: Check: If p is the error probability, then probability of identity = 1-p. Probability of each non-identity = p/15. So total probability for being in S = (1-p) + 3(p/15) = 1-p + p/5 = 1 - (4/5)p. So for each CNOT gate, the probability that the error is “good” (i.e., it gives overlap 1) is 1 - (4/5)p. And these are independent for the two CNOT gates. So then the fidelity F = [1 - (4/5)p] * [1 - (4/5)p] = (1 - 4p/5)^2. But wait: Is that the physical state fidelity? That would be the probability that the final state is exactly |ψ_ideal>? But careful: The fidelity is the overlap squared between the actual state and the ideal state. And if the error is such that E2 and E1 are both in S, then (E2 E1)|ψ_ideal> will have overlap 1 with |ψ_ideal>? But is it exactly |ψ_ideal>? Not necessarily: For example, if E2 = X⊗X and E1 = I, then we computed earlier that <ψ_ideal| P|ψ_ideal> was 1/2, but then squared is 1/4? Wait, we computed for P = X_0 X_3 (and E1=I) earlier and got overlap = 1/2, so squared = 1/4. But then according to our factorization, for E2 = X⊗X, we would have <Φ^+| X⊗X |Φ^+> = 1, so that gives 1. There is a discrepancy. Let’s recheck that example: E2 = X⊗X on qubits (0,3) and E1 = I on qubits (1,2). Then our computed final state from that error was: |ψ_final> = (X_0 X_3) |ψ_ideal>. And we computed earlier: |ψ_ideal> = 1/2 (|0000>+|1001>+|0110>+|1111>). Apply X_0 X_3: X_0 X_3|0000> = |1001> X_0 X_3|1001> = |0000> X_0 X_3|0110> = |1110>? Because: |0110>: qubit0=0 -> becomes 1, qubit3=0 -> becomes 1, so becomes |1110>. X_0 X_3|1111> = |0111>. So |ψ_final> = 1/2 (|1001>+|0000>+|1110>+|0111>) = 1/2 (|0000>+|1001>+|0111>+|1110>). Then <ψ_ideal|ψ_final> = 1/4 (1+1+0+0) = 1/2. So the overlap squared is 1/4. But according to our factorization, we said that |ψ_ideal> = |Φ^+>{03} ⊗ |Φ^+>{12}. And then if E2 = X⊗X, then <Φ^+|{03}| X⊗X |Φ^+>{03} = 1 (as computed: (X⊗X)|Φ^+> = |Φ^+>, so expectation = 1). And E1 = I gives expectation 1. So then the product is 1. So there is a discrepancy. Why? Because careful: When we factor, we have to be careful with the ordering of qubits. Our factorization: |ψ_ideal> = (|00>{03}+|11>{03})/√2 ⊗ (|00>{12}+|11>{12})/√2. But is that true? Let’s check: (|00>{03}+|11>{03})/√2 ⊗ (|00>{12}+|11>{12})/√2 = 1/2 (|00,00> + |00,11> + |11,00> + |11,11>) but careful: |00,00> means qubit0=0, qubit3=0, qubit1=0, qubit2=0 -> that is |0000> = A. Good. |00,11> means qubit0=0, qubit3=0, qubit1=1, qubit2=1 -> that is |0110> = C. Good. |11,00> means qubit0=1, qubit3=1, qubit1=0, qubit2=0 -> that is |1001> = B. Good. |11,11> means qubit0=1, qubit3=1, qubit1=1, qubit2=1 -> that is |1111> = D. Good. So indeed, |ψ_ideal> = |Φ^+>{03} ⊗ |Φ^+>{12} if we define |Φ^+> on qubits (0,3) as (|00>+|11>)/√2 and on qubits (1,2) as (|00>+|11>)/√2. So that is correct. Now, if we apply E2 = X⊗X on qubits (0,3), then what is (X⊗X) on the pair (0,3)? But careful: Our E2 is defined as acting on qubits 0 and 3. But in the tensor product ordering, qubit0 is the first qubit of that pair and qubit3 is the second. So if we write |Φ^+>{03} = (|0_0,0_3> + |1_0,1_3>)/√2, then applying X⊗X means applying X on qubit0 and X on qubit3. Then (X⊗X)|Φ^+>{03} = (|1_0,1_3> + |0_0,0_3>)/√2 = |Φ^+>{03}. So indeed, the expectation is 1. So then why did our direct computation give a different result? Because in our direct computation, we applied X_0 X_3 to |ψ_ideal>. But note: X_0 X_3 is exactly the operator that acts as X on qubit0 and X on qubit3. So it should give the same. Let’s recompute carefully with the factored form: |ψ_ideal> = |Φ^+>{03} ⊗ |Φ^+>{12}, where |Φ^+>{03} = (|0_0,0_3> + |1_0,1_3>)/√2, |Φ^+>{12} = (|0_1,0_2> + |1_1,1_2>)/√2. Now, apply X_0 X_3 to |ψ_ideal>. That acts only on the first factor. So we get: (X_0 X_3|Φ^+>{03}) ⊗ |Φ^+>{12} = |Φ^+>{03} ⊗ |Φ^+>{12} = |ψ_ideal>. So it should leave the state invariant. But our earlier computation gave: X_0 X_3|0000> = |1001> (that is correct: |0000> -> X on qubit0: 0->1, X on qubit3: 0->1, so becomes |1001>) X_0 X_3|1001> = |0000> (since qubit0:1->0, qubit3:1->0) X_0 X_3|0110> = ? For |0110>: qubit0=0 -> becomes 1, qubit3=0 -> becomes 1, so becomes |1110>? But wait: |0110> means: qubit0=0, qubit1=1, qubit2=1, qubit3=0. Applying X_0: qubit0 becomes 1. Applying X_3: qubit3 becomes 1. So the result is: |1,1,1,1>? That is |1111> actually. Because qubit1 and qubit2 remain 1. So |0110> should become |1111>? But then what about |1111>? X_0 X_3|1111> = qubit0:1->0, qubit3:1->0, so becomes |0111>? So then if I compute: |ψ_ideal> = 1/2 (|0000>+|1001>+|0110>+|1111>). Apply X_0 X_3: |0000> -> |1001> |1001> -> |0000> |0110> -> |1111>? (Because: 0->1, so qubit0 becomes 1; qubit3:0->1, so qubit3 becomes 1; so indeed, |0110> -> |1,1,1,1> = |1111>) |1111> -> |0111>? (qubit0:1->0, qubit3:1->0, so becomes |0,1,1,1> = |0111>) So then the resulting state is: 1/2 (|1001> + |0000> + |1111> + |0111>) = 1/2 (|0000>+|1001>+|1111>+|0111>). That is different from |ψ_ideal> because |ψ_ideal> has |0110> and |1111>? Actually, |ψ_ideal> has |0110> and |1111>. Here we have |1111> and |0111>. So there is a discrepancy: According to the factored form, X_0 X_3 should leave |Φ^+>{03} invariant. But wait: What is |Φ^+>{03}? It is (|00>{03} + |11>_{03})/√2. That means states where qubit0 and qubit3 are the same. In |0110>, qubit0=0 and qubit3=0, so that is actually |00> in the (0,3) pair? But careful: When we factor, we consider the pair (0,3) as a system. The basis states for that system are: |00>, |01>, |10>, |11>. In |ψ_ideal>, the (0,3) part is either |00> or |11>. So |0000> gives (0,3)=00, (1,2)=00. |1001> gives (0,3)=11, (1,2)=00. |0110> gives (0,3)=00? Because qubit0=0, qubit3=0, so that is |00> actually. But wait: |0110>: qubit0=0, qubit3=0, so that is |00> in the (0,3) pair. And (1,2) part: qubit1=1, qubit2=1, so that is |11>. |1111> gives (0,3)=11, (1,2)=11. So indeed, |ψ_ideal> = 1/2 (|00,00> + |11,00> + |00,11> + |11,11>) in the factored notation where the first tensor factor is (0,3) and the second is (1,2). Now, what is the action of X_0 X_3 on the (0,3) factor? On the (0,3) factor, the operator X_0 X_3 means: apply X to the first qubit of that factor and X to the second qubit of that factor. So if the state of the (0,3) factor is |00>, then X_0 X_3|00> = |11>. If it is |11>, then X_0 X_3|11> = |00>. So actually, X_0 X_3 swaps |00> and |11>. So then applied to |ψ_ideal>, we get: X_0 X_3|ψ_ideal> = 1/2 (|11,00> + |00,00> + |11,11> + |00,11>) = 1/2 (|00,00> + |11,00> + |00,11> + |11,11>) = |ψ_ideal>. So it is invariant. Now, let’s check with the detailed computation using the 4-qubit basis: We have |ψ_ideal> = 1/2 (|0000> + |1001> + |0110> + |1111>). Now, apply X_0 X_3. For |0000>: qubit0=0 -> becomes 1, qubit3=0 -> becomes 1, so becomes |1001>? But careful: |0000> becomes: qubit0: 0->1, so becomes 1; qubit3: 0->1, so becomes 1; so the result is |1,0,0,1> which is |1001>. For |1001>: qubit0=1 -> becomes 0, qubit3=1 -> becomes 0, so becomes |0000>. For |0110>: qubit0=0 -> becomes 1, qubit3=0 -> becomes 1, so becomes |1,1,1,1>? That would be |1111> if qubit1 and qubit2 remain 1. But wait: |0110> has qubit1=1, qubit2=1, so yes, it becomes |1111>. For |1111>: qubit0=1 -> becomes 0, qubit3=1 -> becomes 0, so becomes |0,1,1,1> which is |0111>? That is not |0110>. There’s the error: For |1111>, applying X_0: qubit0: 1->0, so becomes 0; applying X_3: qubit3: 1->0, so becomes 0; so the result is |0,1,1,0> which is |0110> actually. Because |1111> is: qubit0=1, qubit1=1, qubit2=1, qubit3=1. After X_0: becomes 0; after X_3: becomes 0; so we get |0,1,1,0> = |0110>. I earlier mistakenly wrote |0111> but it should be |0110> because qubit2 is 1, so it’s 0,1,1,0. So then: |0000> -> |1001> |1001> -> |0000> |0110> -> |1111> |1111> -> |0110> So indeed, the state becomes: 1/2 (|1001> + |0000> + |1111> + |0110>) which is the same as |ψ_ideal>. So the overlap is 1. My earlier computation for X_0 X_3 had a mistake: I computed |0110> -> |1110> incorrectly. Because |0110>: qubit0=0 -> becomes 1, so qubit0 becomes 1; qubit3=0 -> becomes 1, so qubit3 becomes 1; so the result is |1,1,1,1> = |1111>, not |1110>. And |1111> -> |0,1,1,0> = |0110>. So it’s consistent. So then indeed, for E2 = X⊗X, we get overlap = 1. And similarly, for any E2 that is of the form σ⊗σ (with σ being X, Y, or Z), we get that it actually leaves |Φ^+> invariant up to a sign? But careful: For Y⊗Y, we computed that (Y⊗Y)|Φ^+> = -|Φ^+>, so the expectation is -1, but squared gives 1. So that’s fine. So then our factorization is correct: |<ψ_ideal| P |ψ_ideal>|^2 = 1 if and only if E2 is in S (the set of Paulis on two qubits that are either identity or of the form σ⊗σ) and E1 is in S. And otherwise, it is 0. So then the fidelity F = (Probability that E2 is in S) * (Probability that E1 is in S) = [1 - (4/5)p]^2.

But wait: Is that the physical state fidelity? That would be the probability that the final state is exactly the ideal state. But careful: The actual state is a mixture: ρ = Σ_{E2,E1} p(E2)p(E1) |ψ_{E2,E1}><ψ_{E2,E1}|, where |ψ_{E2,E1}> = (E2 E1) |ψ_ideal>. And then the fidelity is <ψ_ideal| ρ |ψ_ideal> = Σ_{E2,E1} p(E2)p(E1) |<ψ_ideal| E2 E1 |ψ_ideal>|^2. And we’ve computed that this is 1 if both errors are in S, and 0 otherwise. So indeed, F = (1 - 4p/5)^2.

But wait: There is also the possibility that the errors could cancel each other? But they are independent, so no.

So then the answer for Subproblem 1 is: F = (1 - 4p/5)^2.

But is that the physical state fidelity? The problem asks: “What is the physical state fidelity of the final physical 4-qubit state at the end of the circuit as a function of the two-qubit gate error rate p?” That would be F = (1 - 4p/5)^2.

However, let’s double-check with a couple of examples: - If p=0, then F=1. - If p is small, then F ≈ 1 - (8/5)p? Actually, (1 - 4p/5)^2 = 1 - 8p/5 + … So the infidelity is 8p/5 to first order. That seems plausible.

I will check one more: Consider an error on the first CNOT: E1 = X on qubit1 only (so that is I⊗X on qubits (1,2)? Actually, careful: E1 acts on qubits (2,1) but order: qubit2 and qubit1. But we can consider it as a Pauli on qubits (1,2) as well, but careful with order? Actually, the pair is (2,1) but since the Pauli group is the same regardless of order, it’s fine. So take E1 = X on qubit1 only, which is I⊗X on qubits (1,2) if we order as (1,2). That is not in S because for it to be in S, it would need to be something like X⊗X or I⊗I, etc. So then the overlap should be 0. We computed earlier with E1 = X on qubit1 (and E2=I) and got that the final state was orthogonal to |ψ_ideal>. So that gives 0. So that fits.

Thus, the fidelity is indeed the product of the probabilities that each CNOT’s error is in the set that doesn’t disturb the Bell state on that pair. And that probability is: For a given CNOT, the error channel: With probability 1-p, no error (identity). With probability p, an error occurs, and among the 15, only 3 yield a Pauli that is of the form σ⊗σ. So probability that the error is “good” = 1-p + (3/15)p = 1 - p + p/5 = 1 - 4p/5. So then F = (1 - 4p/5)^2.

Therefore, the answer for Subproblem 1 is: Physical state fidelity F = (1 - 4p/5)^2.

I’ll put that in a box.

But wait: The problem might expect an expression in terms of p. So final answer: \boxed{(1 - \frac{4}{5}p)^2}.

I’ll now produce the step-by-step derivation.

We are given: “What is the logical state fidelity of the final 2-qubit logical state at the end of the circuit as a function of the two-qubit gate error rate p, assuming the state is post-selected on all detectable errors in the code?”

Interpretation: We have the same circuit as before. But now we post-select on detectable errors. That means we measure the stabilizers (XXXX and ZZZZ) and only keep outcomes that are +1? But careful: The [[4,2,2]] code is a error detection code. It can detect any single-qubit error? Actually, [[4,2,2]] code has distance 2, so it can detect any single-qubit error. But here errors are two-qubit Pauli errors from the depolarizing channels on the CNOT gates. But these errors act on two qubits. They might be detectable or not.

We want: “logical state fidelity” after post-selection on detectable errors. That means we only consider cases where the error is detected as no error? Or maybe we discard states that are flagged as errors? Actually: “post-selected on all detectable errors” means that if an error occurs that is detectable by the code, we discard that outcome. So we only keep runs where the error is either no error or an error that is not detectable (i.e., a logical error? But wait: The code is [[4,2,2]]. Its stabilizers are XXXX and ZZZZ. The detectable errors are those that anti-commute with at least one stabilizer. The undetectable errors are those that commute with both stabilizers. But note: The undetectable errors include the logical operators. So if an error is undetectable, it might change the logical state. So post-selection: We measure the stabilizers. If we get a non-trivial syndrome (i.e., not +1,+1), we discard the state. So we only keep states that are in the code space. But careful: Even if the state is in the code space, it might not be the ideal state; it could be a logical error. So the logical fidelity after post-selection is the probability that the state is the correct logical state given that it is in the code space.

So we need to compute: F_logical = Probability( state is |ψ_ideal> ) given that the state is in the code space (i.e., no detectable error syndrome). But careful: The errors come from the two CNOT gates. And we already derived that the final state with errors is: |ψ_final> = (E2 E1) |ψ_ideal>, where E1 is on qubits (1,2) and E2 on qubits (0,3). And these are independent Pauli errors with distribution: P(I) = 1-p, and for non-identity, p/15 each.

Now, what errors are detectable? The code stabilizers are S1 = XXXX, S2 = ZZZZ. An error E (a Pauli on 4 qubits) is detectable if it anti-commutes with at least one stabilizer. It is undetectable if it commutes with both stabilizers. The undetectable errors are exactly the logical operators and the stabilizers themselves. But note: Our errors are of the form E = E2 E1, where E2 acts on qubits (0,3) and E1 on qubits (1,2). But careful: Not every Pauli on 4 qubits can be written as such a product? Actually, any Pauli on 4 qubits can be written as a product of a Pauli on qubits (0,3) and a Pauli on qubits (1,2) if we allow? But wait: The Pauli group on 4 qubits is generated by X0, X1, X2, X3, etc. But here our errors come from two independent sources: E2 on {0,3} and E1 on {1,2}. So the overall error is indeed of the form: E = (σ_a ⊗ σ_b) on qubits (0,3) times (σ_c ⊗ σ_d) on qubits (1,2). So it is a product of two two-qubit Paulis on these disjoint sets. So we can analyze separately.

Now, when do we get a detectable error? The syndrome is determined by measuring S1 and S2. For a Pauli error E, the syndrome is (+1,+1) if E commutes with both S1 and S2. Otherwise, it is detectable.

We want to post-select on no detectable error. That means we only keep cases where the error E commutes with both S1 and S2. In those cases, the state is in the code space. But it might be a logical error. So then the logical fidelity is: Among those errors that are undetectable (i.e., commute with stabilizers), what is the probability that the state is exactly |ψ_ideal> (i.e., that the error is actually a stabilizer element? Because if E is a stabilizer, then E|ψ_ideal> = |ψ_ideal> since |ψ_ideal> is a +1 eigenstate of the stabilizers? But careful: The stabilizers are XXXX and ZZZZ. And we already computed that |ψ_ideal> is an eigenstate of both with +1. So if E is a stabilizer (i.e., E = XXXX or ZZZZ or their products? Actually, the stabilizer group is generated by XXXX and ZZZZ. So the stabilizers are: I, XXXX, ZZZZ, XXXX ZZZZ? But note: XXXX ZZZZ = (XXXX)(ZZZZ) but careful: XXXX and ZZZZ commute? They commute because they act on different Pauli bases? Actually, X and Z anticommute on the same qubit, but here it’s product over all qubits: XXXX and ZZZZ: On each qubit, X and Z anticommute, so overall, the commutator is product of (-1) for each qubit? Actually, check: XXXX and ZZZZ: For each qubit, X and Z anticommute, so the overall commutator is (-1)^4 = 1. So they commute. So the stabilizer group has 4 elements: {I, XXXX, ZZZZ, XXXX ZZZZ}. But wait: Are all these stabilizers? They all commute with the logical operators? Yes. So if E is any of these, then E|ψ_ideal> = |ψ_ideal> because |ψ_ideal> is in the code space and is a +1 eigenstate of these stabilizers? But careful: Is |ψ_ideal> a +1 eigenstate of XXXX? Yes, we computed that. For ZZZZ, yes. For XXXX ZZZZ, that is the product. Let’s check: XXXX ZZZZ|ψ_ideal> =? Since XXXX and ZZZZ both give +1, so product gives +1. So indeed, if E is in the stabilizer group, then the state remains |ψ_ideal>. What if E is a logical operator? Then it will map |ψ_ideal> to some other logical state. For example, if E = X_A (which is XIXI), then that will flip the first logical qubit. So then the state becomes (|10>{AB}+|01>{AB})/√2, which is orthogonal to |ψ_ideal>. So then the fidelity would be 0. But wait: There might be errors that are undetectable but are not exactly stabilizers or logical operators? Actually, the set of undetectable errors is exactly the normalizer of the stabilizer group within the Pauli group. That set includes the stabilizers themselves and the logical operators. So indeed, if an error is undetectable, it is either a stabilizer (which leaves the state invariant) or a logical operator (which changes the logical state). But careful: There could be products of stabilizers and logical operators. But that is just a logical operator up to a stabilizer. So essentially, if the error is in the stabilizer group, fidelity remains 1. If it is a logical operator (not including stabilizers), then fidelity becomes 0. But wait: Could there be errors that are undetectable but give a fidelity that is not 0 or 1? For example, if the error is a logical operator that is not a Pauli? But here errors are Pauli. So yes.

So then, after post-selection, we only keep runs where the error E = E2 E1 commutes with both stabilizers. Among those, the probability that the state is correct is the probability that E is actually a stabilizer (i.e., in the set S_stab = {I, XXXX, ZZZZ, XXXX ZZZZ}). But careful: Is it exactly that? Because if E is a stabilizer, then E|ψ_ideal> = |ψ_ideal>. But if E is a logical operator (like X_A), then E|ψ_ideal> is orthogonal to |ψ_ideal>. So indeed, the fidelity conditional on undetectable error is: F_logical = P(E is in stabilizer group | E is undetectable).

But wait: There is also possibility that the error might be such that it is undetectable but still gives some overlap? But as argued, if E is a logical operator (not a stabilizer), then it will map the logical state to a different logical state. And since the logical GHZ state is not an eigenstate of logical operators (except maybe some? Actually, check: X_A on |ψ_ideal> gives? X_A = XIXI. Apply to |ψ_ideal> = 1/2 (|0000>+|1001>+|0110>+|1111>). X_A flips qubit0 and qubit2. X_A|0000> = |1010> (not in support) X_A|1001> = |0011> (not in support) X_A|0110> = |1110> (not in support) X_A|1111> = |0101> (not in support) So indeed, it becomes orthogonal. So fidelity 0.)

So then we need to compute: Total probability that error is undetectable = Sum_{E such that [E, S1]=0 and [E, S2]=0} P(E) And among these, the probability that error is a stabilizer = Sum_{E in stabilizer group} P(E) Then F_logical = (Sum_{E in stabilizer} P(E)) / (Sum_{E undetectable} P(E))

But careful: The errors E are products E = E2 E1, and E2 and E1 are independent with the distribution from the depolarizing channels on the two CNOT gates. So we need to determine for which pairs (E2, E1) the overall error E is undetectable, and for which it is a stabilizer.

Let’s denote: E2 is a Pauli on qubits {0,3}. There are 16 possibilities with probabilities: P(E2=I) = 1-p, and for each non-identity, p/15. Similarly, E1 on qubits {1,2}.

And overall E = E2 ⊗ E1 (but careful: The order: E acts on qubits 0,1,2,3. But since E2 acts on qubits 0 and 3, and E1 acts on qubits 1 and 2, we can write E = (E2)_0 ⊗ (E1)_1 ⊗ (E1)_2 ⊗ (E2)_3? Actually, careful: E2 is a two-qubit Pauli on qubits 0 and 3. So if we write E2 = σ_a ⊗ σ_b (where the first factor acts on qubit0 and the second on qubit3), and E1 = σ_c ⊗ σ_d (where the first factor acts on qubit1 and the second on qubit2), then overall E = σ_a on qubit0, σ_c on qubit1, σ_d on qubit2, σ_b on qubit3. So E = (σ_a, σ_c, σ_d, σ_b).

We want to know when E commutes with S1 = XXXX and S2 = ZZZZ. Let’s compute commutativity conditions.

For a Pauli operator P on 4 qubits, it commutes with XXXX if and only if the number of qubits where P has a Pauli that anticommutes with X is even? Actually, more precisely: XXXX is product of X on each qubit. A Pauli on a single qubit anticommutes with X if it is Y or Z? Actually, X commutes with X, anticommutes with Y and Z? Check: X and Y: XY = iZ, YX = -iZ, so they anticommute. X and Z: anticommute. So indeed, on each qubit, if the Pauli is X, it commutes with X; if it is Y or Z, it anticommutes with X. But careful: There is also phase factors? But for commutativity, we care about the sign. So, P commutes with XXXX if the number of qubits where P has a Pauli that is Y or Z is even. Because each such qubit gives a minus sign, and overall product gives (-1)^(number). So condition: The number of qubits among {0,1,2,3} where the Pauli is Y or Z must be even.

Similarly, P commutes with ZZZZ if the number of qubits where the Pauli is X or Y (since Z anticommutes with X and Y) is even.

So for E = (σ_a on 0, σ_c on 1, σ_d on 2, σ_b on 3), the conditions are: For S1: The number of qubits among {0,1,2,3} such that the Pauli is Y or Z is even. For S2: The number of qubits among {0,1,2,3} such that the Pauli is X or Y is even.

But maybe we can use the factorization into the two pairs? Alternatively, we can determine the undetectable errors in terms of E2 and E1. Note that the code’s undetectable errors are exactly the logical operators and stabilizers. And we already know that the logical operators are: X_A = XIXI, so that corresponds to: on qubit0: X, qubit1: I, qubit2: X, qubit3: I. X_B = XXII, so: qubit0: X, qubit1: X, qubit2: I, qubit3: I. Z_A = ZZII, so: qubit0: Z, qubit1: Z, qubit2: I, qubit3: I. Z_B = ZIZI, so: qubit0: Z, qubit1: I, qubit2: Z, qubit3: I. And the stabilizers: S1 = XXXX: so: all X. S2 = ZZZZ: all Z. And also products like S1 S2 = XXXX ZZZZ = (XZ on each qubit) but careful: XZ = -iY, so that is essentially Y up to a phase? But since we consider Pauli operators, we consider them up to phase? Actually, the stabilizer group is {I, XXXX, ZZZZ, XXXX ZZZZ}. But note: XXXX ZZZZ on each qubit gives XZ which is not a Pauli matrix? Actually, XZ = -iY, so if we take the product over all qubits, we get (-i)^4 * YYYY = YYYY. So actually, S1 S2 is equivalent to YYYY up to a phase? But careful: (-i)^4 = 1, so indeed, XXXX ZZZZ = YYYY. So the stabilizer group can also be written as {I, XXXX, ZZZZ, YYYY}. So the stabilizers are: all X, all Z, all Y. So the undetectable errors are those that are either: - I - XXXX, ZZZZ, YYYY (the stabilizers) - The logical operators: XIXI, XXII, ZZII, ZIZI, and their products with stabilizers? But note: For example, XIXI multiplied by XXXX gives? XIXI * XXXX = (XX on qubit0 gives I, on qubit1: IX = X, on qubit2: XX = I, on qubit3: IX = X) so that gives IXXI? That is not one of the listed ones? Actually, the logical operators are defined up to stabilizers. So the set of undetectable errors is actually the normalizer of the stabilizer group. That set has size 44=16? Actually, the normalizer of a stabilizer group of size 4 in 4 qubits: The code has 2 logical qubits, so the logical Pauli group has 4^2=16 elements. And the stabilizer group has 4 elements. So total undetectable Paulis are 416=64? Wait, careful: The Pauli group on 4 qubits has 4^4=256 elements. The stabilizer group has 4 elements. The normalizer of the stabilizer group in the Pauli group is the set of Paulis that commute with all stabilizers. That set has size: |N(S)| = |S| * |Logical Pauli| = 4 * 16 = 64. So there are 64 undetectable Paulis. But our errors E come from a specific subset: They are of the form E = E2 ⊗ E1, where E2 acts on qubits {0,3} and E1 acts on qubits {1,2}. How many such Paulis are there? There are 16*16=256. So not all Paulis are of that form? Actually, any Pauli on 4 qubits can be written as a product of a Pauli on qubits {0,3} and a Pauli on qubits {1,2}? That would be if the Pauli on qubit1 and qubit2 are independent? But careful: A general Pauli on 4 qubits is of the form: σ_a on qubit0, σ_b on qubit1, σ_c on qubit2, σ_d on qubit3. That can be written as (σ_a on 0 ⊗ σ_d on 3) * (σ_b on 1 ⊗ σ_c on 2) if these commute? They do commute because they act on disjoint sets. So indeed, every Pauli on 4 qubits can be factored as a product of a Pauli on {0,3} and a Pauli on {1,2}. So our errors cover all Paulis. But wait: Our errors are not arbitrary Paulis? They come from two independent distributions on the two pairs. But that distribution is such that the probability of a given Pauli on {0,3} is: if it is identity, probability = 1-p; if it is any of the 15 non-identity, probability = p/15. So that distribution is not uniform over the 16 Paulis? It is: P(E2) = (1-p) for E2=I, and p/15 for each non-identity. So it is not uniform over all 16. Similarly for E1. So we need to compute the probability that E is undetectable. That is: Sum over E2, E1 such that E = E2 E1 commutes with S1 and S2 of P(E2)P(E1).

And then the probability that E is a stabilizer (i.e., leaves the state invariant) is: Sum over E2, E1 such that E is in the stabilizer group {I, XXXX, ZZZZ, YYYY} of P(E2)P(E1). But careful: Is it exactly that? If E is a stabilizer, then E|ψ_ideal> = |ψ_ideal>. But what if E is a stabilizer times a phase? But Pauli operators are defined up to ±1, ±i. But we take them as the matrices themselves. So we need to check: For example, XXXX: That means on qubits 0,1,2,3 all X. Can that be written as E2 E1? XXXX = (X on 0 and X on 3) * (X on 1 and X on 2). So that corresponds to E2 = X⊗X and E1 = X⊗X. So that is included. ZZZZ: E2 = Z⊗Z, E1 = Z⊗Z. YYYY: E2 = Y⊗Y, E1 = Y⊗Y. And also I itself: E2=I, E1=I. So these are the only stabilizers? But wait, what about the product of stabilizers? For example, XXXX * ZZZZ = YYYY (up to phase? Actually, XXXX * ZZZZ = (XZ on each qubit) = (-iY)^4? Actually, careful: On each qubit, XZ = -iY. So product over 4 qubits gives (-i)^4 YYYY = YYYY. So yes.) So indeed, the stabilizer group elements are exactly: I, (X⊗X on both pairs), (Z⊗Z on both pairs), (Y⊗Y on both pairs). But wait: Is that all? What about, say, E2 = X⊗X and E1 = I? That gives overall: X on qubit0, I on qubit1, I on qubit2, X on qubit3. That is not a stabilizer because stabilizers are all same on all qubits. So that is a logical operator? Actually, check: XIXI is a logical operator (X_A). So indeed, the pattern: For the error to be a stabilizer, we need that E2 and E1 are both of the form σ⊗σ and moreover, the σ must be the same for both pairs? Actually, if E2 = X⊗X and E1 = I⊗I, then overall is X on qubit0, I on qubit1, I on qubit2, X on qubit3. That is not all the same. So for it to be a stabilizer, we require that E2 and E1 are both either I⊗I, or both X⊗X, or both Y⊗Y, or both Z⊗Z. But wait: What about E2 = X⊗X and E1 = Z⊗Z? That gives: on qubit0: X, qubit1: Z, qubit2: Z, qubit3: X. That is not a stabilizer because stabilizers are all X or all Z or all Y. So indeed, the condition for E to be a stabilizer is: E2 and E1 are both in the set S (the set we defined earlier: {I⊗I, X⊗X, Y⊗Y, Z⊗Z}) and moreover, they are the same type? Actually, check: If E2 = I⊗I and E1 = X⊗X, then overall: I on qubit0, X on qubit1, X on qubit2, I on qubit3. That is not a stabilizer because stabilizers are all same on all qubits. So indeed, for E to be a stabilizer, we need E2 = E1? But wait: Consider E2 = X⊗X and E1 = X⊗X gives XXXX, so that’s stabilizer. What about E2 = I⊗I and E1 = I⊗I gives I. What about E2 = Y⊗Y and E1 = Y⊗Y gives YYYY. What about E2 = Z⊗Z and E1 = Z⊗Z gives ZZZZ. So indeed, the stabilizer group elements that can occur from our error model are exactly those where E2 and E1 are both the same Pauli from the set {I, X⊗X, Y⊗Y, Z⊗Z}? But careful: Is it necessary that they are the same? What if E2 = X⊗X and E1 = -X⊗X? That would be a phase, but we don’t consider phases? Actually, Pauli operators are defined up to ±1. But typically, we take the representatives as the matrices themselves. So we consider X⊗X as the operator. So yes, they must be exactly the same type. So, the probability that the error is a stabilizer (and thus fidelity remains 1) is: P(stab) = P(E2 = I and E1 = I) + P(E2 = X⊗X and E1 = X⊗X) + P(E2 = Y⊗Y and E1 = Y⊗Y) + P(E2 = Z⊗Z and E1 = Z⊗Z). Now, P(E2 = I) = 1-p. P(E2 = X⊗X) = p/15 (since X⊗X is one of the 15 non-identity Paulis? But careful: There are 15 non-identity Paulis on two qubits. Among them, the ones that are of the form σ⊗σ are: X⊗X, Y⊗Y, Z⊗Z. So each has probability p/15.) Similarly, P(E2 = Y⊗Y) = p/15, P(E2 = Z⊗Z) = p/15. And these are independent for E1. So, P(stab) = (1-p)^2 + 3*(p/15)^2 = (1-p)^2 + 3p^2/225 = (1-p)^2 + p^2/75.

Now, we need the total probability that the error is undetectable. That is: Sum over all E2, E1 such that E = E2 E1 commutes with both stabilizers. We need to determine conditions on E2 and E1 for E to be undetectable. E is undetectable if it commutes with S1 and S2. As derived, E = (σ_a on 0, σ_c on 1, σ_d on 2, σ_b on 3) where E2 gives σ_a and σ_b, and E1 gives σ_c and σ_d. The condition for commuting with S1 (XXXX) is: The number of qubits where the Pauli is Y or Z is even. That is: For qubit0: if σ_a is Y or Z, then count 1; if X or I, then 0. Qubit1: if σ_c is Y or Z, count 1. Qubit2: if σ_d is Y or Z, count 1. Qubit3: if σ_b is Y or Z, count 1. So condition: (indicator(σ_a ∈ {Y,Z}) + (indicator(σ_c ∈ {Y,Z}) + (indicator(σ_d ∈ {Y,Z}) + (indicator(σ_b ∈ {Y,Z}) is even. Similarly, condition for S2 (ZZZZ): The number of qubits where the Pauli is X or Y is even. So: (indicator(σ_a ∈ {X,Y}) + (indicator(σ_c ∈ {X,Y}) + (indicator(σ_d ∈ {X,Y}) + (indicator(σ_b ∈ {X,Y}) is even.

Now, since E2 and E1 are independent, we can compute the probability that a given two-qubit Pauli (on a pair) yields a certain pattern of these indicators? But careful: The conditions involve both pairs together. So we need to sum over all combinations of E2 and E1 that satisfy these parity conditions. Maybe we can compute the probability that a random two-qubit Pauli (with our distribution) has a certain property regarding these indicators? But note: The distribution on each pair is not uniform over the 16 Paulis. It is: P(I⊗I) = 1-p, and for each of the 15 non-identity, probability = p/15. Let’s denote for a two-qubit Pauli Q on a pair (say on qubits i and j). We want to know: What is the probability that Q has a certain pattern? But careful: The conditions are on the individual qubits. For a two-qubit Pauli Q, it is a product of Paulis on two qubits. Let’s denote the two qubits as, say, qubit A and qubit B. For Q, we have two Paulis: one on A and one on B. They can be I, X, Y, Z. We want to compute: For a given Q, define: a = indicator( Pauli on A is Y or Z ) b = indicator( Pauli on B is Y or Z ) Then the contribution to the first condition from this pair is a+b mod 2. Similarly, for the second condition: c = indicator( Pauli on A is X or Y ) d = indicator( Pauli on B is X or Y ) Then contribution mod 2 is c+d mod 2.

So if we have two independent pairs: pair from E2 and pair from E1. Let for E2, let these be: A: qubit0, B: qubit3. For E1: qubit1 and qubit2. Then the overall condition for S1 is: (a2 + b2) + (a1 + b1) is even mod 2, i.e., (a2+b2) and (a1+b1) have the same parity. For S2: (c2+d2) and (c1+d1) have the same parity.

So if we can compute for a given two-qubit Pauli Q (with our distribution) the probability that (a+b) mod 2 is 0 or 1, and similarly for (c+d) mod 2. And since E2 and E1 are independent, the probability that both conditions hold is: P(undetectable) = [P((a+b) mod 2 = 0 for E2) * P((a+b) mod 2 = 0 for E1) + P((a+b)=1 for E2)*P((a+b)=1 for E1)] * [similar for the second condition? Actually, careful: The conditions are independent? They must both hold. So we need: P(undetectable) = Sum over i,j,k,l in {0,1} such that i = k and j = l? Actually, we need: Condition1: (from E2, let X = (a+b) mod 2) and from E1, let Y = (a+b) mod 2) must be equal. Condition2: from E2, let U = (c+d) mod 2) and from E1, let V = (c+d) mod 2) must be equal. And these two conditions are independent? Not exactly independent because for a given Pauli, the values of (a+b) and (c+d) are related. So we need to compute the joint distribution for a two-qubit Pauli Q from our distribution: We want: For a given Q (on two qubits), what are the probabilities: P( (a+b) mod 2 = 0, (c+d) mod 2 = 0 ) P( (a+b)=0, (c+d)=1 ) P( (a+b)=1, (c+d)=0 ) P( (a+b)=1, (c+d)=1 ) And then since E2 and E1 are independent, the overall probability that E is undetectable is: P(undetectable) = Sum_{r,s in {0,1}} [P(E2 gives (r,s)) * P(E1 gives (r,s))] because we need both pairs to have the same (a+b) mod 2 and the same (c+d) mod 2. So if we denote: q_{00} = probability that for a given pair (from our distribution), we have (a+b)=0 and (c+d)=0. q_{01} = probability for (0,1) q_{10} = probability for (1,0) q_{11} = probability for (1,1) Then P(undetectable) = q_{00}^2 + q_{01}^2 + q_{10}^2 + q_{11}^2.

And we already computed that if the error is a stabilizer (i.e., leaves state invariant), that corresponds to? For a stabilizer, we need E2 and E1 to be both of the form: either I⊗I, or X⊗X, or Y⊗Y, or Z⊗Z. And what are the (a+b) and (c+d) for these? For I⊗I: On each qubit: I. So a: indicator(Y or Z) = 0, so a+b=0. c: indicator(X or Y)=0, so c+d=0. So (0,0). For X⊗X: On each qubit: X. So a: X is not Y or Z? Actually, X is not Y or Z, so a=0. So a+b=0. c: X is X or Y? Yes, so c=1. So c+d=1+1=2 mod2 = 0? Actually, careful: c+d mod2: For each qubit, if it is X, then indicator(X or Y)=1. So for two qubits, that sum is 2, which is even, so mod2 = 0. So (0,0) actually? Wait: Check: For X⊗X, we have: qubit0: X -> indicator for S1: is it Y or Z? No, so 0. qubit3: X -> 0. So total mod2 = 0. For S2: qubit0: X -> indicator for X or Y = 1; qubit3: X -> 1; sum=2, mod2=0. So indeed, X⊗X gives (0,0). For Y⊗Y: On each qubit: Y. Then a: Y is Y or Z? Yes, so a=1. So a+b=1+1=2 mod2=0. c: Y is X or Y? Yes, so c=1. So c+d=2 mod2=0. So Y⊗Y gives (0,0). For Z⊗Z: On each qubit: Z. Then a: Z is Y or Z? Yes, so a=1. So a+b=2 mod2=0. c: Z is X or Y? No, so c=0. So c+d=0 mod2=0. So all these stabilizer type errors (I, X⊗X, Y⊗Y, Z⊗Z) give (0,0). But wait, are there other errors that give (0,0) but are not stabilizers? For example, consider E2 = I⊗X? That is: on qubit0: I, on qubit3: X. Then a: I gives 0, X gives 0, so a+b=0. c: I gives 0, X gives 1, so c+d=1 mod2=1. So that gives (0,1) actually. What about E2 = X⊗I? That gives: qubit0: X -> 0; qubit3: I -> 0; so a+b=0. c: X -> 1; I -> 0; so c+d=1 mod2=1. So (0,1). What about E2 = I⊗Z? Then a: I=0, Z=1 so sum=1 mod2=1; c: I=0, Z=0 so sum=0; so (1,0). What about E2 = Z⊗I? Then (1,0) as well. What about E2 = X⊗Z? Then: qubit0: X -> 0; qubit3: Z -> 1; so a+b=1 mod2=1; c: X -> 1; Z -> 0; so c+d=1 mod2=1; so (1,1). What about E2 = Z⊗X? That gives: qubit0: Z -> 1; qubit3: X -> 0; so a+b=1; c: Z -> 0; X -> 1; so c+d=1; so (1,1). What about E2 = Y⊗I? Then: qubit0: Y -> 1; qubit3: I -> 0; so a+b=1; c: Y -> 1; I -> 0; so c+d=1; so (1,1)? Wait, check Y⊗I: Y on qubit0: indicator for S1: Y is Y or Z -> 1; I: 0; so sum=1 mod2=1. For S2: Y is X or Y -> 1; I: 0; so sum=1 mod2=1. So (1,1). What about E2 = I⊗Y? Then: (1,1) as well. What about E2 = X⊗Y? Then: qubit0: X -> 0; qubit3: Y -> 1; so a+b=1; c: X -> 1; Y -> 1; so c+d=2 mod2=0; so (1,0)? Check X⊗Y: Actually, careful: X⊗Y: On qubit0: X -> for S1: not Y or Z so 0; on qubit3: Y -> yes, so 1; so total=1 mod2=1. For S2: qubit0: X -> yes (X or Y) so 1; qubit3: Y -> yes so 1; total=2 mod2=0. So (1,0). What about E2 = Y⊗X? Then: qubit0: Y -> 1; qubit3: X -> 0; so a+b=1; c: Y -> 1; X -> 1; so c+d=2 mod2=0; so (1,0). What about E2 = Y⊗Z? Then: qubit0: Y -> 1; qubit3: Z -> 1; so a+b=2 mod2=0; c: Y -> 1; Z -> 0; so c+d=1 mod2=1; so (0,1). What about E2 = Z⊗Y? Then: (0,1) as well. So indeed, for a two-qubit Pauli, the pair ( (a+b) mod2, (c+d) mod2 ) can be (0,0), (0,1), (1,0), (1,1). And we want to compute the probabilities of these from our distribution on a pair. Our distribution on a pair (say E2) is: P(I⊗I) = 1-p. For each of the 15 non-identity Paulis, probability = p/15. So we need to list all 16 Paulis on two qubits and determine their (a+b, c+d) values. Let’s do that systematically. Label the two qubits as A and B. The Pauli is σ_A ⊗ σ_B, where σ_A, σ_B ∈ {I, X, Y, Z}. But careful: The order: For E2, the qubits are 0 and 3. But that doesn’t matter; it’s just two qubits. We want: a = 1 if σ is Y or Z, else 0. c = 1 if σ is X or Y, else 0. So for each Pauli, we compute: For I: a=0, c=0. For X: a=0, c=1. For Y: a=1, c=1. For Z: a=1, c=0. Now, for a two-qubit Pauli Q = σ ⊗ τ, we have: (a+b) mod2: actually, we want the sum of the a values from both qubits modulo 2. But careful: The condition for commuting with S1 is that the total number of qubits with a=1 is even. That is exactly (a_A + a_B) mod 2 = 0. Similarly, (c_A + c_B) mod 2 = 0 for commuting with S2. But wait: That is for the pair itself? Actually, for a single pair, if we consider it alone, the condition for that pair to contribute an even number to the overall count? But we are going to combine two pairs. So we want the values: Let u = (a_A + a_B) mod 2. Let v = (c_A + c_B) mod 2. So for each two-qubit Pauli, we compute (u,v).

List all 16:

1. I⊗I: (0+0 mod2=0, 0+0=0) -> (0,0)
2. I⊗X: (0+0=0, 0+1=1 mod2=1) -> (0,1)
3. I⊗Y: (0+1=1, 0+1=1) -> (1,1)
4. I⊗Z: (0+1=1, 0+0=0) -> (1,0)
5. X⊗I: (0+0=0, 1+0=1) -> (0,1)
6. X⊗X: (0+0=0, 1+1=2 mod2=0) -> (0,0)
7. X⊗Y: (0+1=1, 1+1=2 mod2=0) -> (1,0)? Check X⊗Y: a: X=0, Y=1 so sum=1 mod2=1; c: X=1, Y=1 so sum=2 mod2=0 -> (1,0)
8. X⊗Z: (0+1=1, 1+0=1) -> (1,1)
9. Y⊗I: (1+0=1, 1+0=1) -> (1,1)
10. Y⊗X: (1+0=1, 1+1=2 mod2=0) -> (1,0)
11. Y⊗Y: (1+1=2 mod2=0, 1+1=2 mod2=0) -> (0,0)
12. Y⊗Z: (1+1=0, 1+0=1) -> (0,1)
13. Z⊗I: (1+0=1, 0+0=0) -> (1,0)
14. Z⊗X: (1+0=1, 0+1=1) -> (1,1)
15. Z⊗Y: (1+1=0, 0+1=1) -> (0,1)
16. Z⊗Z: (1+1=0, 0+0=0) -> (0,0)

So, summary: (0,0): I⊗I, X⊗X, Y⊗Y, Z⊗Z. That is 4 Paulis. (0,1): I⊗X, X⊗I, Y⊗Z, Z⊗Y. That is 4 Paulis. (1,0): I⊗Z, X⊗Y, Y⊗X, Z⊗I. That is 4 Paulis. (1,1): I⊗Y, X⊗Z, Y⊗I, Z⊗X. That is 4 Paulis.

So indeed, uniformly over the 16 Paulis, each category has probability 4/16 = 1/4. But our distribution is not uniform: It gives probability 1-p for I⊗I, and for each of the other 15, probability p/15. So then, for a given pair (say E2), the probabilities: q_{00} = P((u,v)=(0,0)) = P(I⊗I) + P(X⊗X) + P(Y⊗Y) + P(Z⊗Z) but careful: I⊗I has probability 1-p, while the others have probability p/15 each. So q_{00} = (1-p) + 3*(p/15) = 1-p + p/5 = 1 - (4/5)p. q_{01} = There are 4 Paulis in (0,1): each probability p/15, so q_{01} = 4p/15. q_{10} = 4p/15. q_{11} = 4p/15.

And note that these sum to: (1-p + 3p/15) + 4p/15+4p/15+4p/15 = 1-p + 15p/15 = 1-p+p = 1. Good.

So for a single pair, the probability distribution over (u,v) is: (0,0): 1 - 4p/5 (0,1): 4p/15 (1,0): 4p/15 (1,1): 4p/15

Now, for the overall error E to be undetectable, we need that from E2 and E1, the (u,v) values are the same. Because then the total number of qubits with a=1 will be even (since even+even=even, odd+odd=even) and similarly for c. So, P(undetectable) = = [q_{00}^2 + q_{01}^2 + q_{10}^2 + q_{11}^2] = (1 - 4p/5)^2 + 3(4p/15)^2 = (1 - 4p/5)^2 + 3(16p^2/225) = (1 - 4p/5)^2 + 48p^2/225 = (1 - 8p/5 + 16p^2/25) + 48p^2/225 Note that 16/25 = (144/225) so: = 1 - 8p/5 + (144p^2/225 + 48p^2/225) = 1 - 8p/5 + 192p^2/225 Simplify 192/225 = divide by 3: 192/3=64, 225/3=75, so = 64p^2/75. So P(undetectable) = 1 - (8/5)p + (64/75)p^2.

But wait, is that the total probability that the error is undetectable? Yes, because we require that both conditions hold: the (u,v) from E2 equals the (u,v) from E1.

Now, what is the probability that the error is a stabilizer? That is when E2 and E1 are both in the set that gives (0,0) and moreover, they must be the same type? Actually, careful: For the error to be a stabilizer, we need that E itself is exactly I, XXXX, ZZZZ, or YYYY. But as we computed, that requires that E2 and E1 are both from the set {I⊗I, X⊗X, Y⊗Y, Z⊗Z} and moreover, they are the same specific type? But wait: If E2 = I⊗I and E1 = X⊗X, then overall E = (I⊗I)(X⊗X) = X⊗X on the pair? But that gives: on qubit0: IX = X, on qubit3: IX = X, on qubit1: X⊗? Actually, careful: E2 acts on qubits 0 and 3, E1 acts on qubits 1 and 2. So if E2 = I⊗I and E1 = X⊗X, then overall: qubit0: I, qubit1: X, qubit2: X, qubit3: I. That is not XXXX. That is actually the logical operator X_B? Actually, X_B = XXII, so that would be X on qubit0 and qubit1, but here we have I on qubit0 and X on qubit1? Actually, wait: X_B is defined as XXII. That means X on qubit0 and qubit1, and I on qubit2 and qubit3. So our E here is: qubit0: I, qubit1: X, qubit2: X, qubit3: I. That is not X_B. It is actually? It commutes with stabilizers? Let’s check: This E = I, X, X, I. For S1: XXXX, the number of Y/Z? qubit0: I ->0; qubit1: X ->0; qubit2: X ->0; qubit3: I ->0; so even, so okay. For S2: ZZZZ, qubit0: I->0; qubit1: X ->1; qubit2: X ->1; qubit3: I->0; sum=2 mod2=0. So it is undetectable. But is it a stabilizer? Stabilizers are all same on all qubits. Here we have I on qubits 0 and 3, and X on qubits 1 and 2. That is not all the same. So it is a logical operator. So indeed, for the error to be a stabilizer, we need that E2 and E1 are not only both in the (0,0) category, but they must be exactly the same Pauli? Let’s check: If E2 = I⊗I and E1 = I⊗I, then E = I, stabilizer. If E2 = X⊗X and E1 = X⊗X, then E = XXXX, stabilizer. If E2 = Y⊗Y and E1 = Y⊗Y, then E = YYYY, stabilizer. If E2 = Z⊗Z and E1 = Z⊗Z, then E = ZZZZ, stabilizer. What if E2 = I⊗I and E1 = X⊗X? Then as above, E = (I on 0,3) and (X on 1,2) gives: qubit0: I, qubit1: X, qubit2: X, qubit3: I. That is not a stabilizer because stabilizers are all X, all Y, all Z, or all I. So that is a logical operator. So indeed, for the error to be a stabilizer, we need that E2 and E1 are both in the set {I⊗I, X⊗X, Y⊗Y, Z⊗Z} and that they are the same element. But wait: What about E2 = X⊗X and E1 = I⊗I? That is the same as above, so no. So then, the probability that the error is a stabilizer is: P(stab) = P(E2 = I⊗I and E1 = I⊗I) + P(E2 = X⊗X and E1 = X⊗X) + P(E2 = Y⊗Y and E1 = Y⊗Y) + P(E2 = Z⊗Z and E1 = Z⊗Z) = (1-p)^2 + 3(p/15)^2 = (1-p)^2 + 3p^2/225 = (1-p)^2 + p^2/75.

So then the logical fidelity after post-selection is: F_logical = P(stab) / P(undetectable) = [ (1-p)^2 + p^2/75 ] / [ 1 - (8/5)p + (64/75)p^2 ].

We can simplify that expression: Numerator: (1-p)^2 + p^2/75 = 1 - 2p + p^2 + p^2/75 = 1 - 2p + (76/75)p^2? Actually, p^2 + p^2/75 = (75/75 + 1/75)p^2 = 76/75 p^2. So numerator = 1 - 2p + (76/75)p^2. Denom: 1 - (8/5)p + (64/75)p^2 = 1 - (40/25)p? Actually, better: Write denominator with denominator 75: = (75/75) - (8/5)p = (8/5) = (120/75) so denominator = (75/75) - (120/75)p + (64/75)p^2 = [75 - 120p + 64p^2] / 75. So F_logical = [1 - 2p + (76/75)p^2] / ( [75 - 120p + 64p^2] / 75 ) = 75 * [1 - 2p + (76/75)p^2] / [75 - 120p + 64p^2] = [75 - 150p + 76p^2] / [75 - 120p + 64p^2].

So F_logical = (75 - 150p + 76p^2) / (75 - 120p + 64p^2).

That is the logical fidelity after post-selection on detectable errors.

But wait: Is that the final answer? It might be simplified further? Divide numerator and denominator by? They don’t have a common factor. So that is the expression.

We can also write it as: F_logical = 1 - (30p - 12p^2) / (75 - 120p + 64p^2)? Not simpler.

So the answer for Subproblem 2 is: F_logical = (75 - 150p + 76p^2) / (75 - 120p + 64p^2).

But let’s check for small p: When p=0, F_logical = 75/75 = 1. First order: Expand: Numerator: 75 - 150p + … Denom: 75 - 120p + … So F_logical ≈ 1 - (150-120)/75 p = 1 - (30/75)p = 1 - (2/5)p = 1 - 0.4p. So the infidelity to first order is 0.4p. That makes sense.

So final answer for Subproblem 2.

I’ll now produce the step-by-step derivation.

We are given the expression for the logical state fidelity after post-selection: F_logical = (75 - 150p + 76p^2) / (75 - 120p + 64p^2)

We will represent this in SymPy. Note: The expression is already simplified.

Let’s write it in the function.