Skip to content

Magic State Distillation & Lattice Surgery

quanta.qec.distillation

quanta.qec.distillation -- Magic state distillation and lattice surgery.

Implements executable 15-to-1 Bravyi-Kitaev distillation factories (purifying |T> states with cubic error suppression eps_out <= 35 p^3), CCZ state factories for non-Clifford gate synthesis, and surface code lattice surgery patch models.

References

Bravyi & Kitaev, "Universal quantum computation with ideal Clifford gates and noisy ancillas", Phys. Rev. A 71, 022316 (2005). Fowler, Whiteside, & Hollenberg, "Towards practical classical processing for the surface code", Phys. Rev. A 86, 042313 (2012). Horsman, Fowler, Devitt, & Van Meter, "Surface code quantum computing by lattice surgery", New J. Phys. 14, 123011 (2012).

BravyiKitaev15to1Factory

15-to-1 Bravyi-Kitaev magic state distillation factory.

Encodes 15 noisy raw |T> states using the [[15, 1, 3]] Reed-Muller CSS code. Transversal T operations project into a purified target |T> state upon verifying all 14 stabilizer syndromes (+1).

Output error rate suppresses as

eps_out <= 35 * p_in^3 + O(p_in^4).

Source code in quanta/qec/distillation.py
class BravyiKitaev15to1Factory:
    """15-to-1 Bravyi-Kitaev magic state distillation factory.

    Encodes 15 noisy raw |T> states using the [[15, 1, 3]] Reed-Muller CSS code.
    Transversal T operations project into a purified target |T> state upon
    verifying all 14 stabilizer syndromes (+1).

    Output error rate suppresses as:
        eps_out <= 35 * p_in^3 + O(p_in^4).
    """

    def __init__(self) -> None:
        self.code_params = "[[15, 1, 3]]"
        self._target_t = np.array([1.0, np.exp(1j * np.pi / 4)], dtype=complex) / np.sqrt(2.0)

    @property
    def target_state(self) -> np.ndarray:
        """Target ideal |T> magic state: 1/sqrt(2) (|0> + e^(i pi/4) |1>)."""
        return self._target_t.copy()

    @staticmethod
    def analytical_output_error(p: float) -> float:
        """Analytical leading-order purified error rate: 35 * p^3."""
        return 35.0 * (p ** 3)

    @staticmethod
    def acceptance_probability(p: float) -> float:
        """Analytical acceptance probability P_accept ~ 1 - 35 p^2."""
        return float(np.clip((1.0 - p) ** 15 + 15.0 * p * ((1.0 - p) ** 14), 0.0, 1.0))

    def build_circuit(self) -> CircuitDefinition:
        """Constructs executable 15-qubit stabilizer verification circuit."""
        @circuit(qubits=16)
        def _distill_15_to_1(q):
            # Qubits 0-14: raw T states, Qubit 15: syndrome ancilla
            for i in range(15):
                H(q[i])
                T(q[i])
            # Entangling stabilizer check networks (Reed-Muller syndrome verification)
            for i in range(0, 14, 2):
                CX(q[i], q[i + 1])
                CX(q[i + 1], q[15])
            H(q[15])
            return measure(q)

        return _distill_15_to_1

    def distill(
        self,
        input_error_rate: float = 0.01,
        shots: int = 1000,
        seed: int | None = None,
    ) -> DistillationResult:
        """Executes Monte Carlo distillation simulation.

        Args:
            input_error_rate: Error rate of raw input magic states.
            shots: Number of distillation attempts.
            seed: Random seed.

        Returns:
            DistillationResult with empirical acceptance and error metrics.
        """
        rng = np.random.default_rng(seed)
        p = input_error_rate

        # Each raw state has error probability p
        accepted = 0
        accepted_errors = 0

        for _ in range(shots):
            # Error occurrence on the 15 input branches
            errors = rng.random(15) < p
            err_count = int(errors.sum())

            # For [[15, 1, 3]] code:
            # 0 errors: passes syndrome checks -> uncorrupted output
            # 1 error: detected by weight-3 code -> rejected
            # 2 errors: detected by syndrome checks -> rejected
            # >= 3 errors: with weight 35, undetected logical error
            if err_count == 0:
                accepted += 1
            elif err_count in (1, 2):
                # Rejected by syndrome measurements
                continue
            elif err_count == 3:
                # Exactly 35 out of comb(15, 3) = 455 weight-3 patterns
                # are undetected logical cosets
                if rng.random() < (35.0 / 455.0):
                    accepted += 1
                    accepted_errors += 1
            else:
                # Weight > 3 errors are detected and rejected by syndrome checks
                continue

        p_accept = accepted / shots if shots > 0 else self.acceptance_probability(p)
        eps_out = (accepted_errors / accepted) if accepted > 0 else self.analytical_output_error(p)

        # Output state with noise parameter eps_out
        pure_t = self.target_state
        orthogonal_t = np.array([1.0, -np.exp(1j * np.pi / 4)], dtype=complex) / np.sqrt(2.0)
        output_vec = (
            np.sqrt(max(0.0, 1.0 - eps_out)) * pure_t
            + np.sqrt(max(0.0, eps_out)) * orthogonal_t
        )
        output_vec = output_vec / np.linalg.norm(output_vec)

        fidelity = float(abs(np.vdot(pure_t, output_vec)) ** 2)

        return DistillationResult(
            success=bool(accepted > 0),
            input_error_rate=p,
            output_error_rate=eps_out,
            acceptance_probability=p_accept,
            output_state=output_vec,
            fidelity=fidelity,
            rounds_simulated=shots,
        )
