errors.py

"""
errors.py — Pauli error models for D(S3).

Each operator is an independent Bernoulli channel firing with
probability p. Errors accumulate in the Pauli frame; the decorated
syndromes in lattice.measure_syndromes() read the combined state.

Per edge (r, c, d) per call to apply():

  Qutrit (all edges, including ungauged phase=1):
    X   fires with prob p  → qutrit_X += 1  (mod 3)
    X†  fires with prob p  → qutrit_X += 2  (mod 3)
    Z   fires with prob p  → qutrit_Z += 1  (mod 3)
    Z†  fires with prob p  → qutrit_Z += 2  (mod 3)

  Qubit (D(S3) edges only, phase=0):
    sigma^X  fires with prob p  → qubit_X ^= 1
    sigma^Z  fires with prob p  → qubit_Z ^= 1

All six channels are independent, so multiple operators can fire on
the same edge in one step; the Pauli frame accumulates them correctly
mod 3 (qutrit) or mod 2 (qubit). In particular:

  X then X†  → net qutrit_X += 3 ≡ 0   (identity)    ✓
  X then X   → net qutrit_X += 2 = X†                 ✓
  sigma^X and sigma^Z both fire → effective sigma^Y    ✓

The per-operator probability p is the same for all six channels,
matching the model discussed for the partition function construction.
"""

import numpy as np
from lattice import DS3Lattice


class DepolarizingError:
    """
    Independent Pauli channels on each edge, each firing with prob p.

    Parameters
    ----------
    p      : single-operator error probability (same for all 6 channels)
    p_meas : measurement error probability (defaults to p if not given)
    """

    def __init__(self, p: float, p_meas: float = None):
        self.p = p
        self.p_meas = p if p_meas is None else p_meas

    def apply(self, lattice: DS3Lattice, rng: np.random.Generator):
        L = lattice.L
        p = self.p
        if p <= 0:
            return

        for d in range(2):

            # ----------------------------------------------------------------
            # Qutrit channels — applied to ALL edges regardless of phase.
            # Ungauged (phase=1) edges still carry a live qutrit; errors on
            # them are picked up by compute_z3_syndrome_arrays() in the patch.
            # ----------------------------------------------------------------

            # X:  qutrit_X += 1 mod 3
            mask = rng.random((L, L)) < p
            lattice.qutrit_X[:, :, d] = (
                lattice.qutrit_X[:, :, d] + mask) % 3

            # X†: qutrit_X += 2 mod 3
            mask = rng.random((L, L)) < p
            lattice.qutrit_X[:, :, d] = (
                lattice.qutrit_X[:, :, d] + 2 * mask) % 3

            # Z:  qutrit_Z += 1 mod 3
            mask = rng.random((L, L)) < p
            lattice.qutrit_Z[:, :, d] = (
                lattice.qutrit_Z[:, :, d] + mask) % 3

            # Z†: qutrit_Z += 2 mod 3
            mask = rng.random((L, L)) < p
            lattice.qutrit_Z[:, :, d] = (
                lattice.qutrit_Z[:, :, d] + 2 * mask) % 3

            # ----------------------------------------------------------------
            # Qubit channels — only on D(S3) edges (phase=0).
            # Ungauged edges (phase=1) have no qubit subsystem.
            # ----------------------------------------------------------------
            gauged = (lattice.phase[:, :, d] == 0)

            # sigma^X: qubit_X ^= 1
            mask = (rng.random((L, L)) < p) & gauged
            lattice.qubit_X[:, :, d] ^= mask.astype(np.int8)

            # sigma^Z: qubit_Z ^= 1
            mask = (rng.random((L, L)) < p) & gauged
            lattice.qubit_Z[:, :, d] ^= mask.astype(np.int8)