REPEAT block handling in is_clifford (#128)

## Summary
- `Circuit.is_clifford` and the half-π parametric-rotation expansion now
recurse into `REPEAT` blocks. Non-Clifford gates inside loops were
previously misclassified, and half-π rotations inside loops were not
expanded.
- Internal cleanup: rename `_is_clifford` → `is_clifford`,
`_expand_clifford` → `expand_clifford_rotations`, `clifford_expansion` →
`_try_clifford_expansion` to align underscore prefixes with intended
visibility.

Unrelated change:
- Constant phase values are reduced modulo 8 during compilation to keep
phase representations consistent.

Author: rafaelha

CHANGELOG.md

@@ -15,6 +15,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### Added
 - `TPP` and `TPP_DAG` instructions — applies exp(-i pi/8 P) or exp(+i pi/8 P) (up to global phase) for a Pauli product P, i.e., phases the -1 eigenspace of P by exp(i pi/4) or exp(-i pi/4).
+- `Circuit.is_clifford` now supports `REPEAT` blocks.
 
 ## [0.1.3] - 2026-04-13
 

src/tsim/circuit.py

@@ -2,7 +2,6 @@
 
 from __future__ import annotations
 
-from fractions import Fraction
 from typing import (
     TYPE_CHECKING,
     Any,
@@ -23,7 +22,7 @@
 from tsim.core.graph import build_sampling_graph
 from tsim.core.parse import parse_parametric_tag, parse_stim_circuit
 from tsim.noise.dem import get_detector_error_model
-from tsim.utils.clifford import clifford_expansion
+from tsim.utils.clifford import expand_clifford_rotations, is_clifford
 from tsim.utils.diagram import render_pyzx_d3, render_svg
 from tsim.utils.program_text import (
     enriched_stim_error,
@@ -359,20 +358,10 @@ def stim_circuit(self) -> stim.Circuit:
         """Return the underlying stim circuit.
 
         Parametric rotation instructions whose angles are all half-π multiples
-        are expanded into their equivalent Clifford gates.
+        are expanded into their equivalent Clifford gates. ``REPEAT`` blocks are
+        preserved structurally; their bodies are expanded recursively.
         """
-        circ = stim.Circuit()
-        for instr in self._stim_circ:
-            assert not isinstance(instr, stim.CircuitRepeatBlock)
-
-            expansion = clifford_expansion(instr)
-            if expansion is not None:
-                gates, targets = expansion
-                for gate in gates:
-                    circ.append(gate, targets, [])
-            else:
-                circ.append(instr)
-        return circ
+        return expand_clifford_rotations(self._stim_circ)
 
     @property
     def is_clifford(self) -> bool:
@@ -385,35 +374,7 @@ def is_clifford(self) -> bool:
             True if the circuit is a Clifford circuit, otherwise False.
 
         """
-
-        def is_half_pi_multiple(phase: Fraction) -> bool:
-            return phase.denominator <= 2
-
-        for instr in self._stim_circ:
-            assert not isinstance(instr, stim.CircuitRepeatBlock)
-
-            if instr.name in ["S", "S_DAG", "SPP", "SPP_DAG"] and instr.tag == "T":
-                return False
-
-            if instr.name == "I" and instr.tag:
-                result = parse_parametric_tag(instr.tag)
-                if result is None:
-                    return False
-
-                gate_name, params = result
-                if gate_name in ["R_X", "R_Y", "R_Z"]:
-                    if not is_half_pi_multiple(params["theta"]):
-                        return False
-                elif gate_name == "U3":
-                    if not all(
-                        is_half_pi_multiple(params[name])
-                        for name in ("theta", "phi", "lambda")
-                    ):
-                        return False
-                else:
-                    return False
-
-        return True
+        return is_clifford(self._stim_circ)
 
     @property
     def num_measurements(self) -> int:
@@ -823,7 +784,7 @@ def fix_tags(circuit: stim.Circuit) -> stim.Circuit:
                 args = instr.gate_args_copy()
 
                 if name == "I" and tag:
-                    parsed = parse_parametric_tag(tag)
+                    parsed = parse_parametric_tag(instr)
                     if parsed is not None:
                         gate_name, params = parsed
                         if gate_name == "U3":

src/tsim/compile/compile.py

@@ -75,7 +75,7 @@ def _compile_node_phases(
 
     for i, terms in enumerate(terms_per_graph):
         for j, (const_phase, param_bit) in enumerate(terms):
-            phases[i, j] = const_phase
+            phases[i, j] = const_phase % 8
             params[i, j] = param_bit
 
     return NodePhases(
@@ -254,8 +254,8 @@ def _compile_phase_pairs(
 
     for i, terms in enumerate(terms_per_graph):
         for j, (ca, cb, pa, pb) in enumerate(terms):
-            alpha[i, j] = ca
-            beta[i, j] = cb
+            alpha[i, j] = ca % 8
+            beta[i, j] = cb % 8
             alpha_params_arr[i, j] = pa
             beta_params_arr[i, j] = pb
 
@@ -303,7 +303,7 @@ def _compile_prefactor(g_list: list[BaseGraph]) -> ScalarPrefactor:
     )
 
     phase_indices = jnp.array(
-        [int(float(g.scalar.phase) * 4) for g in g_list], dtype=jnp.uint8
+        [int(float(g.scalar.phase) * 4) % 8 for g in g_list], dtype=jnp.uint8
     )
 
     exact_floatfactor = []

src/tsim/core/parse.py

@@ -25,20 +25,38 @@
     u3,
 )
 
+_PARAMETRIC_GATE_PARAMS: dict[str, frozenset[str]] = {
+    "R_X": frozenset({"theta"}),
+    "R_Y": frozenset({"theta"}),
+    "R_Z": frozenset({"theta"}),
+    "U3": frozenset({"theta", "phi", "lambda"}),
+}
+
 
-def parse_parametric_tag(tag: str) -> tuple[str, dict[str, Fraction]] | None:
-    """Parse a parametric gate tag like R_Z(theta=0.3*pi).
+def parse_parametric_tag(
+    instruction: stim.CircuitInstruction,
+) -> tuple[str, dict[str, Fraction]] | None:
+    """Parse the parametric tag on an instruction (e.g. ``I[R_Z(theta=0.3*pi)]``).
 
     Supports gates: R_Z, R_X, R_Y, U3.
 
     Args:
-        tag: The instruction tag to parse, e.g. "R_Z(theta=0.3*pi)" or
-             "U3(theta=0.3*pi, phi=0.24*pi, lambda=0.49*pi)".
+        instruction: The stim instruction whose tag will be parsed.
 
     Returns:
-        Tuple of (gate_name, params_dict) or None if not a valid parametric tag.
+        Tuple of (gate_name, params_dict) when the instruction's tag is a
+        well-formed parametric tag, or ``None`` when the tag is not
+        parametric-looking (no ``name(...)`` shape, or empty).
+
+    Raises:
+        ValueError: When the tag looks parametric (matches ``name(...)``) but is
+            malformed: a parameter value does not parse, the gate name is unknown,
+            or the parameter keys do not match the expected set for the gate.
 
     """
+    tag = instruction.tag
+    err_prefix = f"Could not parse instruction {str(instruction)!r}"
+
     match = re.match(r"^(\w+)\((.*)\)$", tag)
     if not match:
         return None
@@ -56,11 +74,20 @@ def parse_parametric_tag(tag: str) -> tuple[str, dict[str, Fraction]] | None:
             r"^(\w+)=([-+]?(?:\d+\.?\d*|\.\d+)(?:[eE][-+]?\d+)?)\*pi$", param
         )
         if not param_match:
-            return None
+            raise ValueError(f"{err_prefix}. Malformed parametric tag {tag!r}")
         param_name = param_match.group(1)
         value = Fraction(param_match.group(2))
         params[param_name] = value
 
+    expected = _PARAMETRIC_GATE_PARAMS.get(gate_name)
+    if expected is None:
+        raise ValueError(f"{err_prefix}. Unknown parametric gate {gate_name!r}")
+    if params.keys() != expected:
+        raise ValueError(
+            f"{err_prefix}. Parametric tag {tag!r} has parameters "
+            f"{sorted(params)}, expected {sorted(expected)}"
+        )
+
     return gate_name, params
 
 
@@ -166,7 +193,7 @@ def parse_stim_circuit(
 
         # Handle parametric gates via tags (e.g., I with tag "R_Z(theta=0.3*pi)")
         if name == "I" and instruction.tag:
-            result = parse_parametric_tag(instruction.tag)
+            result = parse_parametric_tag(instruction)
             if result is not None:
                 gate_name, params = result
                 targets = [t.value for t in instruction.targets_copy()]

src/tsim/external/vec_sim/vec_sampler.py

@@ -81,7 +81,7 @@ def sample_circuit_with_vec_sim_return_data(
             for q in inst.targets_copy():
                 sim.do_t_dag(q.qubit_value)
         elif inst.name == "I" and inst.tag:
-            result = parse_parametric_tag(inst.tag)
+            result = parse_parametric_tag(inst)
             if result is not None:
                 gate_name, params = result
                 for q in inst.targets_copy():

src/tsim/utils/clifford.py

@@ -101,7 +101,70 @@ def parametric_to_clifford_gates(
     return None
 
 
-def clifford_expansion(
+def is_clifford(source: stim.Circuit) -> bool:
+    """Return True iff every instruction in ``source`` is Clifford.
+
+    Recurses into ``REPEAT`` block bodies.
+    """
+
+    def is_half_pi_multiple(phase: Fraction) -> bool:
+        return phase.denominator <= 2
+
+    for instr in source:
+        if isinstance(instr, stim.CircuitRepeatBlock):
+            if not is_clifford(instr.body_copy()):
+                return False
+            continue
+
+        if instr.name in ["S", "S_DAG", "SPP", "SPP_DAG"] and instr.tag == "T":
+            return False
+
+        if instr.name == "I" and instr.tag:
+            result = parse_parametric_tag(instr)
+            if result is None:
+                return False
+
+            gate_name, params = result
+            if gate_name in ["R_X", "R_Y", "R_Z"]:
+                if not is_half_pi_multiple(params["theta"]):
+                    return False
+            elif gate_name == "U3":
+                if not all(
+                    is_half_pi_multiple(params[name])
+                    for name in ("theta", "phi", "lambda")
+                ):
+                    return False
+            else:
+                return False
+
+    return True
+
+
+def expand_clifford_rotations(source: stim.Circuit) -> stim.Circuit:
+    """Return ``source`` with half-π parametric rotations expanded to Clifford gates.
+
+    ``REPEAT`` blocks are preserved structurally and expanded recursively.
+    """
+    out = stim.Circuit()
+    for instr in source:
+        if isinstance(instr, stim.CircuitRepeatBlock):
+            out.append(
+                stim.CircuitRepeatBlock(
+                    instr.repeat_count, expand_clifford_rotations(instr.body_copy())
+                )
+            )
+            continue
+        expansion = _try_clifford_expansion(instr)
+        if expansion is not None:
+            gates, targets = expansion
+            for gate in gates:
+                out.append(gate, targets, [])
+        else:
+            out.append(instr)
+    return out
+
+
+def _try_clifford_expansion(
     instr: stim.CircuitInstruction,
 ) -> tuple[list[str], list[int]] | None:
     """Try to expand a tagged ``I`` instruction into equivalent Clifford gates.
@@ -115,7 +178,7 @@ def clifford_expansion(
     if instr.name != "I" or not instr.tag:
         return None
 
-    parsed = parse_parametric_tag(instr.tag)
+    parsed = parse_parametric_tag(instr)
     if parsed is None:
         return None
 

test/unit/compile/test_compile.py

@@ -1,5 +1,8 @@
+from fractions import Fraction
+
 import jax.numpy as jnp
 import numpy as np
+from pyzx_param.graph.base import BaseGraph
 from pyzx_param.graph.graph_s import GraphS
 
 from tsim.compile.compile import compile_scalar_graphs
@@ -40,3 +43,36 @@ def test_evaluate_empty_with_params(self):
         result = evaluate(compiled, param_vals)
         assert result.shape == (3,)
         np.testing.assert_array_equal(np.asarray(result), np.zeros(3, dtype=complex))
+
+
+class TestCompilePrefactorPhase:
+    """Phase indices must be reduced modulo 8 so the uint8 storage and
+    UNIT_PHASES lookup in compile/terms.py see canonical values 0..7."""
+
+    def _scalar_graph(self, phase: Fraction) -> GraphS:
+        g = GraphS()
+        g.scalar.phase = phase
+        return g
+
+    def test_negative_phase_reduced_into_range(self):
+        # phase = -1/4 -> int(-1.0) -> -1 -> -1 % 8 == 7
+        compiled = compile_scalar_graphs([self._scalar_graph(Fraction(-1, 4))], [])
+        assert int(compiled.prefactor.phase_indices[0]) == 7
+
+    def test_phase_above_two_reduced_into_range(self):
+        # phase = 9/4 -> int(9.0) -> 9 -> 9 % 8 == 1
+        compiled = compile_scalar_graphs([self._scalar_graph(Fraction(9, 4))], [])
+        assert int(compiled.prefactor.phase_indices[0]) == 1
+
+    def test_in_range_phase_unchanged(self):
+        # phase = 5/4 -> int(5.0) -> 5 -> 5 % 8 == 5
+        compiled = compile_scalar_graphs([self._scalar_graph(Fraction(5, 4))], [])
+        assert int(compiled.prefactor.phase_indices[0]) == 5
+
+    def test_phase_indices_stay_within_unit_phase_table(self):
+        graphs: list[BaseGraph] = [
+            self._scalar_graph(Fraction(p, 4)) for p in (-7, -1, 0, 3, 7, 11)
+        ]
+        compiled = compile_scalar_graphs(graphs, [])
+        indices = np.asarray(compiled.prefactor.phase_indices)
+        assert np.all((indices >= 0) & (indices < 8))