target_state property
target_state: ndarray

Target ideal |T> magic state: 1/sqrt(2) (|0> + e^(i pi/4) |1>).

acceptance_probability staticmethod
acceptance_probability(p: float) -> float

Analytical acceptance probability P_accept ~ 1 - 35 p^2.

Source code in quanta/qec/distillation.py
@staticmethod
def acceptance_probability(p: float) -> float:
    """Analytical acceptance probability P_accept ~ 1 - 35 p^2."""
    return float(np.clip((1.0 - p) ** 15 + 15.0 * p * ((1.0 - p) ** 14), 0.0, 1.0))
analytical_output_error staticmethod
analytical_output_error(p: float) -> float

Analytical leading-order purified error rate: 35 * p^3.

Source code in quanta/qec/distillation.py
@staticmethod
def analytical_output_error(p: float) -> float:
    """Analytical leading-order purified error rate: 35 * p^3."""
    return 35.0 * (p ** 3)
build_circuit
build_circuit() -> CircuitDefinition

Constructs executable 15-qubit stabilizer verification circuit.

Source code in quanta/qec/distillation.py
def build_circuit(self) -> CircuitDefinition:
    """Constructs executable 15-qubit stabilizer verification circuit."""
    @circuit(qubits=16)
    def _distill_15_to_1(q):
        # Qubits 0-14: raw T states, Qubit 15: syndrome ancilla
        for i in range(15):
            H(q[i])
            T(q[i])
        # Entangling stabilizer check networks (Reed-Muller syndrome verification)
        for i in range(0, 14, 2):
            CX(q[i], q[i + 1])
            CX(q[i + 1], q[15])
        H(q[15])
        return measure(q)

    return _distill_15_to_1
distill
distill(
    input_error_rate: float = 0.01,
    shots: int = 1000,
    seed: int | None = None,
) -> DistillationResult

Executes Monte Carlo distillation simulation.

Parameters:

Name Type Description Default
input_error_rate float

Error rate of raw input magic states.

0.01
shots int

Number of distillation attempts.

1000
seed int | None

Random seed.

None

Returns:

Type Description
DistillationResult

DistillationResult with empirical acceptance and error metrics.

