-
Notifications
You must be signed in to change notification settings - Fork 22
Sliding Window Decoding Support #398
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: main
Are you sure you want to change the base?
Conversation
|
This looks great, thanks! Acknowledging receipt and intent to review 🙂 |
|
Regarding the interaction between
(Note: I tentatively renamed Regarding the construction of windows in |
|
I suppose we could also call it a @beitom do you have any thoughts? Context: we're deciding between |
|
One weird thing I noticed in the Any idea what is going on? |
I would actually go with |
|
Thanks for the input @beitom! I've been thinking of dropping Thinking out loud, though, one possibility is to remedy the aforementioned awkward fact by just making every |
I like this idea, keeps the API more clean and succinct |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@beitom sounds good. I'll make a separate issue (#399) to clean up the API and naming of SinterDecoders. No need to bog down this PR with it.
Getting back to code review (@jviszlai): this PR looks to be in pretty good shape! One final idea for substantive changes:
Using a detector coordinate works for the memory experiments produced by qldpc, but it would be valuable to make the SlidingWindowDecoder usable for externally generated circuits as well (such as those produced by tqec). Maybe instead of passing a time_idx to SlidingWindowDecoder at initialization time, we can pass a get_detector_time: Callable[[int], int] | None (naming tentative) that maps a detector to a time coordinate. If None, this would default to using the first detector coordinate. Thoughts?
This sounds good to me. I do think there might be a slightly different interpretation of
I think this is because for the
I like this idea. In the future, it could also be nice to support spatial window decoding where windows are split along a spatial coordinate instead of a time coordinate. |
|
Ok just finished some of the small changes. I might add another example to |
|
Taking a bit to run is okay. We'll just save outputs to the notebook in the repo, and add a warning/comment for the user. Anyways sorry I've been taking so long to get to this! Lots of November deadlines for me 🙂. I promise to get back to this soon (likely next week). |
|
Just added the new plot which doesn't take too long to run actually. I also had to do some debugging and added some edge case handling in case the number of rounds in the experiment is smaller than the window size. And no worries at all! I've also been a bit swamped with deadlines. |
|
@jviszlai I did a cleanup + simplify pass through the A test script to compare implementations: from qldpc import circuits, codes, decoders
circuit = circuits.get_memory_experiment(codes.SurfaceCode(2), num_rounds=9)
decoder = decoders.SlidingWindowDecoder(4, 2)
decoder.compile_decoder_for_dem(circuit.detector_error_model())
for window in decoder.windows:
print(window)Previous implementation (inside # create windows based on a time coordinate of the detectors in the DetectorErrorModel
self.windows = []
for detectors in self.detector_subsets or [range(dem.num_detectors)]:
# collect detectors according to their time index
time_to_detectors: dict[int, list[int]] = collections.defaultdict(list)
for detector in detectors:
time_to_detectors[self.detector_to_time(detector)].append(detector)
# special case: only one window that covers all detectors
if max(time_to_detectors) < self.window_size:
window_dets = [det for dets in time_to_detectors.values() for det in dets]
commit_dets = window_dets
self.windows.append((window_dets, commit_dets))
continue
# add windows one at a time
for t in range(0, max(time_to_detectors) - self.window_size + 1, self.stride):
window_dets = []
commit_dets = []
for i in range(t, t + self.window_size):
window_dets.extend(time_to_detectors[i])
if i < t + self.stride:
commit_dets.extend(time_to_detectors[i])
self.windows.append((window_dets, commit_dets))
# add trailing detectors to the last window
for i in range(t + self.stride, t + self.window_size):
self.windows[-1][1].extend(time_to_detectors[i])
for i in range(t + self.window_size, max(time_to_detectors) + 1):
self.windows[-1][0].extend(time_to_detectors[i])
self.windows[-1][1].extend(time_to_detectors[i])for which the test script prints New implementation: # construct windows defined by "detection" and "commit" regions
self.windows = []
for detectors in self.detector_subsets or [range(dem.num_detectors)]:
# collect detectors according to their time index
time_to_dets: dict[int, list[int]] = collections.defaultdict(list)
for detector in detectors:
time_to_dets[self.detector_to_time(detector)].append(detector)
# add one window at a time
for start in range(0, max(time_to_dets) + 1, self.stride):
window_time_to_dets = [time_to_dets[start + dt] for dt in range(self.window_size)]
window = ( # defined by (detection, commit) regions
[det for dets in window_time_to_dets for det in dets],
[det for dets in window_time_to_dets[: self.stride] for det in dets],
)
self.windows.append(window)which prints Which is the desired behavior? And is it easy to understand why the remainder is treated the way it is? If so, I'd like to explain this in the docstring. Might be worth mentioning that there is another possible treatment of the remainder, which would look like Also, thanks for the example notebook addition! |
|
Hi All, not sure what the timeline is for wanting to be be merged but if it can wait 1-2 weeks I'm happy to also help review. |
|
@beitom we would be happy to wait for your review 🙂 |
| " max_shots: int = 10**5,\n", | ||
| " max_errors: int = 100,\n", | ||
| " distance_trials: int = 100,\n", | ||
| " **decoding_kwargs: object,\n", |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could remove object annotations for kwargs
| " tasks = []\n", | ||
| " custom_decoders = {} # each code is going to need its own decoder\n", | ||
| " for code_index, code in enumerate(codes_to_simulate):\n", | ||
| " distance = code.get_distance(bound=distance_trials)\n", |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I noticed code.get_distance has return type int | float which causes type checking errors for other functions that expect int args, is there any reason distance will be a float?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is because distance is undefined for trivial codes (with only one code word), in which case get_distance returns np.nan, which the type checker interprets as a float. I'm open to other ideas for dealing with this edge case.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
For the notebook example it's fairly harmless, Ill take a look if this comes up anywhere in the source code and we can add a value error if distance comes in as nan.
if isinstance(distance, float):
raise ValueError(f"Code {code} has unknown distance")I'll let you decide if you want to add this to the notebook or ignore the minor warnings
| subgraph_decoders = [] | ||
| for detectors, observables in zip(self.subgraph_detectors, subgraph_observables): | ||
| # identify the error mechanisms that flip these detectors | ||
| errors = dem_arrays.detector_flip_matrix.getnnz(axis=0) != 0 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@jviszlai should this be errors = dem_arrays.detector_flip_matrix[detectors].getnnz(axis=0) != 0 only keep errors that flip at least one detector in this subgraph
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
| errors = dem_arrays.detector_flip_matrix[detectors].getnnz(axis=0) != 0 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes that sounds correct, thanks!
| } | ||
| ], | ||
| "source": [ | ||
| "codes_to_simulate = [codes.SurfaceCode(dist, rotated=True) for dist in (3, 5, 7)]\n", |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
An error comes up here
89 # collect statistics for Z and X logicals (which are corrupted, respectively, by X and Z errors)
---> 90 stats_z = run_memory_experiments(
91 codes_to_simulate, Pauli.Z, error_rates, max_shots, max_errors, distance_trials, **decoding_kwargs
92 )
93 stats_x = run_memory_experiments(
94 codes_to_simulate, Pauli.X, error_rates, max_shots, max_errors, distance_trials, **decoding_kwargs
95 )
...
File "/home/tom/qLDPC/.venv/lib/python3.11/site-packages/sinter/_data/_anon_task_stats.py", line 37, in __post_init__
assert isinstance(self.errors, int)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^After some investigation I think this is a sinter/numpy 2.x issues where numpy 2.x np.count_nonzero() returns numpy.int64 and isinstance(numpy.int64, int) is false in numpy 2.x so sinter's assert isinstance(self.errors, int) returns false. @perlinm do you know if we can force numpy<2.0
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Is this by chance related to this: quantumlib/Stim#961 ?
If so, we can update Stim requirements to 1.16 when it's released
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yep, it's exactly that. Is there an expected timeline for 1.16? If not, we can ask @Strilanc for a version bump.
| decoders. Default: False (unless decoding with MWPM). | ||
| **decoder_kwargs: Arguments to pass to qldpc.decoders.get_decoder when compiling a | ||
| custom decoder from a detector error model. | ||
| """ |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Suggest adding parameter validation
| """ | |
| """ | |
| if window_size <= 0: | |
| raise ValueError("window_size must be positive") | |
| if stride <= 0: | |
| raise ValueError("stride must be positive") | |
| if stride > window_size: | |
| raise ValueError("stride must be less than or equal to window_size") |
| if not self.detector_to_time: | ||
| dem_coords = dem.get_detector_coordinates() | ||
| self.detector_to_time = lambda det: int(dem_coords[det][0]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
error prone if called with with different DEMs will use the first DEM's coordinates since self.detector_to_time is set.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can you clarify what you mean?
Or: are you referring to the possibility that a detector_to_time passed to SlidingWindowDecoder.__init__ is incompatible with a dem passed to SlidingWindowDecoder.compile_decoder_for_dem? If so, perhaps we could validate with something like
if self.detector_to_time is not None and not all(
isinstance(self.detector_to_time(det), int) for det in range(dem.num_detectors)
):
raise ValueError(...)| qubit_ids: QubitIDs | None = None, | ||
| syndrome_measurement_strategy: SyndromeMeasurementStrategy = EdgeColoring(), | ||
| ) -> stim.Circuit: | ||
| f"""Construct a circuit for testing the performance of a code as a quantum memory. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
unused f
| """Construct a circuit for testing the performance of a code as a quantum memory. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think the string formatting is used later in the comment block for DEFAULT_IMMUNE_OP_TAG.
| *, | ||
| qubit_ids: QubitIDs | None = None, | ||
| syndrome_measurement_strategy: SyndromeMeasurementStrategy = EdgeColoring(), | ||
| ) -> tuple[stim.Circuit, stim.Circuit, stim.Circuit, MeasurementRecord, DetectorRecord, QubitIDs]: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Suggest using a data class for such complex return types e.g.
@dataclasses.dataclass
class MemoryExperimentParts:
initialization: stim.Circuit
qec_cycle: stim.Circuit
readout: stim.Circuit
measurement_record: MeasurementRecord
detector_record: DetectorRecord
qubit_ids: QubitIDs| construct the syndrome for segment S_2, namely s_2 = s_1 + H @ e_1. More generally, the | ||
| syndrome for segment S_(j+1) is s_(j+1) = s_j + H @ e_j = s_1 + H @ sum_(k=1)^j e_k. After | ||
| decoding all segments, the net error sum_(j=1)^n e_j is used to predict observable flips. | ||
| A SequentialWindowDecoder splits a detector error model into (possibly overlapping) "windows". |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could be a nice place to draw a diagram in the docstring as is done in other parts of qldpc, here's a sample I made with claude if you guys like it
Detectors: -------------------------------------------------------------->
Window i: [ ........... Detection Region ........... ]
[ Commit ]
|
+---> 1. Decode errors in Detection Region.
2. Commit to errors in Commit Region.
3. Update syndromes in future windows
based on committed errors.
4. Slide window forward.
|
v
Window i+1: [ ........... Detection Region ........... ]
[ Commit ]
|
v
...
| Returns: | ||
| A Stim circuit of OBSERVABLE_INCLUDE instructions. | ||
| """ | ||
| assert basis is Pauli.X or basis is Pauli.Z or basis is None |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Minor comment: this assert and most others on input validation should be a ValueErrors
| matrix_x = code.field.Zeros((0, len(code))) if basis is Pauli.Z else code.matrix | ||
| code = codes.CSSCode(matrix_x, matrix_z) | ||
|
|
||
| if code.is_subsystem_code: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Should add an issue to work on this. Would be great to try and run SHYPS codes with these new features. @jviszlai probably has the best idea of what would need to be done, something on the order of constructing gauge check detectors and finding the classical mapping from those to the stabilizer detectors.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
There is an issue for this: #348. I would certainly appreciate help with it! The problem may be trickers than it seems at a glance.
(Thinking out loud:) The subsystem code case requires careful choice of both:
- The measurement schedule, which is tricky because simply measuring all gauge operators in a random order may not provide enough information to construct a complete set of detectors (which can diagnose all detectable/correctable errors).
- The definition of detectors, which generally depends on the measurement schedule. This choice has no analog for "ordinary" (Abelian) stabilizer codes.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think the most straightforward subsystem schedule is alternating all Z gauges in one round and then all X gauges in the next round. The detectors would then be the parity of all gauges in each stabilizer between the two rounds they were measured. I agree though that other measurement schedules could perform better and would admit different detectors. Perhaps this is a similar situation to the different options for SyndromeMeasurementStrategy?
In terms of implementing, I think we would need to either generalize the way detectors are constructed in qldpc.circuits.memory and create new subsystem-specific SyndromeMeasurementStrategy or just create separate functionality for subsystem codes entirely.
I think we can choose either behavior. If choosing the smaller remainder preserves the logical error rate I'd vote for that. I think my original choice of the larger remainder was to cautiously keep > d rounds per window, but this might be unnecessary given the last round is reconstructed syndromes from data measurements. |
Addresses #299 and adds window decoding support.
SequentialSinterDecoderinto a more generalWindowSinterDecoderwhich allows only a subset of detectors in the segment/window to correspond to errors whose decoded results are committed.SequentialSinterDecoderis preserved as a subclass ofWindowSinterDecoderwhere all detectors in the window are part of the commit region.SlidingWindowSinterDecoderthat creates a sequence of overlapping windows from aDetectorErrorModelbased on time coordinates.num_roundsbut this was additive to theSHIFT_COORDSin the repeat block)I'd appreciate thoughts on the interaction between
WindowSinterDecoder,SequentialSinterDecoder, andSlidingWindowSinterDecodercurrently. Specifically, constructing the windows forSlidingWindowSinterDecoderincompile_decoder_for_demfeels awkward whenWindowSinterDecoderandSequentialSinterDecoderboth construct it in__init__, but I couldn't think of a cleaner solution.Also I think more general window decoding (e.g. parallel window decoding) should be easy to implement, but I'll save that for a separate PR.