Surface-Memory Reference Comparison
gadgets/04_surface_memory_reference.py
extends the physical-memory introduction with
individual gadget-model exports, explicit data depolarization, and a comparison
against Stim’s generated rotated-memory reference. It also offers physical
circuit and independent-DEM sampling.
Start with gadget examples 00 and 01 for the basic verification and model workflow. Optional examples 02 and 03, covered by the Stim-emission and custom-noise tutorials, instead use ideal MPP stabilizer measurements with explicit data noise. This advanced example uses physical ancillas, four interleaved N/Z CNOT layers, reset and measurement noise, and final data-parity detectors.
Run the example
With CUDA-Q Logical and cudaq-edem installed, install the optional reference
simulator and run from the repository root:
python3 -m pip install 'stim>=1.16,<2'
python3 preview/logical/examples/gadgets/04_surface_memory_reference.py
The default output is:
Verified prep, SEC, readout: ideal logical actions and detector contracts
d=3, rounds=3: 24 detectors, 1 observable, 1307 original fault indicators
Artifacts: surface-memory-circuit
219 DEM error terms match the circuit and Stim reference (atol=1e-12)
Edit the parameter cell, or run the # %% cells interactively:
distance = 3
rounds = 3 # Includes one extraction in preparation and one in readout; >= 2.
after_clifford_depolarization = 0.001
before_measure_flip_probability = 0.001
after_reset_flip_probability = 0.001
data_depolarization = 0.001
output = Path("surface-memory-circuit")
compare_with_stim = is_stim_installed() # pip install 'stim>=1.16,<2'.
shots = 0
seed = 7
sampling_mode = "circuit" # Or `dem`: sample the exactly factorized DEM.
Set compare_with_stim = False and shots = 0 to export the circuit and packed
models without importing Stim. Change rounds to 10 to exercise model reuse
across a longer memory experiment.
Select the concrete gadgets
Select preparation, syndrome extraction (SEC), and destructive Z readout, then connect their code boundaries:
implementation = surface.interleaved_nz_definitions(distance)
prep = surface.gadget_contract(implementation.prepare_zero)
sec = surface.gadget_contract(implementation.syndrome_extraction)
readout = surface.gadget_contract(implementation.measure_z_with_syndrome)
# Both endpoints include an SEC. The example supplies the remaining rounds.
contracts = (prep, *((sec,) * (rounds - 2)), readout)
The interleaved definition set contains prepare_zero, syndrome_extraction,
and measure_z_with_syndrome, together with their code and physical-qubit
coordinates. surface.gadget_contract(gadget) supplies each detector contract.
The separate surface.definitions(distance) collection uses ideal-MPP
syndrome rounds.
Each physical gadget includes one extraction. At rounds = 3, the sequence is:
prepare_zero[SEC] -> syndrome_extraction -> measure_z_with_syndrome[SEC + data MZ]
This experiment requires at least two rounds because its preparation and
readout both contain an SEC. The ideal-MPP examples allow one round because
their direct logical readout contains no SEC. Interleaved readout measures all
data qubits, compares their Z-check parities with its extraction records, and
forms logical Z from the data-measurement parity. At distance three, the circuit
uses nine data qubits and eight check ancillas. With r extraction rounds, it
declares r * (distance**2 - 1) detectors and one logical observable.
Verify the ideal behavior
Before adding noise, verify each reusable gadget once:
verification = {
name: cql.gadgets.verify_gadget(contract)
for name, contract in (("prep", prep), ("SEC", sec), ("readout", readout))
}
assert all(receipt.verified for receipt in verification.values())
print("Verified prep, SEC, readout: "
"ideal logical actions and detector contracts")
The receipts establish logical +Z preparation, an identity syndrome round, and destructive logical Z readout, including its returned measurement parity. They also establish deterministic and complete detector contracts. These are ideal stabilizer checks; noisy circuit distance is checked separately. See Verify a contract and the logical action for the supported claims and failure receipts.
Project the realization and choose physical noise
Compile and project each reusable gadget to choose physical fault sites:
builds = {
contract.gadget: cql.compile(contract) for contract in (prep, sec, readout)
}
projections = {
contract.gadget:
cql.compiler.project_clifford(builds[contract.gadget],
gadget=contract.gadget.name)
for contract in (prep, sec, readout)
}
StandardDepolarizingNoise applies depolarization after Clifford gates,
measurement-flip errors before measurements, and reset errors after resets.
This example also places data depolarization at the recipe’s declared
before_syndrome site:
noise = cql.compiler.StandardDepolarizingNoise(
after_clifford_depolarization=after_clifford_depolarization,
before_measure_flip_probability=before_measure_flip_probability,
after_reset_flip_probability=after_reset_flip_probability)
data_sites = {
contract.gadget:
surface.data_depol(contract, projections[contract.gadget],
data_depolarization)
for contract in (prep, sec, readout)
}
placements = {
index: data_sites[contract.gadget]
for index, contract in enumerate(contracts)
}
surface.data_depol(contract, projection, probability) returns ordinary
Depolarize1Placement objects. The helper checks the supplied projection
against the registered realization without compiling it again. It rejects
unknown sites, mismatched projections, and gadgets without the requested site.
The example chooses the probabilities and noisy occurrences.
These data placements follow initialization resets and precede syndrome
extraction. Readout’s data noise precedes its SEC, rather than its final
measurements. A direct readout without SEC has no before_syndrome site; the
compiler does not infer extraction boundaries from gadget names.
Emit the circuit and reuse gadget models
The same contracts and physical noise placements feed circuit emission and EDEM derivation:
coordinates = "\n".join(f"QUBIT_COORDS({x}, {y}) {q}"
for q, (x, y) in enumerate(implementation.coordinates))
text = coordinates + "\n" + cql.compiler.emit_stim(
contracts, noise=noise, depolarizations=placements)
models = {}
for contract in (prep, sec, readout):
projection = projections[contract.gadget]
faults = cql.analysis.standard_depolarizing_faults(
projection, noise, depolarizations=data_sites[contract.gadget])
models[contract.gadget] = cql.analysis.edem_from_gadget(contract.gadget,
contract=contract,
faults=faults)
occurrences = tuple(models[contract.gadget] for contract in contracts)
composed = edem.compose_chain(occurrences)
projected = edem.project_DEM(composed)
edem_from_gadget accepts each gadget and its explicit contract directly. The
projections above are used to choose faults; callers do not need a compiled
Build for EDEM derivation.
An EDEM retains fault effects on an open gadget’s code boundaries, allowing a
fault to influence later gadgets’ detectors. The example derives preparation,
SEC, and readout models once, then reuses the SEC model at every intermediate
round. compose_chain preserves distinct fault indicators and exact
categorical laws for each occurrence. Preparation and destructive readout
close the chain’s code interfaces, allowing project_DEM to produce a DEM.
At rounds = 10, the eight intermediate SEC occurrences still share one
derived SEC model. The composed model gains each occurrence’s detectors and
fault indicators, without repeating the gadget derivation. The
EDEM tutorial explains these boundary coordinates and the
shorter derive_edem_chain API used in example 01.
Export the original models
Save the circuit, projected DEM, and reusable gadget EDEMs before any noise conversion:
output = Path(output)
output.mkdir(parents=True, exist_ok=True)
(output / "memory.stim").write_text(text)
(output / "memory.edemb").write_bytes(projected.to_bytes())
for name, contract in (("prep", prep), ("sec", sec), ("mz", readout)):
(output / f"{name}.edemb").write_bytes(models[contract.gadget].to_bytes())
print(f"d={distance}, rounds={rounds}: "
f"{sum(projected.detectors)} detectors, "
f"{sum(projected.observables)} observable, "
f"{sum(f.size for f in projected.faults)} original fault indicators")
print(f"Artifacts: {output}")
The output directory contains:
File |
Content |
|---|---|
|
Emitted physical circuit with noise, detectors, and logical Z |
|
Closed projected DEM with original categorical fault laws |
|
Reusable gadget EDEMs with their code interfaces |
|
Stim’s generated rotated-memory circuit |
|
DEM obtained from the emitted circuit |
|
Formatted, exactly converted and coalesced gadget-derived DEM |
The last three files are written when compare_with_stim = True. The packed
models retain their original fault laws and response columns. The Surface
recipe supplies the gadgets and contracts directly, so another composer can
reuse them without importing this example.
Compare with Stim and choose a sampler
The reference comparison explicitly converts the categorical DEM to supported independent mechanisms, coalesces equal complete responses, and formats the result for Stim:
if compare_with_stim:
import stim
circuit = stim.Circuit(text)
formatted = to_independent_stim_dem(projected)
circuit_dem = circuit.detector_error_model()
reference = stim.Circuit.generated(
"surface_code:rotated_memory_z",
distance=distance,
rounds=rounds,
before_round_data_depolarization=data_depolarization,
after_clifford_depolarization=noise.after_clifford_depolarization,
before_measure_flip_probability=noise.before_measure_flip_probability,
after_reset_flip_probability=noise.after_reset_flip_probability)
def error_terms(model):
terms = {}
for instruction in model.flattened():
if instruction.type != "error":
continue
support = tuple(sorted(str(t) for t in instruction.targets_copy()))
if support in terms:
raise ValueError(
"reference comparison requires coalesced terms")
terms[support] = instruction.args_copy()[0]
return terms
actual = error_terms(formatted)
for expected in (error_terms(circuit_dem),
error_terms(reference.detector_error_model())):
if actual.keys() != expected.keys():
raise AssertionError("circuit and composed-EDEM supports differ")
for support, probability in expected.items():
if not isclose(
actual[support], probability, rel_tol=0, abs_tol=1e-12):
raise AssertionError(
"circuit and composed-EDEM probabilities differ")
(output / "reference.stim").write_text(str(reference))
(output / "circuit.dem").write_text(str(circuit_dem))
(output / "composed.dem").write_text(str(formatted))
print(f"{len(actual)} DEM error terms match the circuit and Stim reference "
"(atol=1e-12)")
factorize_dem converts the uniform Pauli channels exactly, with transformed
mechanism probabilities. coalesce_dem_columns merges equal complete
responses. to_stim_dem formats that transformed model; it performs no
probability conversion or column merging.
In the original EDEM, a depolarizing channel’s X/Y/Z or two-qubit alternatives are mutually exclusive, with the total probability split among three or fifteen nonidentity Paulis. Reusing those alternative probabilities as independent rates would change the model. The exact conversion above preserves the model’s distribution, and the saved packed files precede that conversion. The custom-noise tutorial explains independently specified Pauli errors as a separate modeling choice.
To sample, set shots to a positive number and choose sampling_mode:
"circuit"samples the physical circuit with its depolarizing channels."dem"samples the exactly factorized independent DEM derived from the composed EDEMs. Factorization runs even when the reference comparison is disabled; unsupported exact conversions raise an error.
if shots:
import stim
if sampling_mode == "circuit":
detectors, observables = stim.Circuit(text).compile_detector_sampler(
seed=seed).sample(shots, separate_observables=True)
else:
dem = to_independent_stim_dem(projected)
detectors, observables, _ = dem.compile_sampler(seed=seed).sample(shots)
print(f"Sampled {shots} shots from the {sampling_mode} model")
print(f"Shots with detection events: {detectors.any(axis=1).sum()}/{shots}")
print("Shots with uncorrected observable flips: "
f"{observables.any(axis=1).sum()}/{shots}")
Both samplers use the same physical noise model. The reported detector-event and observable-flip counts are uncorrected; the example does not run a decoder.
Physical schedule and reference equivalence
The standard-library calls correspond to Stim’s generated reference as follows:
Stim |
Example configuration |
|---|---|
|
|
|
Preparation, |
|
Same field of |
|
Same field of |
|
Same field of |
|
|
The final two CNOT interactions put X hooks horizontally and Z hooks vertically, transverse to the same-Pauli logical strings. Regression tests check the circuit’s graphlike fault distance at distances three, five, and seven. This evidence is separate from the ideal verification receipts.
Readout resets the check ancillas after their final measurements. These terminal resets are part of this concrete schedule and match the reference circuit’s reset fault sites. Their errors have zero response because the ancillas no longer interact with data.
Physical operations, noise probabilities, and detector/observable parities match the reference. Text layout, qubit numbering, ticks, and annotation coordinates can differ. Tests allow disjoint operations to commute and normalize these representation choices before comparing per-qubit operation order and parity rows. They also compare exact categorical fault-response laws before independent conversion.