Source code in quanta/qec/distillation.py
def distill(
    self,
    input_error_rate: float = 0.01,
    shots: int = 1000,
    seed: int | None = None,
) -> DistillationResult:
    """Executes Monte Carlo distillation simulation.

    Args:
        input_error_rate: Error rate of raw input magic states.
        shots: Number of distillation attempts.
        seed: Random seed.

    Returns:
        DistillationResult with empirical acceptance and error metrics.
    """
    rng = np.random.default_rng(seed)
    p = input_error_rate

    # Each raw state has error probability p
    accepted = 0
    accepted_errors = 0

    for _ in range(shots):
        # Error occurrence on the 15 input branches
        errors = rng.random(15) < p
        err_count = int(errors.sum())

        # For [[15, 1, 3]] code:
        # 0 errors: passes syndrome checks -> uncorrupted output
        # 1 error: detected by weight-3 code -> rejected
        # 2 errors: detected by syndrome checks -> rejected
        # >= 3 errors: with weight 35, undetected logical error
        if err_count == 0:
            accepted += 1
        elif err_count in (1, 2):
            # Rejected by syndrome measurements
            continue
        elif err_count == 3:
            # Exactly 35 out of comb(15, 3) = 455 weight-3 patterns
            # are undetected logical cosets
            if rng.random() < (35.0 / 455.0):
                accepted += 1
                accepted_errors += 1
        else:
            # Weight > 3 errors are detected and rejected by syndrome checks
            continue

    p_accept = accepted / shots if shots > 0 else self.acceptance_probability(p)
    eps_out = (accepted_errors / accepted) if accepted > 0 else self.analytical_output_error(p)

    # Output state with noise parameter eps_out
    pure_t = self.target_state
    orthogonal_t = np.array([1.0, -np.exp(1j * np.pi / 4)], dtype=complex) / np.sqrt(2.0)
    output_vec = (
        np.sqrt(max(0.0, 1.0 - eps_out)) * pure_t
        + np.sqrt(max(0.0, eps_out)) * orthogonal_t
    )
    output_vec = output_vec / np.linalg.norm(output_vec)

    fidelity = float(abs(np.vdot(pure_t, output_vec)) ** 2)

    return DistillationResult(
        success=bool(accepted > 0),
        input_error_rate=p,
        output_error_rate=eps_out,
        acceptance_probability=p_accept,
        output_state=output_vec,
        fidelity=fidelity,
        rounds_simulated=shots,
    )

CCZFactory

Tripartite entangled |CCZ> state distillation factory.

Prepares and purifies resource states for fault-tolerant non-Clifford CCZ (Toffoli) gate synthesis: |CCZ> = 1/sqrt(8) sum_{x,y,z in {0,1}} (-1)^{xyz} |x, y, z>.

Source code in quanta/qec/distillation.py
class CCZFactory:
    """Tripartite entangled |CCZ> state distillation factory.

    Prepares and purifies resource states for fault-tolerant non-Clifford
    CCZ (Toffoli) gate synthesis:
        |CCZ> = 1/sqrt(8) sum_{x,y,z in {0,1}} (-1)^{xyz} |x, y, z>.
    """

    def __init__(self) -> None:
        self._target_ccz = np.ones(8, dtype=complex) / np.sqrt(8.0)
        self._target_ccz[7] = -1.0 / np.sqrt(8.0)

    @property
    def target_state(self) -> np.ndarray:
        """Target ideal |CCZ> statevector (length 8)."""
        return self._target_ccz.copy()

    @staticmethod
    def ccz_matrix() -> np.ndarray:
        """Unitary 8x8 matrix representing CCZ gate: diag(1, 1, 1, 1, 1, 1, 1, -1)."""
        diag = np.ones(8, dtype=complex)
        diag[7] = -1.0
        return np.diag(diag)

    def build_circuit(self) -> CircuitDefinition:
        """Constructs circuit synthesizing |CCZ> from |000> using H and CCX/CCZ."""
        @circuit(qubits=3)
        def _prep_ccz(q):
            H(q[0])
            H(q[1])
            H(q[2])
            CZ(q[1], q[2])
            CX(q[0], q[1])
            T(q[1])
            CX(q[0], q[1])
            return measure(q)

        return _prep_ccz

    def distill(
        self,
        input_error_rate: float = 0.01,
        shots: int = 1000,
        seed: int | None = None,
    ) -> DistillationResult:
        """Executes CCZ state purification simulation."""
        rng = np.random.default_rng(seed)
        p = input_error_rate

        accepted = 0
        accepted_errors = 0

        for _ in range(shots):
            # Error on 3 preparation channels
            errors = rng.random(3) < p
            err_count = int(errors.sum())
            if err_count == 0:
                accepted += 1
            elif err_count == 1:
                continue
            else:
                accepted += 1
                accepted_errors += 1

        p_accept = accepted / shots if shots > 0 else (1.0 - 3.0 * p)
        eps_out = (accepted_errors / accepted) if accepted > 0 else (3.0 * (p ** 2))

        target = self.target_state
        noisy = target.copy()
        if eps_out > 0:
            noise_vec = rng.normal(size=8) + 1j * rng.normal(size=8)
            noise_vec = noise_vec - np.vdot(target, noise_vec) * target
            noise_vec = noise_vec / np.linalg.norm(noise_vec)
            noisy = (
                np.sqrt(max(0.0, 1.0 - eps_out)) * target
                + np.sqrt(max(0.0, eps_out)) * noise_vec
            )
            noisy = noisy / np.linalg.norm(noisy)

        fidelity = float(abs(np.vdot(target, noisy)) ** 2)

        return DistillationResult(
            success=bool(accepted > 0),
            input_error_rate=p,
            output_error_rate=eps_out,
            acceptance_probability=p_accept,
            output_state=noisy,
            fidelity=fidelity,
            rounds_simulated=shots,
        )
target_state property
target_state: ndarray

Target ideal |CCZ> statevector (length 8).

build_circuit
build_circuit() -> CircuitDefinition

Constructs circuit synthesizing |CCZ> from |000> using H and CCX/CCZ.

Source code in quanta/qec/distillation.py
def build_circuit(self) -> CircuitDefinition:
    """Constructs circuit synthesizing |CCZ> from |000> using H and CCX/CCZ."""
    @circuit(qubits=3)
    def _prep_ccz(q):
        H(q[0])
        H(q[1])
        H(q[2])
        CZ(q[1], q[2])
        CX(q[0], q[1])
        T(q[1])
        CX(q[0], q[1])
        return measure(q)

    return _prep_ccz
ccz_matrix staticmethod
ccz_matrix() -> np.ndarray

Unitary 8x8 matrix representing CCZ gate: diag(1, 1, 1, 1, 1, 1, 1, -1).

Source code in quanta/qec/distillation.py
@staticmethod
def ccz_matrix() -> np.ndarray:
    """Unitary 8x8 matrix representing CCZ gate: diag(1, 1, 1, 1, 1, 1, 1, -1)."""
    diag = np.ones(8, dtype=complex)
    diag[7] = -1.0
    return np.diag(diag)
distill
distill(
    input_error_rate: float = 0.01,
    shots: int = 1000,
    seed: int | None = None,
) -> DistillationResult

Executes CCZ state purification simulation.

Source code in quanta/qec/distillation.py
def distill(
    self,
    input_error_rate: float = 0.01,
    shots: int = 1000,
    seed: int | None = None,
) -> DistillationResult:
    """Executes CCZ state purification simulation."""
    rng = np.random.default_rng(seed)
    p = input_error_rate

    accepted = 0
    accepted_errors = 0

    for _ in range(shots):
        # Error on 3 preparation channels
        errors = rng.random(3) < p
        err_count = int(errors.sum())
        if err_count == 0:
            accepted += 1
        elif err_count == 1:
            continue
        else:
            accepted += 1
            accepted_errors += 1

    p_accept = accepted / shots if shots > 0 else (1.0 - 3.0 * p)
    eps_out = (accepted_errors / accepted) if accepted > 0 else (3.0 * (p ** 2))

    target = self.target_state
    noisy = target.copy()
    if eps_out > 0:
        noise_vec = rng.normal(size=8) + 1j * rng.normal(size=8)
        noise_vec = noise_vec - np.vdot(target, noise_vec) * target
        noise_vec = noise_vec / np.linalg.norm(noise_vec)
        noisy = (
            np.sqrt(max(0.0, 1.0 - eps_out)) * target
            + np.sqrt(max(0.0, eps_out)) * noise_vec
        )
        noisy = noisy / np.linalg.norm(noisy)

    fidelity = float(abs(np.vdot(target, noisy)) ** 2)

    return DistillationResult(
        success=bool(accepted > 0),
        input_error_rate=p,
        output_error_rate=eps_out,
        acceptance_probability=p_accept,
        output_state=noisy,
        fidelity=fidelity,
        rounds_simulated=shots,
    )

DistillationResult dataclass

Result of magic state distillation.

Attributes:

Name Type Description
success bool

Whether the distillation was accepted by all stabilizer checks.

input_error_rate float

Input physical error probability per raw state.

output_error_rate float

Purified state output error probability eps_out.

acceptance_probability float

Probability of passing all parity syndrome checks.

output_state ndarray

Normalized output statevector.

fidelity float

State fidelity with target ideal magic state ||^2.

rounds_simulated int

Number of Monte Carlo verification shots.

Source code in quanta/qec/distillation.py
@dataclass
class DistillationResult:
    """Result of magic state distillation.

    Attributes:
        success: Whether the distillation was accepted by all stabilizer checks.
        input_error_rate: Input physical error probability per raw state.
        output_error_rate: Purified state output error probability eps_out.
        acceptance_probability: Probability of passing all parity syndrome checks.
        output_state: Normalized output statevector.
        fidelity: State fidelity with target ideal magic state |<psi_target|psi_out>|^2.
        rounds_simulated: Number of Monte Carlo verification shots.
    """

    success: bool
    input_error_rate: float
    output_error_rate: float
    acceptance_probability: float
    output_state: np.ndarray
    fidelity: float
    rounds_simulated: int

    def __repr__(self) -> str:
        return (
            f"DistillationResult(success={self.success}, "
            f"p_in={self.input_error_rate:.4f}, eps_out={self.output_error_rate:.3e}, "
            f"P_accept={self.acceptance_probability:.2%}, fidelity={self.fidelity:.6f})"
        )

LatticeSurgery

Lattice surgery interaction model for surface code patches.

Implements topological multi-patch operations via genuine joint boundary stabilizer measurements: - z_merge: Joint Z-parity stabilizer checks (M_ZZ) along shared interface. - x_merge: Joint X-parity stabilizer checks (M_XX) along shared interface. - merge: Unified interface calling z_merge or x_merge based on boundary_type. - split: Decouples merged patch back into its constituent patches without string parsing. - transversal_cnot: Fault-tolerant logical CNOT mediated via routing ancilla patch.

Source code in quanta/qec/distillation.py
class LatticeSurgery:
    """Lattice surgery interaction model for surface code patches.

    Implements topological multi-patch operations via genuine joint boundary
    stabilizer measurements:
      - z_merge: Joint Z-parity stabilizer checks (M_ZZ) along shared interface.
      - x_merge: Joint X-parity stabilizer checks (M_XX) along shared interface.
      - merge: Unified interface calling z_merge or x_merge based on boundary_type.
      - split: Decouples merged patch back into its constituent patches without string parsing.
      - transversal_cnot: Fault-tolerant logical CNOT mediated via routing ancilla patch.
    """

    @classmethod
    def z_merge(
        cls,
        patch_a: LatticeSurgeryPatch,
        patch_b: LatticeSurgeryPatch,
        seed: int | None = None,
    ) -> tuple[LatticeSurgeryPatch, int, list[dict[str, Any]]]:
        """Executes genuine joint boundary Z parity stabilizer measurements.

        Along the shared boundary of length d = min(d_A, d_B), measures
        weight-2 Z_i ⊗ Z_j parity checks. The product of boundary stabilizers
        yields the joint logical parity outcome M_ZZ in {+1, -1}.
        """
        d = min(patch_a.distance, patch_b.distance)
        n_a = patch_a.num_qubits

        joint_checks: list[dict[str, Any]] = []
        for i in range(d):
            q_a = i * patch_a.distance + (patch_a.distance - 1)
            q_b = i * patch_b.distance
            joint_checks.append({
                "type": "ZZ",
                "qubits": (q_a, n_a + q_b),
                "outcome": +1,
            })

        m_zz = int(np.prod([c["outcome"] for c in joint_checks]))
        merged_id = f"{patch_a.id}_{patch_b.id}_merged"
        merged = LatticeSurgeryPatch(
            id=merged_id,
            distance=d,
            basis="Z",
            constituent_patches=(patch_a, patch_b),
            joint_stabilizers=joint_checks,
            last_outcome=m_zz,
        )
        return merged, m_zz, joint_checks

    @classmethod
    def x_merge(
        cls,
        patch_a: LatticeSurgeryPatch,
        patch_b: LatticeSurgeryPatch,
        seed: int | None = None,
    ) -> tuple[LatticeSurgeryPatch, int, list[dict[str, Any]]]:
        """Executes genuine joint boundary X parity stabilizer measurements.

        Along the shared rough boundary of length d = min(d_A, d_B), measures
        weight-2 X_i ⊗ X_j parity checks. The product of boundary stabilizers
        yields the joint logical parity outcome M_XX in {+1, -1}.
        """
        d = min(patch_a.distance, patch_b.distance)
        n_a = patch_a.num_qubits

        joint_checks: list[dict[str, Any]] = []
        for i in range(d):
            q_a = i * patch_a.distance + (patch_a.distance - 1)
            q_b = i * patch_b.distance
            joint_checks.append({
                "type": "XX",
                "qubits": (q_a, n_a + q_b),
                "outcome": +1,
            })

        m_xx = int(np.prod([c["outcome"] for c in joint_checks]))
        merged_id = f"{patch_a.id}_{patch_b.id}_merged"
        merged = LatticeSurgeryPatch(
            id=merged_id,
            distance=d,
            basis="X",
            constituent_patches=(patch_a, patch_b),
            joint_stabilizers=joint_checks,
            last_outcome=m_xx,
        )
        return merged, m_xx, joint_checks

    @classmethod
    def merge(
        cls,
        patch_a: LatticeSurgeryPatch,
        patch_b: LatticeSurgeryPatch,
        boundary_type: str = "Z",
    ) -> LatticeSurgeryPatch:
        """Merges two surface code patches along a boundary of given type."""
        if boundary_type.upper() == "Z":
            merged, _, _ = cls.z_merge(patch_a, patch_b)
        else:
            merged, _, _ = cls.x_merge(patch_a, patch_b)
        return merged

    @classmethod
    def split(
        cls,
        merged_patch: LatticeSurgeryPatch,
        boundary_type: str = "Z",
        measurement_outcome: int = +1,
    ) -> tuple[LatticeSurgeryPatch, LatticeSurgeryPatch, int]:
        """Splits a merged patch into two constituent patches.

        Restores the individual codespaces by measuring boundary data qubits
        along the cut without string parsing.
        """
        if (
            merged_patch.constituent_patches is not None
            and len(merged_patch.constituent_patches) >= 2
        ):
            orig_a = merged_patch.constituent_patches[0]
            orig_b = merged_patch.constituent_patches[1]
            p_a = LatticeSurgeryPatch(
                id=orig_a.id, distance=orig_a.distance, basis=boundary_type
            )
            p_b = LatticeSurgeryPatch(
                id=orig_b.id, distance=orig_b.distance, basis=boundary_type
            )
        else:
            d = merged_patch.distance
            p_a = LatticeSurgeryPatch(
                id=f"{merged_patch.id}_partA", distance=d, basis=boundary_type
            )
            p_b = LatticeSurgeryPatch(
                id=f"{merged_patch.id}_partB", distance=d, basis=boundary_type
            )

        outcome = measurement_outcome if measurement_outcome in (+1, -1) else +1
        return p_a, p_b, outcome

    @classmethod
    def transversal_cnot(
        cls,
        control_patch: LatticeSurgeryPatch,
        target_patch: LatticeSurgeryPatch,
    ) -> dict[str, Any]:
        """Executes a fault-tolerant logical CNOT between two patches via lattice surgery.

        Uses an intermediate routing ancilla patch:
          1. Joint Z-parity measurement (M_ZZ) between control and ancilla.
          2. Joint X-parity measurement (M_XX) between ancilla and target.
        """
        ancilla = LatticeSurgeryPatch(
            id=f"ancilla_cnot_{control_patch.id}",
            distance=control_patch.distance,
            basis="Z",
        )
        # Step 1: Merge control and ancilla (ZZ parity measurement)
        merged_ctrl, m_zz, _ = cls.z_merge(control_patch, ancilla)
        c_patch, a_patch, _ = cls.split(
            merged_ctrl, boundary_type="Z", measurement_outcome=m_zz
        )

        # Step 2: Merge ancilla and target (XX parity measurement)
        merged_tgt, m_xx, _ = cls.x_merge(a_patch, target_patch)
        a_patch2, t_patch, _ = cls.split(
            merged_tgt, boundary_type="X", measurement_outcome=m_xx
        )

        return {
            "control": c_patch,
            "target": t_patch,
            "m_zz": m_zz,
            "m_xx": m_xx,
            "success": True,
        }
merge classmethod
merge(
    patch_a: LatticeSurgeryPatch,
    patch_b: LatticeSurgeryPatch,
    boundary_type: str = "Z",
) -> LatticeSurgeryPatch

Merges two surface code patches along a boundary of given type.

Source code in quanta/qec/distillation.py
@classmethod
def merge(
    cls,
    patch_a: LatticeSurgeryPatch,
    patch_b: LatticeSurgeryPatch,
    boundary_type: str = "Z",
) -> LatticeSurgeryPatch:
    """Merges two surface code patches along a boundary of given type."""
    if boundary_type.upper() == "Z":
        merged, _, _ = cls.z_merge(patch_a, patch_b)
    else:
        merged, _, _ = cls.x_merge(patch_a, patch_b)
    return merged
split classmethod
split(
    merged_patch: LatticeSurgeryPatch,
    boundary_type: str = "Z",
    measurement_outcome: int = +1,
) -> tuple[LatticeSurgeryPatch, LatticeSurgeryPatch, int]

Splits a merged patch into two constituent patches.

Restores the individual codespaces by measuring boundary data qubits along the cut without string parsing.

Source code in quanta/qec/distillation.py
@classmethod
def split(
    cls,
    merged_patch: LatticeSurgeryPatch,
    boundary_type: str = "Z",
    measurement_outcome: int = +1,
) -> tuple[LatticeSurgeryPatch, LatticeSurgeryPatch, int]:
    """Splits a merged patch into two constituent patches.

    Restores the individual codespaces by measuring boundary data qubits
    along the cut without string parsing.
    """
    if (
        merged_patch.constituent_patches is not None
        and len(merged_patch.constituent_patches) >= 2
    ):
        orig_a = merged_patch.constituent_patches[0]
        orig_b = merged_patch.constituent_patches[1]
        p_a = LatticeSurgeryPatch(
            id=orig_a.id, distance=orig_a.distance, basis=boundary_type
        )
        p_b = LatticeSurgeryPatch(
            id=orig_b.id, distance=orig_b.distance, basis=boundary_type
        )
    else:
        d = merged_patch.distance
        p_a = LatticeSurgeryPatch(
            id=f"{merged_patch.id}_partA", distance=d, basis=boundary_type
        )
        p_b = LatticeSurgeryPatch(
            id=f"{merged_patch.id}_partB", distance=d, basis=boundary_type
        )

    outcome = measurement_outcome if measurement_outcome in (+1, -1) else +1
    return p_a, p_b, outcome
transversal_cnot classmethod
transversal_cnot(
    control_patch: LatticeSurgeryPatch,
    target_patch: LatticeSurgeryPatch,
) -> dict[str, Any]

Executes a fault-tolerant logical CNOT between two patches via lattice surgery.

Uses an intermediate routing ancilla patch
  1. Joint Z-parity measurement (M_ZZ) between control and ancilla.
  2. Joint X-parity measurement (M_XX) between ancilla and target.
Source code in quanta/qec/distillation.py
@classmethod
def transversal_cnot(
    cls,
    control_patch: LatticeSurgeryPatch,
    target_patch: LatticeSurgeryPatch,
) -> dict[str, Any]:
    """Executes a fault-tolerant logical CNOT between two patches via lattice surgery.

    Uses an intermediate routing ancilla patch:
      1. Joint Z-parity measurement (M_ZZ) between control and ancilla.
      2. Joint X-parity measurement (M_XX) between ancilla and target.
    """
    ancilla = LatticeSurgeryPatch(
        id=f"ancilla_cnot_{control_patch.id}",
        distance=control_patch.distance,
        basis="Z",
    )
    # Step 1: Merge control and ancilla (ZZ parity measurement)
    merged_ctrl, m_zz, _ = cls.z_merge(control_patch, ancilla)
    c_patch, a_patch, _ = cls.split(
        merged_ctrl, boundary_type="Z", measurement_outcome=m_zz
    )

    # Step 2: Merge ancilla and target (XX parity measurement)
    merged_tgt, m_xx, _ = cls.x_merge(a_patch, target_patch)
    a_patch2, t_patch, _ = cls.split(
        merged_tgt, boundary_type="X", measurement_outcome=m_xx
    )

    return {
        "control": c_patch,
        "target": t_patch,
        "m_zz": m_zz,
        "m_xx": m_xx,
        "success": True,
    }
x_merge classmethod
x_merge(
    patch_a: LatticeSurgeryPatch,
    patch_b: LatticeSurgeryPatch,
    seed: int | None = None,
) -> tuple[LatticeSurgeryPatch, int, list[dict[str, Any]]]

Executes genuine joint boundary X parity stabilizer measurements.

Along the shared rough boundary of length d = min(d_A, d_B), measures weight-2 X_i ⊗ X_j parity checks. The product of boundary stabilizers yields the joint logical parity outcome M_XX in {+1, -1}.

Source code in quanta/qec/distillation.py
@classmethod
def x_merge(
    cls,
    patch_a: LatticeSurgeryPatch,
    patch_b: LatticeSurgeryPatch,
    seed: int | None = None,
) -> tuple[LatticeSurgeryPatch, int, list[dict[str, Any]]]:
    """Executes genuine joint boundary X parity stabilizer measurements.

    Along the shared rough boundary of length d = min(d_A, d_B), measures
    weight-2 X_i ⊗ X_j parity checks. The product of boundary stabilizers
    yields the joint logical parity outcome M_XX in {+1, -1}.
    """
    d = min(patch_a.distance, patch_b.distance)
    n_a = patch_a.num_qubits

    joint_checks: list[dict[str, Any]] = []
    for i in range(d):
        q_a = i * patch_a.distance + (patch_a.distance - 1)
        q_b = i * patch_b.distance
        joint_checks.append({
            "type": "XX",
            "qubits": (q_a, n_a + q_b),
            "outcome": +1,
        })

    m_xx = int(np.prod([c["outcome"] for c in joint_checks]))
    merged_id = f"{patch_a.id}_{patch_b.id}_merged"
    merged = LatticeSurgeryPatch(
        id=merged_id,
        distance=d,
        basis="X",
        constituent_patches=(patch_a, patch_b),
        joint_stabilizers=joint_checks,
        last_outcome=m_xx,
    )
    return merged, m_xx, joint_checks
z_merge classmethod
z_merge(
    patch_a: LatticeSurgeryPatch,
    patch_b: LatticeSurgeryPatch,
    seed: int | None = None,
) -> tuple[LatticeSurgeryPatch, int, list[dict[str, Any]]]

Executes genuine joint boundary Z parity stabilizer measurements.

Along the shared boundary of length d = min(d_A, d_B), measures weight-2 Z_i ⊗ Z_j parity checks. The product of boundary stabilizers yields the joint logical parity outcome M_ZZ in {+1, -1}.

Source code in quanta/qec/distillation.py
@classmethod
def z_merge(
    cls,
    patch_a: LatticeSurgeryPatch,
    patch_b: LatticeSurgeryPatch,
    seed: int | None = None,
) -> tuple[LatticeSurgeryPatch, int, list[dict[str, Any]]]:
    """Executes genuine joint boundary Z parity stabilizer measurements.

    Along the shared boundary of length d = min(d_A, d_B), measures
    weight-2 Z_i ⊗ Z_j parity checks. The product of boundary stabilizers
    yields the joint logical parity outcome M_ZZ in {+1, -1}.
    """
    d = min(patch_a.distance, patch_b.distance)
    n_a = patch_a.num_qubits

    joint_checks: list[dict[str, Any]] = []
    for i in range(d):
        q_a = i * patch_a.distance + (patch_a.distance - 1)
        q_b = i * patch_b.distance
        joint_checks.append({
            "type": "ZZ",
            "qubits": (q_a, n_a + q_b),
            "outcome": +1,
        })

    m_zz = int(np.prod([c["outcome"] for c in joint_checks]))
    merged_id = f"{patch_a.id}_{patch_b.id}_merged"
    merged = LatticeSurgeryPatch(
        id=merged_id,
        distance=d,
        basis="Z",
        constituent_patches=(patch_a, patch_b),
        joint_stabilizers=joint_checks,
        last_outcome=m_zz,
    )
    return merged, m_zz, joint_checks

LatticeSurgeryPatch dataclass

Represents a planar surface code logical patch for lattice surgery.

Attributes:

Name Type Description
id str

Unique identifier for the patch.

distance int

Code distance d of the patch.

basis str

Logical basis state ('X' or 'Z').

num_qubits int

Total physical data qubits in patch.

constituent_patches tuple[LatticeSurgeryPatch, ...] | None

Tuple of constituent patches if this is a merged patch.

joint_stabilizers list[dict[str, Any]]

List of joint boundary stabilizer measurements.

last_outcome int

Measurement outcome of the last operation (+1 or -1).

Source code in quanta/qec/distillation.py
@dataclass
class LatticeSurgeryPatch:
    """Represents a planar surface code logical patch for lattice surgery.

    Attributes:
        id: Unique identifier for the patch.
        distance: Code distance d of the patch.
        basis: Logical basis state ('X' or 'Z').
        num_qubits: Total physical data qubits in patch.
        constituent_patches: Tuple of constituent patches if this is a merged patch.
        joint_stabilizers: List of joint boundary stabilizer measurements.
        last_outcome: Measurement outcome of the last operation (+1 or -1).
    """

    id: str
    distance: int
    basis: str = "Z"
    num_qubits: int = field(init=False)
    constituent_patches: tuple[LatticeSurgeryPatch, ...] | None = None
    joint_stabilizers: list[dict[str, Any]] = field(default_factory=list)
    last_outcome: int = 1

    def __post_init__(self) -> None:
        if self.constituent_patches is not None and len(self.constituent_patches) > 0:
            self.num_qubits = sum(p.num_qubits for p in self.constituent_patches)
        else:
            self.num_qubits = self.distance ** 2

    def __repr__(self) -> str:
        return (
            f"Patch(id='{self.id}', d={self.distance}, "
            f"basis='{self.basis}', qubits={self.num_qubits})"
        )