diff --git a/Cargo.lock b/Cargo.lock index caf09cea5..5dbae7b93 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -4363,6 +4363,17 @@ dependencies = [ "pecos-simulators", ] +[[package]] +name = "pecos-frontier" +version = "0.2.0-dev.0" +dependencies = [ + "pecos-decoder-core", + "rand 0.10.2", + "rand_xoshiro 0.8.1", + "serde", + "serde_json", +] + [[package]] name = "pecos-fusion-blossom" version = "0.2.0-dev.0" @@ -4782,6 +4793,7 @@ dependencies = [ "rand 0.10.2", "rayon", "serde_json", + "sha2 0.11.0", "tempfile", ] @@ -4805,6 +4817,7 @@ dependencies = [ "num-complex 0.4.6", "pecos-core", "pecos-eeg", + "pecos-frontier", "pecos-neo", "pecos-qec", "pecos-quantum", diff --git a/Cargo.toml b/Cargo.toml index d3078c698..3642aae22 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -202,6 +202,7 @@ pecos-eeg = { version = "0.2.0-dev.0", path = "exp/pecos-eeg" } pecos-stab-tn = { version = "0.2.0-dev.0", path = "exp/pecos-stab-tn" } pecos-experimental = { version = "0.2.0-dev.0", path = "exp/pecos-experimental" } pecos-foreign = { version = "0.2.0-dev.0", path = "crates/pecos-foreign" } +pecos-frontier = { version = "0.2.0-dev.0", path = "exp/pecos-frontier" } pecos-fusion-blossom = { version = "0.2.0-dev.0", path = "crates/pecos-fusion-blossom" } pecos-gpu-sims = { version = "0.2.0-dev.0", path = "crates/pecos-gpu-sims" } pecos-hugr = { version = "0.2.0-dev.0", path = "crates/pecos-hugr" } diff --git a/crates/pecos-decoder-core/src/dem.rs b/crates/pecos-decoder-core/src/dem.rs index 99e3fe233..9884c6ddf 100644 --- a/crates/pecos-decoder-core/src/dem.rs +++ b/crates/pecos-decoder-core/src/dem.rs @@ -5,6 +5,17 @@ use crate::errors::DecoderError; +fn dimension_count(max_index: Option, kind: &str) -> Result { + max_index.map_or(Ok(0), |index| { + let count = u64::from(index) + 1; + usize::try_from(count).map_err(|_| { + DecoderError::InvalidConfiguration(format!( + "{kind} count for index {index} does not fit usize on this platform" + )) + }) + }) +} + /// Trait for decoders that can be constructed from detector error models pub trait DemDecoder: super::Decoder { /// Configuration type for DEM construction @@ -173,8 +184,20 @@ pub mod utils { } } - let detector_count = max_detector.map_or(0, |m| m + 1); - let observable_count = max_observable.map_or(0, |m| m + 1); + let detector_count = max_detector.map_or(Ok(0), |index| { + index.checked_add(1).ok_or_else(|| { + DecoderError::InvalidConfiguration(format!( + "detector count for index {index} does not fit usize on this platform" + )) + }) + })?; + let observable_count = max_observable.map_or(Ok(0), |index| { + index.checked_add(1).ok_or_else(|| { + DecoderError::InvalidConfiguration(format!( + "observable count for index {index} does not fit usize on this platform" + )) + }) + })?; Ok((detector_count, observable_count)) } @@ -398,8 +421,8 @@ impl SparseDem { Ok(Self { mechanisms, detector_coords, - num_detectors: max_detector.map_or(0, |m| m as usize + 1), - num_observables: max_observable.map_or(0, |m| m as usize + 1), + num_detectors: dimension_count(max_detector, "detector")?, + num_observables: dimension_count(max_observable, "observable")?, }) } } @@ -552,8 +575,8 @@ impl DemCheckMatrix { mechanisms.push((probability, detectors, observables)); } - let num_detectors = max_detector.map_or(0, |m| m as usize + 1); - let num_observables = max_observable.map_or(0, |m| m as usize + 1); + let num_detectors = dimension_count(max_detector, "detector")?; + let num_observables = dimension_count(max_observable, "observable")?; let num_mechanisms = mechanisms.len(); // Build matrices. @@ -839,8 +862,8 @@ impl DemMatchingGraph { fault_id += 1; } - let num_detectors = max_detector.map_or(0, |m| m as usize + 1); - let num_observables = max_observable.map_or(0, |m| m as usize + 1); + let num_detectors = dimension_count(max_detector, "detector")?; + let num_observables = dimension_count(max_observable, "observable")?; let edges = Self::merge_parallel_edges(edges); @@ -1262,6 +1285,30 @@ mod tests { assert_eq!(dets, 8, "parse_dem_metadata must count bare detector D7"); } + #[test] + fn sparse_dem_max_u32_detector_index_is_platform_checked() { + let dem = "error(0.01) D4294967295\n"; + let parsed = SparseDem::from_dem_str(dem); + + #[cfg(target_pointer_width = "64")] + { + assert_eq!(parsed.unwrap().num_detectors, 4_294_967_296); + assert_eq!(utils::parse_dem_metadata(dem).unwrap().0, 4_294_967_296); + } + + #[cfg(target_pointer_width = "32")] + { + // On 32-bit targets, the promoted u64 count reaches the fallible + // usize::try_from branch instead of wrapping or panicking. + let error = parsed.unwrap_err(); + assert!(matches!(error, DecoderError::InvalidConfiguration(_))); + assert!(error.to_string().contains("4294967295")); + + let metadata_error = utils::parse_dem_metadata(dem).unwrap_err(); + assert!(metadata_error.to_string().contains("4294967295")); + } + } + #[test] fn test_parsers_reject_malformed_detector_token() { // A `D` / `L` token in an error line is malformed. All three diff --git a/exp/pecos-frontier/Cargo.toml b/exp/pecos-frontier/Cargo.toml new file mode 100644 index 000000000..0953489f7 --- /dev/null +++ b/exp/pecos-frontier/Cargo.toml @@ -0,0 +1,28 @@ +[package] +name = "pecos-frontier" +version.workspace = true +edition.workspace = true +readme = "README.md" +authors.workspace = true +homepage.workspace = true +repository.workspace = true +license.workspace = true +keywords.workspace = true +categories.workspace = true +description = "Frontier approximate logical maximum-likelihood decoder for PECOS" +publish = false + +[dependencies] +pecos-decoder-core.workspace = true + +[lib] +name = "pecos_frontier" + +[dev-dependencies] +rand.workspace = true +rand_xoshiro.workspace = true +serde.workspace = true +serde_json.workspace = true + +[lints] +workspace = true diff --git a/exp/pecos-frontier/README.md b/exp/pecos-frontier/README.md new file mode 100644 index 000000000..398bf41fc --- /dev/null +++ b/exp/pecos-frontier/README.md @@ -0,0 +1,18 @@ +# PECOS Frontier Decoder + +Native Rust implementation of the Frontier approximate logical maximum-likelihood +decoder (Leverrier & Urbanke, arXiv:2606.20513). Not a wrap of the upstream +`frontier` package; the upstream implementation is used as a verification oracle. + +**Experimental** (`exp/`): the algorithm core is enumeration- and upstream-verified +(per-shot parity on matched models), but the crate has not yet accumulated real-user +mileage. Graduation to `crates/` and registration in the `pecos-decoders` meta-crate +are planned once it has been exercised more broadly (larger code families, Python +bindings, human users). + +Pruning ranks accumulated prefix log mass plus a `score_alpha`-weighted +suffix-compatibility estimate. Unpruned results are exact and upstream-verified. + +Deterministic ordering and tie-breaking are bitwise reproducible for a fixed +build and platform. The platform's `ln` and `exp` implementations may differ +across platforms. diff --git a/exp/pecos-frontier/examples/bridge_ab.rs b/exp/pecos-frontier/examples/bridge_ab.rs new file mode 100644 index 000000000..0a39bc660 --- /dev/null +++ b/exp/pecos-frontier/examples/bridge_ab.rs @@ -0,0 +1,113 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +//! Cross-implementation A/B harness: decode upstream-frontier sample shots +//! with `FrontierDecoder` on the identical model and column order. +//! +//! Input JSON (produced by an external extraction script from the upstream +//! `frontier` package): `{num_detectors, num_observables, mechanisms: +//! [[p, [detectors], [observables]], ...], shots: [{syndrome, truth_logical}]}` +//! where mechanism order IS the processing order and `syndrome` packs detector +//! `i` into bit `i`. +//! +//! Usage: `bridge_ab ` +//! Prints one `shot,predicted,truth,status` line per shot plus a summary line. + +use pecos_decoder_core::dem::SparseDem; +use pecos_frontier::{FrontierConfig, FrontierDecoder}; +use serde::Deserialize; +use std::collections::BTreeMap; + +#[derive(Deserialize)] +struct BridgeModel { + num_detectors: usize, + num_observables: usize, + mechanisms: Vec<(f64, Vec, Vec)>, + shots: Vec, +} + +#[derive(Deserialize)] +struct Shot { + /// Fired detector indices (supports arbitrary detector counts). + fired: Vec, + truth_logical: u128, +} + +fn main() { + let mut args = std::env::args().skip(1); + let path = args + .next() + .expect("usage: bridge_ab "); + let k: usize = args.next().expect("missing k").parse().expect("k"); + let delta: f64 = args.next().expect("missing delta").parse().expect("delta"); + let score_alpha: f64 = args + .next() + .expect("missing score_alpha") + .parse() + .expect("score_alpha"); + + let model: BridgeModel = + serde_json::from_str(&std::fs::read_to_string(&path).expect("read model json")) + .expect("parse model json"); + let dem = SparseDem { + mechanisms: model.mechanisms, + detector_coords: BTreeMap::new(), + num_detectors: model.num_detectors, + num_observables: model.num_observables, + }; + let config = FrontierConfig { + k, + delta, + score_alpha, + column_order: None, + }; + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, config).expect("build decoder"); + + let mut failures = 0_u32; + let mut no_path = 0_u32; + let started = std::time::Instant::now(); + assert!( + model.num_observables <= 128, + "bridge truth_logical is u128; wider observables need a format change" + ); + let mut syndrome = vec![0_u8; model.num_detectors]; + for (shot, entry) in model.shots.iter().enumerate() { + syndrome.fill(0); + for &fired in &entry.fired { + syndrome[fired as usize] = 1; + } + if let Ok(result) = decoder.decode(&syndrome) { + let words = result.predicted.words(); + assert!(words.iter().skip(2).all(|&w| w == 0), "label fits u128"); + let predicted = u128::from(words.first().copied().unwrap_or(0)) + | (u128::from(words.get(1).copied().unwrap_or(0)) << 64); + let status = if predicted == entry.truth_logical { + "ok" + } else { + failures += 1; + "logical_fail" + }; + println!("{shot},{predicted},{},{status}", entry.truth_logical); + } else { + failures += 1; + no_path += 1; + println!("{shot},,{},no_path", entry.truth_logical); + } + } + let elapsed = started.elapsed().as_secs_f64(); + let trials = u32::try_from(model.shots.len()).expect("shot count fits u32"); + println!( + "SUMMARY trials={trials} fail={failures} no_path={no_path} fer={} k={k} delta={delta} alpha={score_alpha} decode_s_mean={}", + f64::from(failures) / f64::from(trials), + elapsed / f64::from(trials), + ); +} diff --git a/exp/pecos-frontier/src/lib.rs b/exp/pecos-frontier/src/lib.rs new file mode 100644 index 000000000..e7ab9f8a4 --- /dev/null +++ b/exp/pecos-frontier/src/lib.rs @@ -0,0 +1,999 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +//! Frontier approximate logical maximum-likelihood decoding. +//! +//! The decoder performs ordered dynamic programming over independent binary +//! fault mechanisms. Prefixes with identical active detector boundary and +//! logical labels are merged by log-sum-exp, preserving degeneracy mass. The +//! configured frontier width and log-mass window provide deterministic pruning +//! for a fixed build and platform; underlying `ln`/`exp` implementations may +//! differ across platforms. +//! +//! Pruning ranks accumulated prefix log mass plus a `score_alpha`-weighted +//! suffix-compatibility estimate. Unpruned results are exact and +//! upstream-verified. + +use pecos_decoder_core::ObservableDecoder; +pub use pecos_decoder_core::dem::SparseDem; +pub use pecos_decoder_core::errors::DecoderError; +pub use pecos_decoder_core::obs_mask::ObsMask; +use std::cmp::Ordering; +use std::collections::BTreeMap; + +const WORD_BITS: usize = u64::BITS as usize; + +/// Frontier pruning and column-order configuration. +/// +/// The [`Default`] pruning values are provisional pending benchmarking. +/// Pruning ranks accumulated prefix log mass plus a `score_alpha`-weighted +/// suffix-compatibility estimate. Unpruned results are exact and +/// upstream-verified. +#[derive(Clone, Debug, PartialEq)] +pub struct FrontierConfig { + /// Maximum number of boundary states retained after each column. + pub k: usize, + /// Log-mass window below the best boundary state retained after each column. + pub delta: f64, + /// Weight applied to the suffix-compatibility score during pruning. + /// Defaults to `0.8`, matching upstream Frontier. + pub score_alpha: f64, + /// Optional permutation of the DEM mechanism indices. + pub column_order: Option>, +} + +impl Default for FrontierConfig { + fn default() -> Self { + // Provisional defaults pending benchmarking. + Self { + k: 64, + delta: 50.0, + score_alpha: 0.8, + column_order: None, + } + } +} + +/// Generate the deadline-optimized processing order for a sparse DEM. +/// +/// The input mechanism order is treated as time order. Mechanisms that can +/// close detectors earlier are placed first; detector-free mechanisms sort +/// last. The returned permutation maps target positions to source mechanism +/// indices and can be assigned directly to [`FrontierConfig::column_order`]. +/// +/// # Errors +/// +/// Returns [`DecoderError::InvalidConfiguration`] if a mechanism contains an +/// out-of-range or duplicate detector index. +pub fn deadline_column_order(dem: &SparseDem) -> Result, DecoderError> { + let time_order: Vec = (0..dem.mechanisms.len()).collect(); + deadline_order_for_sequence(dem, &time_order) +} + +/// Generate the backward deadline-optimized processing order for a sparse DEM. +/// +/// This first computes the forward deadline order, reverses that ordered +/// sequence, reruns deadline optimization in the reversed time coordinates, +/// and composes the result back to original mechanism indices. +/// +/// # Errors +/// +/// Returns [`DecoderError::InvalidConfiguration`] if a mechanism contains an +/// out-of-range or duplicate detector index. +pub fn backward_deadline_column_order(dem: &SparseDem) -> Result, DecoderError> { + let mut reversed_forward = deadline_column_order(dem)?; + reversed_forward.reverse(); + deadline_order_for_sequence(dem, &reversed_forward) +} + +/// Retained unnormalized joint log mass for one logical label. +#[derive(Clone, Debug, PartialEq)] +pub struct FrontierLogicalMass { + /// Logical-observable flip label. + pub logical: ObsMask, + /// Unnormalized joint mass `ln P(logical class, observed syndrome)`. + /// Subtract [`FrontierResult::log_evidence`] to obtain the label's log + /// posterior probability within the retained terminal mass. + pub log_mass: f64, +} + +/// Result of one Frontier decode. +#[derive(Clone, Debug, PartialEq)] +pub struct FrontierResult { + /// Predicted logical-observable flip mask. + pub predicted: ObsMask, + /// Log evidence: the logarithm of the total retained joint mass over all + /// terminal logical labels, approximating `ln P(observed syndrome)` when + /// pruning is enabled. + /// + /// The winning label's own log mass is [`Self::logical_masses`]'s first + /// entry. + pub log_evidence: f64, + /// Difference between the winning and runner-up unnormalized joint log + /// masses, if a runner-up exists. + pub runner_up_gap: Option, + /// Largest retained frontier size, including the initial boundary state. + pub peak_retained_states: usize, + /// Number of probabilistic columns processed (`0 < p < 1`). + pub processed_columns: usize, + /// Retained unnormalized joint terminal masses, ordered by mass descending + /// and numeric label ascending. The first entry is the winning label and + /// its retained log mass. + pub logical_masses: Vec, +} + +/// Direction selected by [`FrontierCommittee`]. +#[derive(Clone, Copy, Debug, Eq, PartialEq)] +pub enum CommitteeDirection { + /// The configured processing order. + Forward, + /// The plain reverse of the configured processing order. + Backward, +} + +/// Decode status for one committee leg. +#[derive(Clone, Copy, Debug, Eq, PartialEq)] +pub enum CommitteeStatus { + /// The leg retained at least one terminal state. + Ok, + /// The leg found no retained path for the syndrome. + NoPath, +} + +/// Summary of one committee leg. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct CommitteeMember { + /// Decode status for this leg. + pub status: CommitteeStatus, + /// Total retained log evidence, or negative infinity for no path. + pub log_evidence: f64, +} + +/// Result selected from the forward/backward committee. +#[derive(Clone, Debug, PartialEq)] +pub struct FrontierCommitteeResult { + /// Full result from the selected leg. + pub selected: FrontierResult, + /// Direction of the selected leg. + pub direction: CommitteeDirection, + /// Forward-leg status and evidence. + pub forward: CommitteeMember, + /// Backward-leg status and evidence. + pub backward: CommitteeMember, +} + +#[derive(Clone, Debug)] +struct Column { + detector_toggle: Vec, + logical_toggle: Vec, + close_mask: Vec, + active_mask: Vec, + suffix_compatibility: Vec, + log_odds: f64, + log_one_minus_probability: f64, +} + +#[derive(Clone, Debug)] +struct SuffixCompatibility { + word_index: usize, + bit_mask: u64, + log_probability_zero: f64, + log_probability_one: f64, +} + +#[derive(Clone, Debug, Eq, PartialEq)] +struct StateKey { + active_syndrome: Vec, + logical: Vec, +} + +impl Ord for StateKey { + fn cmp(&self, other: &Self) -> Ordering { + compare_words_as_unsigned(&self.active_syndrome, &other.active_syndrome) + .then_with(|| compare_words_as_unsigned(&self.logical, &other.logical)) + } +} + +impl PartialOrd for StateKey { + fn partial_cmp(&self, other: &Self) -> Option { + Some(self.cmp(other)) + } +} + +#[derive(Clone, Debug)] +struct Candidate { + key: StateKey, + log_mass: f64, +} + +#[derive(Clone, Debug)] +struct ScoredCandidate { + candidate: Candidate, + score: f64, +} + +/// Ordered, pruned dynamic-programming decoder for sparse detector error models. +#[derive(Clone, Debug)] +pub struct FrontierDecoder { + config: FrontierConfig, + columns: Vec, + num_detectors: usize, + detector_words: usize, + logical_words: usize, + touched_detectors: Vec, + forced_syndrome: Vec, + forced_logical: Vec, +} + +/// Two-leg Frontier decoder using a processing order and its plain reverse. +#[derive(Clone, Debug)] +pub struct FrontierCommittee { + forward: FrontierDecoder, + backward: FrontierDecoder, +} + +impl FrontierDecoder { + /// Construct a decoder from a sparse detector error model. + /// + /// Zero-probability mechanisms are discarded and probability-one mechanisms + /// are folded into the initial state after validating the optional ordering + /// permutation. All indices and probabilities are checked. + /// + /// # Errors + /// + /// Returns [`DecoderError::InvalidConfiguration`] for invalid pruning + /// parameters, probabilities, indices, or column order. + pub fn from_sparse_dem(dem: &SparseDem, config: FrontierConfig) -> Result { + validate_config(&config, dem.mechanisms.len())?; + + let detector_words = words_for(dem.num_detectors); + let logical_words = words_for(dem.num_observables); + let order = config + .column_order + .clone() + .unwrap_or_else(|| (0..dem.mechanisms.len()).collect()); + let mut raw_columns = Vec::with_capacity(dem.mechanisms.len()); + let mut forced_syndrome = vec![0; detector_words]; + let mut forced_logical = vec![0; logical_words]; + + for mechanism_index in order { + let (probability, detectors, observables) = &dem.mechanisms[mechanism_index]; + validate_probability(*probability, mechanism_index)?; + validate_indices(detectors, dem.num_detectors, "detector", mechanism_index)?; + validate_indices( + observables, + dem.num_observables, + "observable", + mechanism_index, + )?; + if *probability == 0.0 { + continue; + } + + let detector_toggle = indices_to_words(detectors, detector_words); + let logical_toggle = indices_to_words(observables, logical_words); + if probability.to_bits() == 1.0_f64.to_bits() { + xor_assign(&mut forced_syndrome, &detector_toggle); + xor_assign(&mut forced_logical, &logical_toggle); + continue; + } + + raw_columns.push((detector_toggle, logical_toggle, *probability)); + } + + let mut touched_detectors = vec![0; detector_words]; + let mut last_touch = vec![None; dem.num_detectors]; + for (column_index, (detectors, _, _)) in raw_columns.iter().enumerate() { + or_assign(&mut touched_detectors, detectors); + for detector in set_bits(detectors) { + last_touch[detector] = Some(column_index); + } + } + + // Seed with the forced contribution: detectors carrying a forced bit + // must stay in every active mask until their closing column, or the + // per-step projection would erase the bit before that column arrives. + let mut open_detectors = forced_syndrome.clone(); + let mut columns = Vec::with_capacity(raw_columns.len()); + let mut column_moments = Vec::with_capacity(raw_columns.len()); + for (column_index, (detector_toggle, logical_toggle, probability)) in + raw_columns.into_iter().enumerate() + { + or_assign(&mut open_detectors, &detector_toggle); + let mut close_mask = vec![0; detector_words]; + for (detector, &last) in last_touch.iter().enumerate() { + if last == Some(column_index) { + set_bit(&mut close_mask, detector); + } + } + and_not_assign(&mut open_detectors, &close_mask); + + column_moments.push(1.0 - 2.0 * probability); + columns.push(Column { + detector_toggle, + logical_toggle, + close_mask, + active_mask: open_detectors.clone(), + suffix_compatibility: Vec::new(), + log_odds: (probability / (1.0 - probability)).ln(), + log_one_minus_probability: (1.0 - probability).ln(), + }); + } + + let mut row_moments: Vec = vec![1.0; dem.num_detectors]; + for (column, moment) in columns + .iter_mut() + .rev() + .zip(column_moments.into_iter().rev()) + { + column.suffix_compatibility = set_bits(&column.active_mask) + .map(|detector| { + let eta = row_moments[detector]; + SuffixCompatibility { + word_index: detector / WORD_BITS, + bit_mask: 1 << (detector % WORD_BITS), + log_probability_zero: 1.0_f64.midpoint(eta).ln(), + log_probability_one: 1.0_f64.midpoint(-eta).ln(), + } + }) + .collect(); + for detector in set_bits(&column.detector_toggle) { + row_moments[detector] *= moment; + } + } + + Ok(Self { + config, + columns, + num_detectors: dem.num_detectors, + detector_words, + logical_words, + touched_detectors, + forced_syndrome, + forced_logical, + }) + } + + /// Parse a Stim-format detector error model and construct a decoder. + /// + /// # Errors + /// + /// Returns [`DecoderError`] if parsing or decoder validation fails. + pub fn from_dem_str(dem_str: &str, config: FrontierConfig) -> Result { + let dem = SparseDem::from_dem_str(dem_str)?; + Self::from_sparse_dem(&dem, config) + } + + /// Decode a dense detector syndrome. + /// + /// Every nonzero byte is treated as a fired detector. + /// + /// # Errors + /// + /// Returns [`DecoderError`] for a dimension mismatch or when the syndrome + /// is unexplainable with the retained frontier. + pub fn decode(&mut self, syndrome: &[u8]) -> Result { + if syndrome.len() != self.num_detectors { + return Err(DecoderError::InvalidDimensions { + expected: self.num_detectors, + actual: syndrome.len(), + }); + } + + let observed = syndrome_to_words(syndrome, self.detector_words); + if observed + .iter() + .zip(&self.forced_syndrome) + .zip(&self.touched_detectors) + .any(|((&seen, &forced), &touched)| (seen ^ forced) & !touched != 0) + { + return Err(unexplainable_error()); + } + + let mut initial_syndrome = self.forced_syndrome.clone(); + and_assign(&mut initial_syndrome, &self.touched_detectors); + let initial = StateKey { + active_syndrome: initial_syndrome, + logical: self.forced_logical.clone(), + }; + let mut frontier = BTreeMap::from([(initial, 0.0)]); + let mut peak_retained_states = frontier.len(); + + for column in &self.columns { + let mut merged = BTreeMap::new(); + for (state, &log_mass) in &frontier { + let branch_base = log_mass + column.log_one_minus_probability; + merge_branch(&mut merged, state.clone(), branch_base, column, &observed); + + let mut taken = state.clone(); + xor_assign(&mut taken.active_syndrome, &column.detector_toggle); + xor_assign(&mut taken.logical, &column.logical_toggle); + merge_branch( + &mut merged, + taken, + branch_base + column.log_odds, + column, + &observed, + ); + } + + if merged.is_empty() { + return Err(unexplainable_error()); + } + frontier = prune( + merged, + self.config.k, + self.config.delta, + self.config.score_alpha, + &column.suffix_compatibility, + &observed, + ); + if frontier.is_empty() { + return Err(unexplainable_error()); + } + peak_retained_states = peak_retained_states.max(frontier.len()); + } + + let mut terminal: Vec = frontier + .into_iter() + .map(|(key, log_mass)| Candidate { key, log_mass }) + .collect(); + sort_candidates(&mut terminal); + let winner = &terminal[0]; + let log_evidence = terminal.iter().fold(f64::NEG_INFINITY, |total, candidate| { + logaddexp(total, candidate.log_mass) + }); + let logical_masses = terminal + .iter() + .map(|candidate| FrontierLogicalMass { + logical: ObsMask::from_words(&candidate.key.logical), + log_mass: candidate.log_mass, + }) + .collect(); + + Ok(FrontierResult { + predicted: ObsMask::from_words(&winner.key.logical), + log_evidence, + runner_up_gap: terminal + .get(1) + .map(|runner_up| winner.log_mass - runner_up.log_mass), + peak_retained_states, + processed_columns: self.columns.len(), + logical_masses, + }) + } +} + +impl FrontierCommittee { + /// Construct a forward/backward committee from a sparse DEM. + /// + /// The forward leg uses `config.column_order` (or DEM order when absent). + /// The backward leg uses the plain reverse of that same sequence. + /// + /// # Errors + /// + /// Returns [`DecoderError::InvalidConfiguration`] when the configuration or + /// DEM is invalid. + pub fn from_sparse_dem(dem: &SparseDem, config: FrontierConfig) -> Result { + let FrontierConfig { + k, + delta, + score_alpha, + column_order, + } = config; + let mut forward_order = column_order.unwrap_or_else(|| (0..dem.mechanisms.len()).collect()); + let forward = FrontierDecoder::from_sparse_dem( + dem, + FrontierConfig { + k, + delta, + score_alpha, + column_order: Some(forward_order.clone()), + }, + )?; + forward_order.reverse(); + let backward_config = FrontierConfig { + k, + delta, + score_alpha, + column_order: Some(forward_order), + }; + let backward = FrontierDecoder::from_sparse_dem(dem, backward_config)?; + Ok(Self { forward, backward }) + } + + /// Parse a Stim-format detector error model and construct a committee. + /// + /// # Errors + /// + /// Returns [`DecoderError`] if parsing or committee construction fails. + pub fn from_dem_str(dem_str: &str, config: FrontierConfig) -> Result { + let dem = SparseDem::from_dem_str(dem_str)?; + Self::from_sparse_dem(&dem, config) + } + + /// Decode with both processing directions and select the stronger result. + /// + /// # Errors + /// + /// Returns the standard unexplainable-syndrome error if both legs find no + /// retained path. + pub fn decode(&mut self, syndrome: &[u8]) -> Result { + let forward_result = self.forward.decode(syndrome); + let backward_result = self.backward.decode(syndrome); + let forward = committee_member(&forward_result); + let backward = committee_member(&backward_result); + let (selected, direction) = match (forward_result, backward_result) { + (Err(_), Err(_)) => return Err(unexplainable_error()), + (Ok(selected), Err(_)) => (selected, CommitteeDirection::Forward), + (Err(_), Ok(selected)) => (selected, CommitteeDirection::Backward), + (Ok(forward_result), Ok(backward_result)) => { + if compare_committee_legs(Some(&forward_result), Some(&backward_result)) + == Ordering::Less + { + (backward_result, CommitteeDirection::Backward) + } else { + (forward_result, CommitteeDirection::Forward) + } + } + }; + + Ok(FrontierCommitteeResult { + selected, + direction, + forward, + backward, + }) + } +} + +impl ObservableDecoder for FrontierDecoder { + fn decode_obs(&mut self, syndrome: &[u8]) -> Result { + Ok(self.decode(syndrome)?.predicted) + } + + fn decode_to_observables(&mut self, syndrome: &[u8]) -> Result { + if self.logical_words > 1 { + return Err(DecoderError::InvalidConfiguration( + "decoder has more than 64 observables; use decode_obs() for the wide mask".into(), + )); + } + let decoded = self.decode(syndrome)?.predicted; + Ok(decoded.words().first().copied().unwrap_or(0)) + } +} + +impl ObservableDecoder for FrontierCommittee { + fn decode_obs(&mut self, syndrome: &[u8]) -> Result { + Ok(self.decode(syndrome)?.selected.predicted) + } +} + +type DeadlineKey = (usize, usize, usize, usize, usize); + +fn deadline_order_for_sequence( + dem: &SparseDem, + sequence: &[usize], +) -> Result, DecoderError> { + let sentinel = dem.mechanisms.len() + 1; + let mut first_touch = vec![sentinel; dem.num_detectors]; + let mut last_touch = vec![sentinel; dem.num_detectors]; + + for (position, &mechanism_index) in sequence.iter().enumerate() { + let detectors = &dem.mechanisms[mechanism_index].1; + validate_indices(detectors, dem.num_detectors, "detector", mechanism_index)?; + for &detector in detectors { + let detector = detector as usize; + first_touch[detector] = first_touch[detector].min(position); + last_touch[detector] = position; + } + } + + let mut positions: Vec = (0..sequence.len()).collect(); + positions.sort_by_key(|&position| -> DeadlineKey { + let mechanism_index = sequence[position]; + let detectors = &dem.mechanisms[mechanism_index].1; + if detectors.is_empty() { + return (sentinel, sentinel, sentinel, mechanism_index, position); + } + + let (earliest_last, latest_last, earliest_first) = detectors.iter().fold( + (sentinel, 0, sentinel), + |(min_last, max_last, min_first), &detector| { + let detector = detector as usize; + ( + min_last.min(last_touch[detector]), + max_last.max(last_touch[detector]), + min_first.min(first_touch[detector]), + ) + }, + ); + ( + earliest_last, + latest_last, + earliest_first, + mechanism_index, + position, + ) + }); + Ok(positions + .into_iter() + .map(|position| sequence[position]) + .collect()) +} + +fn committee_member(result: &Result) -> CommitteeMember { + match result { + Ok(decoded) => CommitteeMember { + status: CommitteeStatus::Ok, + log_evidence: decoded.log_evidence, + }, + Err(_) => CommitteeMember { + status: CommitteeStatus::NoPath, + log_evidence: f64::NEG_INFINITY, + }, + } +} + +fn compare_committee_legs( + forward: Option<&FrontierResult>, + backward: Option<&FrontierResult>, +) -> Ordering { + let forward_rank = committee_rank(forward, true); + let backward_rank = committee_rank(backward, false); + forward_rank + .iter() + .zip(backward_rank) + .find_map(|(forward_component, backward_component)| { + let ordering = forward_component.total_cmp(&backward_component); + (ordering != Ordering::Equal).then_some(ordering) + }) + .unwrap_or(Ordering::Equal) +} + +fn committee_rank(result: Option<&FrontierResult>, is_forward: bool) -> [f64; 6] { + let forward_bonus = if is_forward { 1.0 } else { 0.0 }; + let Some(result) = result else { + return [ + 1.0, + f64::NEG_INFINITY, + f64::NEG_INFINITY, + f64::NEG_INFINITY, + 0.0, + forward_bonus, + ]; + }; + + let terminal_gap = result.runner_up_gap.unwrap_or(f64::INFINITY); + let terminal_gap = if terminal_gap.is_nan() { + f64::NEG_INFINITY + } else { + terminal_gap + }; + let top_one_posterior = result + .logical_masses + .first() + .map_or(f64::NEG_INFINITY, |winner| { + winner.log_mass - result.log_evidence + }); + let top_one_posterior = if top_one_posterior.is_finite() { + top_one_posterior + } else { + f64::NEG_INFINITY + }; + [ + 2.0, + result.log_evidence, + terminal_gap, + top_one_posterior, + 0.0, + forward_bonus, + ] +} + +fn validate_config(config: &FrontierConfig, mechanism_count: usize) -> Result<(), DecoderError> { + if config.k == 0 { + return Err(DecoderError::InvalidConfiguration( + "FrontierConfig.k must be at least 1".into(), + )); + } + if config.delta.is_nan() || config.delta < 0.0 { + return Err(DecoderError::InvalidConfiguration(format!( + "FrontierConfig.delta must be non-negative and not NaN, got {}", + config.delta + ))); + } + if !config.score_alpha.is_finite() || config.score_alpha < 0.0 { + return Err(DecoderError::InvalidConfiguration(format!( + "FrontierConfig.score_alpha must be finite and non-negative, got {}", + config.score_alpha + ))); + } + if let Some(order) = &config.column_order { + if order.len() != mechanism_count { + return Err(DecoderError::InvalidConfiguration(format!( + "column_order must be a permutation of 0..{mechanism_count}" + ))); + } + let mut seen = vec![false; mechanism_count]; + for &index in order { + if index >= mechanism_count || seen[index] { + return Err(DecoderError::InvalidConfiguration(format!( + "column_order must be a permutation of 0..{mechanism_count}" + ))); + } + seen[index] = true; + } + } + Ok(()) +} + +fn validate_probability(probability: f64, index: usize) -> Result<(), DecoderError> { + if !(0.0..=1.0).contains(&probability) { + return Err(DecoderError::InvalidConfiguration(format!( + "mechanism {index} probability must satisfy 0 <= p <= 1, got {probability}" + ))); + } + Ok(()) +} + +fn validate_indices( + indices: &[u32], + upper_bound: usize, + kind: &str, + mechanism_index: usize, +) -> Result<(), DecoderError> { + let mut seen = std::collections::BTreeSet::new(); + for &index in indices { + if index as usize >= upper_bound { + return Err(DecoderError::InvalidConfiguration(format!( + "mechanism {mechanism_index} {kind} index {index} is out of range 0..{upper_bound}" + ))); + } + if !seen.insert(index) { + return Err(DecoderError::InvalidConfiguration(format!( + "mechanism {mechanism_index} repeats {kind} index {index}" + ))); + } + } + Ok(()) +} + +fn compare_words_as_unsigned(left: &[u64], right: &[u64]) -> Ordering { + left.iter().rev().cmp(right.iter().rev()) +} + +fn merge_branch( + merged: &mut BTreeMap, + mut state: StateKey, + log_mass: f64, + column: &Column, + observed: &[u64], +) { + if state + .active_syndrome + .iter() + .zip(observed) + .zip(&column.close_mask) + .any(|((&accumulated, &expected), &closing)| (accumulated ^ expected) & closing != 0) + { + return; + } + and_assign(&mut state.active_syndrome, &column.active_mask); + merged + .entry(state) + .and_modify(|mass| *mass = logaddexp(*mass, log_mass)) + .or_insert(log_mass); +} + +fn prune( + frontier: BTreeMap, + k: usize, + delta: f64, + score_alpha: f64, + suffix_compatibility: &[SuffixCompatibility], + observed: &[u64], +) -> BTreeMap { + if k == usize::MAX && delta.is_infinite() { + return frontier; + } + + let mut candidates: Vec = frontier + .into_iter() + .map(|(key, log_mass)| { + let score = if score_alpha == 0.0 { + log_mass + } else { + log_mass + + score_alpha + * suffix_compatibility_score( + &key.active_syndrome, + observed, + suffix_compatibility, + ) + }; + ScoredCandidate { + candidate: Candidate { key, log_mass }, + score, + } + }) + .collect(); + candidates.sort_by(|left, right| { + right + .score + .total_cmp(&left.score) + .then_with(|| left.candidate.key.cmp(&right.candidate.key)) + }); + let cutoff = candidates[0].score - delta; + candidates + .into_iter() + .take(k) + .take_while(|candidate| candidate.score >= cutoff) + .map(|scored| (scored.candidate.key, scored.candidate.log_mass)) + .collect() +} + +fn suffix_compatibility_score( + active_syndrome: &[u64], + observed: &[u64], + suffix_compatibility: &[SuffixCompatibility], +) -> f64 { + suffix_compatibility + .iter() + .map(|row| { + if (active_syndrome[row.word_index] ^ observed[row.word_index]) & row.bit_mask == 0 { + row.log_probability_zero + } else { + row.log_probability_one + } + }) + .sum() +} + +fn sort_candidates(candidates: &mut [Candidate]) { + candidates.sort_by(|left, right| { + right + .log_mass + .total_cmp(&left.log_mass) + .then_with(|| left.key.cmp(&right.key)) + }); +} + +fn logaddexp(left: f64, right: f64) -> f64 { + if left == f64::NEG_INFINITY { + return right; + } + if right == f64::NEG_INFINITY { + return left; + } + let (high, low) = if left.total_cmp(&right) == Ordering::Less { + (right, left) + } else { + (left, right) + }; + high + (low - high).exp().ln_1p() +} + +fn unexplainable_error() -> DecoderError { + DecoderError::DecodingFailed("syndrome is unexplainable at the given pruning parameters".into()) +} + +const fn words_for(bits: usize) -> usize { + bits.div_ceil(WORD_BITS) +} + +fn indices_to_words(indices: &[u32], word_count: usize) -> Vec { + let mut words = vec![0; word_count]; + for &index in indices { + set_bit(&mut words, index as usize); + } + words +} + +fn syndrome_to_words(syndrome: &[u8], word_count: usize) -> Vec { + let mut words = vec![0; word_count]; + for (index, &value) in syndrome.iter().enumerate() { + if value != 0 { + set_bit(&mut words, index); + } + } + words +} + +fn set_bit(words: &mut [u64], index: usize) { + words[index / WORD_BITS] |= 1 << (index % WORD_BITS); +} + +fn set_bits(words: &[u64]) -> impl Iterator + '_ { + words.iter().enumerate().flat_map(|(word_index, &word)| { + (0..WORD_BITS) + .filter(move |&bit| word & (1 << bit) != 0) + .map(move |bit| word_index * WORD_BITS + bit) + }) +} + +fn xor_assign(left: &mut [u64], right: &[u64]) { + for (left_word, &right_word) in left.iter_mut().zip(right) { + *left_word ^= right_word; + } +} + +fn or_assign(left: &mut [u64], right: &[u64]) { + for (left_word, &right_word) in left.iter_mut().zip(right) { + *left_word |= right_word; + } +} + +fn and_assign(left: &mut [u64], right: &[u64]) { + for (left_word, &right_word) in left.iter_mut().zip(right) { + *left_word &= right_word; + } +} + +fn and_not_assign(left: &mut [u64], right: &[u64]) { + for (left_word, &right_word) in left.iter_mut().zip(right) { + *left_word &= !right_word; + } +} + +#[cfg(test)] +mod tests { + use super::{FrontierLogicalMass, FrontierResult, committee_rank, logaddexp}; + use pecos_decoder_core::obs_mask::ObsMask; + + #[test] + fn logaddexp_handles_negative_infinity_on_either_side() { + assert_eq!( + logaddexp(f64::NEG_INFINITY, -2.5).to_bits(), + (-2.5_f64).to_bits() + ); + assert_eq!( + logaddexp(-2.5, f64::NEG_INFINITY).to_bits(), + (-2.5_f64).to_bits() + ); + } + + #[test] + fn committee_rank_maps_special_terminal_statistics() { + let no_runner_up = FrontierResult { + predicted: ObsMask::new(), + log_evidence: -1.0, + runner_up_gap: None, + peak_retained_states: 1, + processed_columns: 0, + logical_masses: vec![FrontierLogicalMass { + logical: ObsMask::new(), + log_mass: f64::NAN, + }], + }; + let rank = committee_rank(Some(&no_runner_up), true); + assert_eq!(rank[2].to_bits(), f64::INFINITY.to_bits()); + assert_eq!(rank[3].to_bits(), f64::NEG_INFINITY.to_bits()); + assert_eq!(rank[5].to_bits(), 1.0_f64.to_bits()); + + let nan_gap = FrontierResult { + runner_up_gap: Some(f64::NAN), + logical_masses: vec![FrontierLogicalMass { + logical: ObsMask::new(), + log_mass: -1.5, + }], + ..no_runner_up + }; + let rank = committee_rank(Some(&nan_gap), false); + assert_eq!(rank[2].to_bits(), f64::NEG_INFINITY.to_bits()); + assert_eq!(rank[3].to_bits(), (-0.5_f64).to_bits()); + + let no_path_rank = committee_rank(None, false); + assert_eq!(no_path_rank[0].to_bits(), 1.0_f64.to_bits()); + assert_eq!(no_path_rank[1].to_bits(), f64::NEG_INFINITY.to_bits()); + } +} diff --git a/exp/pecos-frontier/tests/fixtures/generate_upstream_fixtures.py b/exp/pecos-frontier/tests/fixtures/generate_upstream_fixtures.py new file mode 100644 index 000000000..1ffda19d2 --- /dev/null +++ b/exp/pecos-frontier/tests/fixtures/generate_upstream_fixtures.py @@ -0,0 +1,200 @@ +"""Golden-fixture generator for pecos-frontier, run against the upstream frontier package. + +Orchestrator-owned oracle: this script and its JSON output are authored and +committed by the reviewer, not by the implementation. The implementation must +never edit them. + +Usage (from a clone of the upstream repo with its venv built): + .venv/bin/python generate_upstream_fixtures.py > upstream_fixtures.json + +Each fixture: a binary-fault model given as mechanisms [p, [detectors], +[observables]] in processing order (identical to the column order PECOS uses), +decoded for every syndrome bitmask in `syndromes`, unpruned (K=10^9, +Delta=inf) and pruned (K/Delta from the fixture). Expected values come from +upstream FrontierResult: status, logical_hat, log_evidence, and +terminal_log_masses (logical label -> log mass, unnormalized). +""" + +from __future__ import annotations + +import json +import math +import random +import sys + +from frontier import FrontierModel, decode_frontier +from frontier.progressive import ( + FactorTransition, + OutcomeTransition, + build_frontier_layout, + columns_from_factor_transitions, +) + +UNPRUNED_K = 10**9 + + +def build_model( + mechanisms: list, + num_detectors: int, + num_observables: int, +) -> FrontierModel: + """Build an upstream FrontierModel from binary (p, detectors, observables) mechanisms.""" + factors = [] + for idx, (p, dets, obs) in enumerate(mechanisms): + det_mask = 0 + for d in dets: + det_mask |= 1 << d + log_mask = 0 + for o in obs: + log_mask |= 1 << o + factors.append( + FactorTransition( + factor_id=idx, + outcomes=( + OutcomeTransition(probability=1.0 - p, detector_mask=0, logical_mask=0), + OutcomeTransition(probability=p, detector_mask=det_mask, logical_mask=log_mask), + ), + instruction_offset=idx, + label=f"f{idx}", + ), + ) + columns = tuple(columns_from_factor_transitions(tuple(factors))) + return FrontierModel( + columns=columns, + layout=build_frontier_layout(list(columns), num_detectors=num_detectors), + num_detectors=num_detectors, + num_observables=num_observables, + ) + + +def decode_all(model: FrontierModel, syndromes: list, k: int, delta: float) -> list: + """Decode each syndrome and record upstream expected results as JSON-safe dicts.""" + out = [] + for syndrome in syndromes: + r = decode_frontier(model, syndrome, K=k, Delta=delta) + # Keep the JSON strictly standard: no-path decodes report -inf + # log_evidence upstream; emit null instead (allow_nan=False enforces). + out.append( + { + "syndrome": syndrome, + "status": r.status, + "logical_hat": r.logical_hat, + "log_evidence": r.log_evidence if math.isfinite(r.log_evidence) else None, + "terminal_log_masses": { + str(k_): v for k_, v in sorted(r.terminal_log_masses.items()) if math.isfinite(v) + }, + "engine": r.engine, + }, + ) + return out + + +def random_mechanisms( + rng: random.Random, + num_mechs: int, + num_detectors: int, + num_observables: int, +) -> list: + """Sample a seeded random binary-mechanism list.""" + mechs = [] + for _ in range(num_mechs): + n_d = rng.choice([1, 1, 2, 2, 3]) + dets = sorted(rng.sample(range(num_detectors), min(n_d, num_detectors))) + obs = sorted(rng.sample(range(num_observables), rng.choice([0, 0, 1, 1, 2]))) + p = rng.uniform(0.01, 0.4) + mechs.append((round(p, 6), dets, obs)) + return mechs + + +def main() -> int: + """Emit the fixture JSON to stdout.""" + fixtures = [] + + # F1: hand-built degeneracy case -- logical-ML differs from most-likely-error. + # Mechanism 0: p=0.20, flips detector 0, flips observable 0. + # Mechanisms 1,2: p=0.15 each, flip detector 0, no observable flip. + # Syndrome {det0}: MLE representative is mech 0 (0.20 > 0.15), but coset mass + # of label 0 (mech1 alone + mech2 alone + all three) exceeds label 1's mass. + fixtures.append( + { + "name": "degeneracy_ml_vs_mle", + "mechanisms": [[0.20, [0], [0]], [0.15, [0], []], [0.15, [0], []]], + "num_detectors": 1, + "num_observables": 1, + "syndromes": [0, 1], + "pruned": {"k": 2, "delta": 100.0}, + }, + ) + + # F2: repetition-code-like chain with hyperedge, two observables. + fixtures.append( + { + "name": "rep_chain_hyperedge", + "mechanisms": [ + [0.05, [0], [0]], + [0.04, [0, 1], []], + [0.03, [1, 2], [1]], + [0.06, [2, 3], []], + [0.02, [3], [0, 1]], + [0.07, [1, 2, 3], []], + ], + "num_detectors": 4, + "num_observables": 2, + "syndromes": list(range(16)), + "pruned": {"k": 4, "delta": 30.0}, + }, + ) + + # F3: wide observables -- winning label flips observable index 70. + fixtures.append( + { + "name": "wide_observable_70", + "mechanisms": [ + [0.10, [0], [70]], + [0.02, [0], [3]], + [0.05, [1], [0, 70]], + [0.03, [0, 1], []], + ], + "num_detectors": 2, + "num_observables": 71, + "syndromes": [0, 1, 2, 3], + "pruned": {"k": 8, "delta": 100.0}, + }, + ) + + # F4/F5: seeded random models. + for seed, num_mechs, num_dets, num_obs in ((11, 8, 5, 3), (23, 12, 6, 2)): + rng = random.Random(seed) + fixtures.append( + { + "name": f"random_seed{seed}", + "mechanisms": random_mechanisms(rng, num_mechs, num_dets, num_obs), + "num_detectors": num_dets, + "num_observables": num_obs, + "syndromes": sorted(rng.sample(range(1 << num_dets), 12)), + "pruned": {"k": 3, "delta": 20.0}, + }, + ) + + for fx in fixtures: + model = build_model(fx["mechanisms"], fx["num_detectors"], fx["num_observables"]) + fx["expected_unpruned"] = decode_all(model, fx["syndromes"], UNPRUNED_K, math.inf) + fx["expected_pruned"] = decode_all( + model, + fx["syndromes"], + fx["pruned"]["k"], + fx["pruned"]["delta"], + ) + + json.dump( + {"generator": "generate_upstream_fixtures.py", "fixtures": fixtures}, + sys.stdout, + indent=1, + allow_nan=False, + ) + sys.stdout.write("\n") + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) diff --git a/exp/pecos-frontier/tests/fixtures/generate_upstream_order_fixtures.py b/exp/pecos-frontier/tests/fixtures/generate_upstream_order_fixtures.py new file mode 100644 index 000000000..edd1eab5e --- /dev/null +++ b/exp/pecos-frontier/tests/fixtures/generate_upstream_order_fixtures.py @@ -0,0 +1,173 @@ +"""Ordering + committee golden-fixture generator, run against the upstream frontier package. + +Orchestrator-owned oracle: this script and its JSON output are authored and +committed by the reviewer, not by the implementation. The implementation must +never edit them. + +Usage (from a clone of the upstream repo with its venv built): + .venv/bin/python generate_upstream_order_fixtures.py > upstream_order_fixtures.json + +Per fixture: mechanisms in TIME order. Expected values from upstream: +- forward_ordering: optimize_column_order permutation over the time-ordered + columns (deadline reorder). ordering[target] = source index. +- backward_ordering: reverse the forward-ordered columns, then + optimize_column_order again (build_backward_deadline_ordered_family + semantics), composed back to original mechanism indices. +- committee: decode_frontier_committee on the forward-ordered model per + syndrome: status, logical_hat, direction, log_evidence. +""" + +from __future__ import annotations + +import json +import math +import random +import sys + +from frontier import FrontierModel, decode_frontier_committee +from frontier.progressive import ( + FactorTransition, + OutcomeTransition, + build_frontier_layout, + columns_from_factor_transitions, + optimize_column_order, +) +from tools.frontier_progressive import _reverse_progressive_columns + + +def build_columns(mechanisms: list) -> list: + """Build upstream progressive columns from binary mechanisms in time order.""" + factors = [] + for idx, (p, dets, obs) in enumerate(mechanisms): + det_mask = 0 + for d in dets: + det_mask |= 1 << d + log_mask = 0 + for o in obs: + log_mask |= 1 << o + factors.append( + FactorTransition( + factor_id=idx, + outcomes=( + OutcomeTransition(probability=1.0 - p, detector_mask=0, logical_mask=0), + OutcomeTransition(probability=p, detector_mask=det_mask, logical_mask=log_mask), + ), + instruction_offset=idx, + label=f"f{idx}", + ), + ) + return list(columns_from_factor_transitions(tuple(factors))) + + +def random_mechanisms( + rng: random.Random, + num_mechs: int, + num_detectors: int, + num_observables: int, +) -> list: + """Sample a seeded random binary-mechanism list (detector-free allowed).""" + mechs = [] + for _ in range(num_mechs): + n_d = rng.choice([0, 1, 1, 2, 2, 3]) + dets = sorted(rng.sample(range(num_detectors), min(n_d, num_detectors))) + obs = sorted(rng.sample(range(num_observables), rng.choice([0, 0, 1, 1, 2]))) + p = rng.uniform(0.01, 0.4) + mechs.append((round(p, 6), dets, obs)) + return mechs + + +def main() -> int: + """Emit the ordering/committee fixture JSON to stdout.""" + fixtures = [] + + # F1: hand-built chain where time order != deadline order. + fixtures.append( + { + "name": "chain_reorder", + "mechanisms": [ + [0.05, [3], [0]], + [0.04, [0, 1], []], + [0.03, [0], [1]], + [0.06, [2, 3], []], + [0.02, [1, 2], [0, 1]], + [0.07, [], [0]], + [0.08, [3], []], + ], + "num_detectors": 4, + "num_observables": 2, + "syndromes": list(range(16)), + "pruned": {"k": 3, "delta": 25.0}, + }, + ) + + # F2/F3: seeded random models (include detector-free mechanisms). + for seed, num_mechs, num_dets, num_obs in ((7, 10, 5, 2), (41, 14, 6, 3)): + rng = random.Random(seed) + fixtures.append( + { + "name": f"order_random_seed{seed}", + "mechanisms": random_mechanisms(rng, num_mechs, num_dets, num_obs), + "num_detectors": num_dets, + "num_observables": num_obs, + "syndromes": sorted(rng.sample(range(1 << num_dets), 10)), + "pruned": {"k": 3, "delta": 25.0}, + }, + ) + + for fx in fixtures: + time_columns = build_columns(fx["mechanisms"]) + forward_columns, forward_ordering = optimize_column_order( + list(time_columns), + num_detectors=fx["num_detectors"], + ) + fx["forward_ordering"] = [int(v) for v in forward_ordering] + + reversed_columns = _reverse_progressive_columns(forward_columns) + _backward_columns, backward_ordering_local = optimize_column_order( + list(reversed_columns), + num_detectors=fx["num_detectors"], + ) + # Compose back to original mechanism indices: reversed[i] came from + # forward position len-1-i, which came from mechanism forward_ordering[...]. + n = len(forward_ordering) + backward_in_original = [int(forward_ordering[n - 1 - int(local)]) for local in backward_ordering_local] + fx["backward_ordering"] = backward_in_original + + model = FrontierModel( + columns=tuple(forward_columns), + layout=build_frontier_layout(list(forward_columns), num_detectors=fx["num_detectors"]), + num_detectors=fx["num_detectors"], + num_observables=fx["num_observables"], + ) + committee = [] + for syndrome in fx["syndromes"]: + r = decode_frontier_committee( + model, + syndrome, + K=fx["pruned"]["k"], + Delta=fx["pruned"]["delta"], + ) + committee.append( + { + "syndrome": syndrome, + "status": r.status, + "logical_hat": r.logical_hat, + "direction": r.direction, + "log_evidence": r.log_evidence if math.isfinite(r.log_evidence) else None, + "engine": r.engine, + }, + ) + fx["expected_committee"] = committee + + json.dump( + {"generator": "generate_upstream_order_fixtures.py", "fixtures": fixtures}, + sys.stdout, + indent=1, + allow_nan=False, + ) + sys.stdout.write("\n") + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) diff --git a/exp/pecos-frontier/tests/fixtures/upstream_fixtures.json b/exp/pecos-frontier/tests/fixtures/upstream_fixtures.json new file mode 100644 index 000000000..82027a65b --- /dev/null +++ b/exp/pecos-frontier/tests/fixtures/upstream_fixtures.json @@ -0,0 +1,1488 @@ +{ + "generator": "generate_upstream_fixtures.py", + "fixtures": [ + { + "name": "degeneracy_ml_vs_mle", + "mechanisms": [ + [ + 0.2, + [ + 0 + ], + [ + 0 + ] + ], + [ + 0.15, + [ + 0 + ], + [] + ], + [ + 0.15, + [ + 0 + ], + [] + ] + ], + "num_detectors": 1, + "num_observables": 1, + "syndromes": [ + 0, + 1 + ], + "pruned": { + "k": 2, + "delta": 100.0 + }, + "expected_unpruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -0.4354089844812365, + "terminal_log_masses": { + "0": -0.5175146119167873, + "1": -2.9759296462578115 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 0, + "log_evidence": -1.0412872220488403, + "terminal_log_masses": { + "0": -1.5896352851379207, + "1": -1.9038089730366778 + }, + "engine": "native_binary" + } + ], + "expected_pruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -0.5050092384445508, + "terminal_log_masses": { + "0": -0.5481814103097595, + "1": -3.6690768268177565 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 1, + "log_evidence": -1.4003932854993923, + "terminal_log_masses": { + "0": -2.282782465697866, + "1": -1.93447577142965 + }, + "engine": "native_binary" + } + ] + }, + { + "name": "rep_chain_hyperedge", + "mechanisms": [ + [ + 0.05, + [ + 0 + ], + [ + 0 + ] + ], + [ + 0.04, + [ + 0, + 1 + ], + [] + ], + [ + 0.03, + [ + 1, + 2 + ], + [ + 1 + ] + ], + [ + 0.06, + [ + 2, + 3 + ], + [] + ], + [ + 0.02, + [ + 3 + ], + [ + 0, + 1 + ] + ], + [ + 0.07, + [ + 1, + 2, + 3 + ], + [] + ] + ], + "num_detectors": 4, + "num_observables": 2, + "syndromes": [ + 0, + 1, + 2, + 3, + 4, + 5, + 6, + 7, + 8, + 9, + 10, + 11, + 12, + 13, + 14, + 15 + ], + "pruned": { + "k": 4, + "delta": 30.0 + }, + "expected_unpruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -0.27716516945765746, + "terminal_log_masses": { + "0": -0.27722321191185384, + "1": -10.031530036817655 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 1, + "log_evidence": -3.217786917942595, + "terminal_log_masses": { + "0": -8.78108842722251, + "1": -3.221630385188459 + }, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 0, + "log_evidence": -5.233687599573494, + "terminal_log_masses": { + "0": -5.615426272349639, + "1": -6.381511448991767 + }, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.4491289645180316, + "terminal_log_masses": { + "0": -3.455226241668234, + "1": -8.552089030335141 + }, + "engine": "native_binary" + }, + { + "syndrome": 4, + "status": "ok", + "logical_hat": 3, + "log_evidence": -6.764670168798568, + "terminal_log_masses": { + "2": -9.069129045682928, + "3": -6.86982270281349 + }, + "engine": "native_binary" + }, + { + "syndrome": 5, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.82797229734218, + "terminal_log_masses": { + "2": -6.879539903005537, + "3": -9.818506902744613 + }, + "engine": "native_binary" + }, + { + "syndrome": 6, + "status": "ok", + "logical_hat": 2, + "log_evidence": -3.7047508205521926, + "terminal_log_masses": { + "2": -3.7532296281899296, + "3": -6.7555208339606345 + }, + "engine": "native_binary" + }, + { + "syndrome": 7, + "status": "ok", + "logical_hat": 3, + "log_evidence": -6.6145097115094025, + "terminal_log_masses": { + "2": -9.626383227722268, + "3": -6.664960641265254 + }, + "engine": "native_binary" + }, + { + "syndrome": 8, + "status": "ok", + "logical_hat": 3, + "log_evidence": -4.0608263079955025, + "terminal_log_masses": { + "2": -6.339918972287872, + "3": -4.168831489862692 + }, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.9357967928880235, + "terminal_log_masses": { + "2": -7.039693883624326, + "3": -9.251649985363196 + }, + "engine": "native_binary" + }, + { + "syndrome": 10, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.432607404086188, + "terminal_log_masses": { + "2": -6.482439701584985, + "3": -9.456512046911433 + }, + "engine": "native_binary" + }, + { + "syndrome": 11, + "status": "ok", + "logical_hat": 3, + "log_evidence": -7.130114138247697, + "terminal_log_masses": { + "2": -9.46622924710348, + "3": -7.231817558646671 + }, + "engine": "native_binary" + }, + { + "syndrome": 12, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.026106948529199, + "terminal_log_masses": { + "0": -3.028736928251697, + "1": -8.968200793089709 + }, + "engine": "native_binary" + }, + { + "syndrome": 13, + "status": "ok", + "logical_hat": 1, + "log_evidence": -5.309778798539666, + "terminal_log_masses": { + "0": -6.041915585766176, + "1": -5.965399686237199 + }, + "engine": "native_binary" + }, + { + "syndrome": 14, + "status": "ok", + "logical_hat": 0, + "log_evidence": -2.8537193014310196, + "terminal_log_masses": { + "0": -2.8639125560097964, + "1": -7.444840692719713 + }, + "engine": "native_binary" + }, + { + "syndrome": 15, + "status": "ok", + "logical_hat": 1, + "log_evidence": -5.289694649257151, + "terminal_log_masses": { + "0": -6.194399083124567, + "1": -5.808319729286401 + }, + "engine": "native_binary" + } + ], + "expected_pruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -0.27722330026295683, + "terminal_log_masses": { + "0": -0.27722330026295683 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 1, + "log_evidence": -3.221662279429397, + "terminal_log_masses": { + "1": -3.221662279429397 + }, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 1, + "log_evidence": -6.381511448991767, + "terminal_log_masses": { + "1": -6.381511448991767 + }, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.4552771306109022, + "terminal_log_masses": { + "0": -3.4552771306109022 + }, + "engine": "native_binary" + }, + { + "syndrome": 4, + "status": "ok", + "logical_hat": 3, + "log_evidence": -6.86982270281349, + "terminal_log_masses": { + "3": -6.86982270281349 + }, + "engine": "native_binary" + }, + { + "syndrome": 5, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.879539903005537, + "terminal_log_masses": { + "2": -6.879539903005537 + }, + "engine": "native_binary" + }, + { + "syndrome": 6, + "status": "ok", + "logical_hat": 2, + "log_evidence": -3.75332199009823, + "terminal_log_masses": { + "2": -3.75332199009823 + }, + "engine": "native_binary" + }, + { + "syndrome": 7, + "status": "ok", + "logical_hat": 3, + "log_evidence": -6.664960641265254, + "terminal_log_masses": { + "3": -6.664960641265254 + }, + "engine": "native_binary" + }, + { + "syndrome": 8, + "status": "ok", + "logical_hat": 3, + "log_evidence": -4.169043598373583, + "terminal_log_masses": { + "3": -4.169043598373583 + }, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 2, + "log_evidence": -7.039693883624326, + "terminal_log_masses": { + "2": -7.039693883624326 + }, + "engine": "native_binary" + }, + { + "syndrome": 10, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.504857303140179, + "terminal_log_masses": { + "2": -6.504857303140179 + }, + "engine": "native_binary" + }, + { + "syndrome": 11, + "status": "ok", + "logical_hat": 3, + "log_evidence": -7.231817558646671, + "terminal_log_masses": { + "3": -7.231817558646671 + }, + "engine": "native_binary" + }, + { + "syndrome": 12, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.028758613304906, + "terminal_log_masses": { + "0": -3.028758613304906 + }, + "engine": "native_binary" + }, + { + "syndrome": 13, + "status": "ok", + "logical_hat": 1, + "log_evidence": -5.965399686237199, + "terminal_log_masses": { + "1": -5.965399686237199 + }, + "engine": "native_binary" + }, + { + "syndrome": 14, + "status": "ok", + "logical_hat": 0, + "log_evidence": -2.855561927630394, + "terminal_log_masses": { + "0": -2.8639126443608993, + "1": -7.645142288208857 + }, + "engine": "native_binary" + }, + { + "syndrome": 15, + "status": "ok", + "logical_hat": 0, + "log_evidence": -6.194399083124567, + "terminal_log_masses": { + "0": -6.194399083124567 + }, + "engine": "native_binary" + } + ] + }, + { + "name": "wide_observable_70", + "mechanisms": [ + [ + 0.1, + [ + 0 + ], + [ + 70 + ] + ], + [ + 0.02, + [ + 0 + ], + [ + 3 + ] + ], + [ + 0.05, + [ + 1 + ], + [ + 0, + 70 + ] + ], + [ + 0.03, + [ + 0, + 1 + ], + [] + ] + ], + "num_detectors": 2, + "num_observables": 71, + "syndromes": [ + 0, + 1, + 2, + 3 + ], + "pruned": { + "k": 8, + "delta": 100.0 + }, + "expected_unpruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -0.20483714083157648, + "terminal_log_masses": { + "0": -0.20731572484760488, + "1": -8.825077971185538, + "1180591620717411303432": -6.2963606002944505, + "1180591620717411303433": -10.519673691959945 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 1180591620717411303424, + "log_evidence": -2.2235890867675474, + "terminal_log_masses": { + "8": -4.099136022958231, + "9": -12.716898269296165, + "1180591620717411303424": -2.404540302183824, + "1180591620717411303425": -6.6278533938493185 + }, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 1180591620717411303425, + "log_evidence": -3.0752084750557875, + "terminal_log_masses": { + "8": -7.575234712793504, + "9": -9.24079957946089, + "1180591620717411303424": -5.880638992019097, + "1180591620717411303425": -3.1517547040140452 + }, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.4795914490901674, + "terminal_log_masses": { + "0": -3.6834144146828782, + "1": -5.348979281350265, + "1180591620717411303432": -9.772459290129724, + "1180591620717411303433": -7.043575002124672 + }, + "engine": "native_binary" + } + ], + "expected_pruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -0.20483714083157648, + "terminal_log_masses": { + "0": -0.20731572484760488, + "1": -8.825077971185538, + "1180591620717411303432": -6.2963606002944505, + "1180591620717411303433": -10.519673691959945 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 1180591620717411303424, + "log_evidence": -2.2235890867675474, + "terminal_log_masses": { + "8": -4.099136022958231, + "9": -12.716898269296165, + "1180591620717411303424": -2.404540302183824, + "1180591620717411303425": -6.6278533938493185 + }, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 1180591620717411303425, + "log_evidence": -3.0752084750557875, + "terminal_log_masses": { + "8": -7.575234712793504, + "9": -9.24079957946089, + "1180591620717411303424": -5.880638992019097, + "1180591620717411303425": -3.1517547040140452 + }, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.4795914490901674, + "terminal_log_masses": { + "0": -3.6834144146828782, + "1": -5.348979281350265, + "1180591620717411303432": -9.772459290129724, + "1180591620717411303433": -7.043575002124672 + }, + "engine": "native_binary" + } + ] + }, + { + "name": "random_seed11", + "mechanisms": [ + [ + 0.343468, + [ + 3, + 4 + ], + [ + 2 + ] + ], + [ + 0.046708, + [ + 1 + ], + [ + 0, + 1 + ] + ], + [ + 0.24222, + [ + 0, + 1 + ], + [ + 0, + 2 + ] + ], + [ + 0.033225, + [ + 1, + 3 + ], + [ + 0, + 2 + ] + ], + [ + 0.137254, + [ + 1 + ], + [ + 0, + 1 + ] + ], + [ + 0.011792, + [ + 1, + 2, + 4 + ], + [ + 1 + ] + ], + [ + 0.399086, + [ + 3 + ], + [ + 1 + ] + ], + [ + 0.210025, + [ + 2 + ], + [ + 0 + ] + ] + ], + "num_detectors": 5, + "num_observables": 3, + "syndromes": [ + 0, + 1, + 2, + 3, + 9, + 12, + 18, + 21, + 24, + 26, + 27, + 30 + ], + "pruned": { + "k": 3, + "delta": 20.0 + }, + "expected_unpruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -1.671574535204913, + "terminal_log_masses": { + "0": -1.6765593657035884, + "2": -11.448346632924705, + "4": -7.034064799488613, + "6": -10.064525695922914 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 6, + "log_evidence": -4.284979715475963, + "terminal_log_masses": { + "0": -9.627505046429114, + "2": -6.597044149994812, + "4": -14.166460686329314, + "6": -4.394673419108197 + }, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 3, + "log_evidence": -3.144433013520559, + "terminal_log_masses": { + "1": -13.02591398437391, + "3": -3.2541267171527934, + "5": -8.48695834447371, + "7": -5.456497448039408 + }, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 5, + "log_evidence": -2.812121237160317, + "terminal_log_masses": { + "1": -8.174611501444016, + "3": -11.205072397878318, + "5": -2.8171060676589925, + "7": -12.588893334880108 + }, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 4, + "log_evidence": -4.570633149175533, + "terminal_log_masses": { + "0": -6.187769253052399, + "2": -9.218230149486702, + "4": -4.80394831605061, + "6": -14.575735583271726 + }, + "engine": "native_binary" + }, + { + "syndrome": 12, + "status": "ok", + "logical_hat": 3, + "log_evidence": -3.391991752771167, + "terminal_log_masses": { + "1": -10.532846801727313, + "3": -3.4106089907858057, + "5": -8.330476070840698, + "7": -7.949564630686005 + }, + "engine": "native_binary" + }, + { + "syndrome": 18, + "status": "ok", + "logical_hat": 5, + "log_evidence": -4.052729957929209, + "terminal_log_masses": { + "1": -5.695100038325779, + "3": -7.429805960302515, + "5": -4.311279101323989, + "7": -12.787311394087538 + }, + "engine": "native_binary" + }, + { + "syndrome": 21, + "status": "ok", + "logical_hat": 1, + "log_evidence": -6.145749443548663, + "terminal_log_masses": { + "1": -6.776600531419197, + "3": -12.603083367903139, + "5": -8.160421468420987, + "7": -7.245577934118114 + }, + "engine": "native_binary" + }, + { + "syndrome": 24, + "status": "ok", + "logical_hat": 4, + "log_evidence": -2.318700728115651, + "terminal_log_masses": { + "0": -7.681942286717397, + "2": -9.41664820869413, + "4": -2.3244368529323713, + "6": -10.800469145695923 + }, + "engine": "native_binary" + }, + { + "syndrome": 26, + "status": "ok", + "logical_hat": 7, + "log_evidence": -3.779561516643309, + "terminal_log_masses": { + "1": -7.839080857244927, + "3": -6.104374935268192, + "5": -12.378036497145127, + "7": -3.9020042043815764 + }, + "engine": "native_binary" + }, + { + "syndrome": 27, + "status": "ok", + "logical_hat": 1, + "log_evidence": -3.459247430071055, + "terminal_log_masses": { + "1": -3.4649835548877754, + "3": -11.941015847651325, + "5": -8.8224889886728, + "7": -10.557194910649535 + }, + "engine": "native_binary" + }, + { + "syndrome": 30, + "status": "ok", + "logical_hat": 6, + "log_evidence": -4.89787704322171, + "terminal_log_masses": { + "0": -6.5143061291051225, + "2": -7.429149663407996, + "4": -11.053261769005323, + "6": -5.2267789325213805 + }, + "engine": "native_binary" + } + ], + "expected_pruned": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "log_evidence": -1.6765593657035884, + "terminal_log_masses": { + "0": -1.6765593657035884 + }, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 6, + "log_evidence": -4.290596505646497, + "terminal_log_masses": { + "2": -6.604808770272863, + "6": -4.394673419108197 + }, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 3, + "log_evidence": -3.150049803691093, + "terminal_log_masses": { + "3": -3.2541267171527934, + "7": -5.464262068317459 + }, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 5, + "log_evidence": -2.8171060676589925, + "terminal_log_masses": { + "5": -2.8171060676589925 + }, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 4, + "log_evidence": -4.581860766597815, + "terminal_log_masses": { + "0": -6.19553387333045, + "4": -4.80394831605061 + }, + "engine": "native_binary" + }, + { + "syndrome": 12, + "status": "ok", + "logical_hat": 3, + "log_evidence": -3.4106089907858057, + "terminal_log_masses": { + "3": -3.4106089907858057 + }, + "engine": "native_binary" + }, + { + "syndrome": 18, + "status": "ok", + "logical_hat": 5, + "log_evidence": -4.0891915518711945, + "terminal_log_masses": { + "1": -5.70286465860383, + "5": -4.311279101323989 + }, + "engine": "native_binary" + }, + { + "syndrome": 21, + "status": "ok", + "logical_hat": 1, + "log_evidence": -6.554512981966402, + "terminal_log_masses": { + "1": -6.776600531419197, + "5": -8.168186088699038 + }, + "engine": "native_binary" + }, + { + "syndrome": 24, + "status": "ok", + "logical_hat": 4, + "log_evidence": -2.3244368529323713, + "terminal_log_masses": { + "4": -2.3244368529323713 + }, + "engine": "native_binary" + }, + { + "syndrome": 26, + "status": "ok", + "logical_hat": 7, + "log_evidence": -3.7979272909198762, + "terminal_log_masses": { + "3": -6.112139555546243, + "7": -3.9020042043815764 + }, + "engine": "native_binary" + }, + { + "syndrome": 27, + "status": "ok", + "logical_hat": 1, + "log_evidence": -3.4649835548877754, + "terminal_log_masses": { + "1": -3.4649835548877754 + }, + "engine": "native_binary" + }, + { + "syndrome": 30, + "status": "ok", + "logical_hat": 6, + "log_evidence": -5.12270201905968, + "terminal_log_masses": { + "2": -7.436914283686047, + "6": -5.2267789325213805 + }, + "engine": "native_binary" + } + ] + }, + { + "name": "random_seed23", + "mechanisms": [ + [ + 0.157882, + [ + 0, + 5 + ], + [ + 0, + 1 + ] + ], + [ + 0.01549, + [ + 1, + 5 + ], + [ + 1 + ] + ], + [ + 0.392697, + [ + 0, + 3, + 4 + ], + [] + ], + [ + 0.310602, + [ + 0, + 3 + ], + [ + 0, + 1 + ] + ], + [ + 0.252187, + [ + 0 + ], + [ + 0, + 1 + ] + ], + [ + 0.233464, + [ + 2, + 4 + ], + [ + 0, + 1 + ] + ], + [ + 0.080991, + [ + 1, + 3 + ], + [ + 0, + 1 + ] + ], + [ + 0.237096, + [ + 4 + ], + [ + 0 + ] + ], + [ + 0.231524, + [ + 0 + ], + [ + 0 + ] + ], + [ + 0.154165, + [ + 2, + 4, + 5 + ], + [] + ], + [ + 0.110783, + [ + 2, + 4 + ], + [ + 0 + ] + ], + [ + 0.285786, + [ + 4 + ], + [] + ] + ], + "num_detectors": 6, + "num_observables": 2, + "syndromes": [ + 6, + 9, + 10, + 15, + 21, + 24, + 28, + 34, + 37, + 47, + 58, + 63 + ], + "pruned": { + "k": 3, + "delta": 20.0 + }, + "expected_unpruned": [ + { + "syndrome": 6, + "status": "ok", + "logical_hat": 3, + "log_evidence": -6.261360364568913, + "terminal_log_masses": { + "0": -8.022112546944442, + "1": -7.381200831217981, + "2": -8.108959728041507, + "3": -7.328074193125911 + }, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 3, + "log_evidence": -2.746146123465367, + "terminal_log_masses": { + "0": -4.161166730740183, + "1": -4.178207023817416, + "2": -5.201547393913861, + "3": -3.5844707493527204 + }, + "engine": "native_binary" + }, + { + "syndrome": 10, + "status": "ok", + "logical_hat": 3, + "log_evidence": -4.825670397095194, + "terminal_log_masses": { + "0": -7.181705090376861, + "1": -6.69623585114158, + "2": -7.093662900641675, + "3": -5.260077866992471 + }, + "engine": "native_binary" + }, + { + "syndrome": 15, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.267659716409388, + "terminal_log_masses": { + "0": -7.536940959893618, + "1": -7.801949941948513, + "2": -7.502800422375863, + "3": -7.816172145726686 + }, + "engine": "native_binary" + }, + { + "syndrome": 21, + "status": "ok", + "logical_hat": 0, + "log_evidence": -3.728576225544649, + "terminal_log_masses": { + "0": -4.705768715947511, + "1": -6.048612948850653, + "2": -4.742088593348227, + "3": -5.546134838095753 + }, + "engine": "native_binary" + }, + { + "syndrome": 24, + "status": "ok", + "logical_hat": 3, + "log_evidence": -3.195049731702389, + "terminal_log_masses": { + "0": -5.198501326106871, + "1": -4.225594619504352, + "2": -5.2657461414212445, + "3": -4.156807045414695 + }, + "engine": "native_binary" + }, + { + "syndrome": 28, + "status": "ok", + "logical_hat": 3, + "log_evidence": -4.048424295677088, + "terminal_log_masses": { + "0": -5.663367805631611, + "1": -5.241754487058921, + "2": -5.818941255520636, + "3": -5.164254047244065 + }, + "engine": "native_binary" + }, + { + "syndrome": 34, + "status": "ok", + "logical_hat": 2, + "log_evidence": -5.919368393148876, + "terminal_log_masses": { + "0": -7.509757989395482, + "1": -7.781149574063392, + "2": -6.8199852871255935, + "3": -7.369926971975649 + }, + "engine": "native_binary" + }, + { + "syndrome": 37, + "status": "ok", + "logical_hat": 0, + "log_evidence": -4.376752960041338, + "terminal_log_masses": { + "0": -5.645368208329156, + "1": -5.805270421594858, + "2": -5.838870051154233, + "3": -5.773622533026615 + }, + "engine": "native_binary" + }, + { + "syndrome": 47, + "status": "ok", + "logical_hat": 1, + "log_evidence": -6.454195312203719, + "terminal_log_masses": { + "0": -8.016944210315367, + "1": -7.677035952817707, + "2": -7.997361477046367, + "3": -7.718822085051365 + }, + "engine": "native_binary" + }, + { + "syndrome": 58, + "status": "ok", + "logical_hat": 1, + "log_evidence": -6.346280908592142, + "terminal_log_masses": { + "0": -7.840804515838283, + "1": -7.5791547012495855, + "2": -7.971308624305485, + "3": -7.593563081590824 + }, + "engine": "native_binary" + }, + { + "syndrome": 63, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.293149510399605, + "terminal_log_masses": { + "0": -7.548704846288836, + "1": -7.991352744805838, + "2": -7.53198768208597, + "3": -7.710522924736072 + }, + "engine": "native_binary" + } + ], + "expected_pruned": [ + { + "syndrome": 6, + "status": "no_path", + "logical_hat": null, + "log_evidence": null, + "terminal_log_masses": {}, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 3, + "log_evidence": -2.9484790481067105, + "terminal_log_masses": { + "0": -4.206361122627545, + "1": -4.459088846121866, + "3": -3.6517386468059914 + }, + "engine": "native_binary" + }, + { + "syndrome": 10, + "status": "no_path", + "logical_hat": null, + "log_evidence": null, + "terminal_log_masses": {}, + "engine": "native_binary" + }, + { + "syndrome": 15, + "status": "no_path", + "logical_hat": null, + "log_evidence": null, + "terminal_log_masses": {}, + "engine": "native_binary" + }, + { + "syndrome": 21, + "status": "ok", + "logical_hat": 2, + "log_evidence": -4.618721109597032, + "terminal_log_masses": { + "0": -6.136926003885628, + "2": -4.866035499364936 + }, + "engine": "native_binary" + }, + { + "syndrome": 24, + "status": "ok", + "logical_hat": 3, + "log_evidence": -3.615451514982718, + "terminal_log_masses": { + "1": -4.490147557262374, + "2": -5.767403798098758, + "3": -4.377403835207338 + }, + "engine": "native_binary" + }, + { + "syndrome": 28, + "status": "ok", + "logical_hat": 2, + "log_evidence": -5.970649988579859, + "terminal_log_masses": { + "0": -7.488854882868454, + "2": -6.217964378347763 + }, + "engine": "native_binary" + }, + { + "syndrome": 34, + "status": "ok", + "logical_hat": 2, + "log_evidence": -6.8891711710562795, + "terminal_log_masses": { + "2": -7.006381644683862, + "3": -9.09098808483612 + }, + "engine": "native_binary" + }, + { + "syndrome": 37, + "status": "ok", + "logical_hat": 0, + "log_evidence": -5.770110042340254, + "terminal_log_masses": { + "0": -6.6332976046584315, + "1": -6.672398242645703, + "2": -7.527211797628485 + }, + "engine": "native_binary" + }, + { + "syndrome": 47, + "status": "no_path", + "logical_hat": null, + "log_evidence": null, + "terminal_log_masses": {}, + "engine": "native_binary" + }, + { + "syndrome": 58, + "status": "no_path", + "logical_hat": null, + "log_evidence": null, + "terminal_log_masses": {}, + "engine": "native_binary" + }, + { + "syndrome": 63, + "status": "no_path", + "logical_hat": null, + "log_evidence": null, + "terminal_log_masses": {}, + "engine": "native_binary" + } + ] + } + ] +} diff --git a/exp/pecos-frontier/tests/fixtures/upstream_order_fixtures.json b/exp/pecos-frontier/tests/fixtures/upstream_order_fixtures.json new file mode 100644 index 000000000..393f52993 --- /dev/null +++ b/exp/pecos-frontier/tests/fixtures/upstream_order_fixtures.json @@ -0,0 +1,721 @@ +{ + "generator": "generate_upstream_order_fixtures.py", + "fixtures": [ + { + "name": "chain_reorder", + "mechanisms": [ + [ + 0.05, + [ + 3 + ], + [ + 0 + ] + ], + [ + 0.04, + [ + 0, + 1 + ], + [] + ], + [ + 0.03, + [ + 0 + ], + [ + 1 + ] + ], + [ + 0.06, + [ + 2, + 3 + ], + [] + ], + [ + 0.02, + [ + 1, + 2 + ], + [ + 0, + 1 + ] + ], + [ + 0.07, + [], + [ + 0 + ] + ], + [ + 0.08, + [ + 3 + ], + [] + ] + ], + "num_detectors": 4, + "num_observables": 2, + "syndromes": [ + 0, + 1, + 2, + 3, + 4, + 5, + 6, + 7, + 8, + 9, + 10, + 11, + 12, + 13, + 14, + 15 + ], + "pruned": { + "k": 3, + "delta": 25.0 + }, + "forward_ordering": [ + 2, + 1, + 4, + 3, + 0, + 6, + 5 + ], + "backward_ordering": [ + 0, + 6, + 3, + 4, + 1, + 2, + 5 + ], + "expected_committee": [ + { + "syndrome": 0, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -0.28346776513264615, + "engine": "native_binary" + }, + { + "syndrome": 1, + "status": "ok", + "logical_hat": 2, + "direction": "forward", + "log_evidence": -3.759322862835824, + "engine": "native_binary" + }, + { + "syndrome": 2, + "status": "ok", + "logical_hat": 2, + "direction": "forward", + "log_evidence": -6.806188349883318, + "engine": "native_binary" + }, + { + "syndrome": 3, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -3.4613874828963334, + "engine": "native_binary" + }, + { + "syndrome": 4, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -5.005668051298697, + "engine": "native_binary" + }, + { + "syndrome": 5, + "status": "ok", + "logical_hat": 3, + "direction": "forward", + "log_evidence": -7.073753194315811, + "engine": "native_binary" + }, + { + "syndrome": 6, + "status": "ok", + "logical_hat": 3, + "direction": "forward", + "log_evidence": -4.174728408076651, + "engine": "native_binary" + }, + { + "syndrome": 7, + "status": "ok", + "logical_hat": 1, + "direction": "forward", + "log_evidence": -7.1904896117099835, + "engine": "native_binary" + }, + { + "syndrome": 8, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -2.2570814664348506, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 2, + "direction": "forward", + "log_evidence": -5.720641390924243, + "engine": "native_binary" + }, + { + "syndrome": 10, + "status": "ok", + "logical_hat": 3, + "direction": "forward", + "log_evidence": -6.7980259206306055, + "engine": "native_binary" + }, + { + "syndrome": 11, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -5.42821286277637, + "engine": "native_binary" + }, + { + "syndrome": 12, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -3.0349460619312487, + "engine": "native_binary" + }, + { + "syndrome": 13, + "status": "ok", + "logical_hat": 2, + "direction": "forward", + "log_evidence": -6.452971892213298, + "engine": "native_binary" + }, + { + "syndrome": 14, + "status": "ok", + "logical_hat": 3, + "direction": "forward", + "log_evidence": -6.120320259211788, + "engine": "native_binary" + }, + { + "syndrome": 15, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -6.180612551164233, + "engine": "native_binary" + } + ] + }, + { + "name": "order_random_seed7", + "mechanisms": [ + [ + 0.03825, + [ + 1 + ], + [ + 0 + ] + ], + [ + 0.024623, + [ + 0, + 2 + ], + [ + 0, + 1 + ] + ], + [ + 0.045378, + [ + 0, + 3 + ], + [] + ], + [ + 0.255944, + [ + 0, + 4 + ], + [] + ], + [ + 0.390739, + [ + 0, + 3 + ], + [] + ], + [ + 0.173464, + [], + [ + 0, + 1 + ] + ], + [ + 0.236824, + [ + 0, + 2 + ], + [ + 0, + 1 + ] + ], + [ + 0.251414, + [ + 0, + 1, + 2 + ], + [ + 0, + 1 + ] + ], + [ + 0.238369, + [ + 3, + 4 + ], + [ + 1 + ] + ], + [ + 0.319808, + [ + 2, + 4 + ], + [] + ] + ], + "num_detectors": 5, + "num_observables": 2, + "syndromes": [ + 2, + 9, + 10, + 14, + 15, + 16, + 18, + 23, + 28, + 31 + ], + "pruned": { + "k": 3, + "delta": 25.0 + }, + "forward_ordering": [ + 0, + 2, + 4, + 7, + 1, + 3, + 6, + 8, + 9, + 5 + ], + "backward_ordering": [ + 9, + 3, + 8, + 1, + 6, + 7, + 2, + 4, + 0, + 5 + ], + "expected_committee": [ + { + "syndrome": 2, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -3.7718451171394514, + "engine": "native_binary" + }, + { + "syndrome": 9, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -2.356365085425516, + "engine": "native_binary" + }, + { + "syndrome": 10, + "status": "no_path", + "logical_hat": null, + "direction": "forward", + "log_evidence": null, + "engine": "native_binary" + }, + { + "syndrome": 14, + "status": "ok", + "logical_hat": 3, + "direction": "forward", + "log_evidence": -3.447450196392692, + "engine": "native_binary" + }, + { + "syndrome": 15, + "status": "no_path", + "logical_hat": null, + "direction": "forward", + "log_evidence": null, + "engine": "native_binary" + }, + { + "syndrome": 16, + "status": "no_path", + "logical_hat": null, + "direction": "forward", + "log_evidence": null, + "engine": "native_binary" + }, + { + "syndrome": 18, + "status": "no_path", + "logical_hat": null, + "direction": "forward", + "log_evidence": null, + "engine": "native_binary" + }, + { + "syndrome": 23, + "status": "no_path", + "logical_hat": null, + "direction": "forward", + "log_evidence": null, + "engine": "native_binary" + }, + { + "syndrome": 28, + "status": "no_path", + "logical_hat": null, + "direction": "forward", + "log_evidence": null, + "engine": "native_binary" + }, + { + "syndrome": 31, + "status": "ok", + "logical_hat": 1, + "direction": "forward", + "log_evidence": -3.7142361581954058, + "engine": "native_binary" + } + ] + }, + { + "name": "order_random_seed41", + "mechanisms": [ + [ + 0.361782, + [ + 1, + 2 + ], + [] + ], + [ + 0.342266, + [ + 2, + 4 + ], + [ + 1, + 2 + ] + ], + [ + 0.181364, + [ + 0, + 1 + ], + [] + ], + [ + 0.11043, + [ + 1 + ], + [ + 0 + ] + ], + [ + 0.315883, + [ + 0, + 4, + 5 + ], + [] + ], + [ + 0.338008, + [], + [ + 1 + ] + ], + [ + 0.195322, + [ + 0 + ], + [ + 2 + ] + ], + [ + 0.048489, + [ + 0, + 5 + ], + [] + ], + [ + 0.030877, + [ + 5 + ], + [] + ], + [ + 0.17961, + [ + 1 + ], + [ + 0, + 2 + ] + ], + [ + 0.066904, + [], + [ + 1 + ] + ], + [ + 0.099857, + [ + 0, + 3, + 5 + ], + [ + 1 + ] + ], + [ + 0.040368, + [], + [] + ], + [ + 0.265496, + [ + 2, + 3, + 5 + ], + [ + 0, + 2 + ] + ] + ], + "num_detectors": 6, + "num_observables": 3, + "syndromes": [ + 2, + 4, + 8, + 11, + 41, + 43, + 45, + 46, + 51, + 55 + ], + "pruned": { + "k": 3, + "delta": 25.0 + }, + "forward_ordering": [ + 1, + 4, + 3, + 9, + 2, + 0, + 6, + 7, + 11, + 13, + 8, + 5, + 10, + 12 + ], + "backward_ordering": [ + 11, + 13, + 3, + 9, + 2, + 0, + 7, + 8, + 6, + 4, + 1, + 5, + 10, + 12 + ], + "expected_committee": [ + { + "syndrome": 2, + "status": "ok", + "logical_hat": 5, + "direction": "forward", + "log_evidence": -3.6936105673212936, + "engine": "native_binary" + }, + { + "syndrome": 4, + "status": "ok", + "logical_hat": 5, + "direction": "forward", + "log_evidence": -4.552330763901034, + "engine": "native_binary" + }, + { + "syndrome": 8, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -8.710128457552232, + "engine": "native_binary" + }, + { + "syndrome": 11, + "status": "ok", + "logical_hat": 6, + "direction": "backward", + "log_evidence": -6.883257276707938, + "engine": "native_binary" + }, + { + "syndrome": 41, + "status": "ok", + "logical_hat": 2, + "direction": "forward", + "log_evidence": -4.658431103478807, + "engine": "native_binary" + }, + { + "syndrome": 43, + "status": "ok", + "logical_hat": 1, + "direction": "forward", + "log_evidence": -5.475933441779121, + "engine": "native_binary" + }, + { + "syndrome": 45, + "status": "ok", + "logical_hat": 1, + "direction": "forward", + "log_evidence": -4.908295344487828, + "engine": "native_binary" + }, + { + "syndrome": 46, + "status": "ok", + "logical_hat": 0, + "direction": "forward", + "log_evidence": -4.898928838303473, + "engine": "native_binary" + }, + { + "syndrome": 51, + "status": "ok", + "logical_hat": 5, + "direction": "forward", + "log_evidence": -4.7528127098070145, + "engine": "native_binary" + }, + { + "syndrome": 55, + "status": "ok", + "logical_hat": 0, + "direction": "backward", + "log_evidence": -3.815302120068774, + "engine": "native_binary" + } + ] + } + ] +} diff --git a/exp/pecos-frontier/tests/frontier.rs b/exp/pecos-frontier/tests/frontier.rs new file mode 100644 index 000000000..4b63947b3 --- /dev/null +++ b/exp/pecos-frontier/tests/frontier.rs @@ -0,0 +1,618 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +use pecos_decoder_core::dem::SparseDem; +use pecos_decoder_core::obs_mask::ObsMask; +use pecos_decoder_core::{DecoderError, ObservableDecoder}; +use pecos_frontier::{FrontierConfig, FrontierDecoder, FrontierResult}; +use rand::{RngExt, SeedableRng}; +use rand_xoshiro::Xoshiro256PlusPlus; +use std::collections::BTreeMap; + +fn sparse_dem( + mechanisms: Vec<(f64, Vec, Vec)>, + num_detectors: usize, + num_observables: usize, +) -> SparseDem { + SparseDem { + mechanisms, + detector_coords: BTreeMap::new(), + num_detectors, + num_observables, + } +} + +fn exact_config() -> FrontierConfig { + FrontierConfig { + k: usize::MAX, + delta: f64::INFINITY, + score_alpha: 0.8, + column_order: None, + } +} + +fn independent_logaddexp(left: f64, right: f64) -> f64 { + if left == f64::NEG_INFINITY { + return right; + } + if right == f64::NEG_INFINITY { + return left; + } + let high = left.max(right); + let low = left.min(right); + high + (low - high).exp().ln_1p() +} + +fn independent_enumeration(dem: &SparseDem, observed: &[u8]) -> BTreeMap, f64> { + let logical_words = dem.num_observables.div_ceil(u64::BITS as usize); + let mut masses = BTreeMap::new(); + for subset in 0..(1usize << dem.mechanisms.len()) { + let mut detectors = vec![0_u8; dem.num_detectors]; + let mut logical = vec![0_u64; logical_words]; + let mut log_mass = 0.0; + for (index, (probability, detector_set, observable_set)) in + dem.mechanisms.iter().enumerate() + { + if subset & (1 << index) == 0 { + log_mass += (1.0 - probability).ln(); + } else { + log_mass += probability.ln(); + for &detector in detector_set { + detectors[detector as usize] ^= 1; + } + for &observable in observable_set { + logical[observable as usize / 64] ^= 1 << (observable % 64); + } + } + } + if log_mass.is_finite() && detectors == observed { + masses + .entry(logical) + .and_modify(|mass| *mass = independent_logaddexp(*mass, log_mass)) + .or_insert(log_mass); + } + } + masses +} + +fn numeric_words_cmp(left: &[u64], right: &[u64]) -> std::cmp::Ordering { + left.iter().rev().cmp(right.iter().rev()) +} + +fn independent_winner(masses: &BTreeMap, f64>) -> Vec { + masses + .iter() + .min_by(|(left_label, left_mass), (right_label, right_mass)| { + right_mass + .total_cmp(left_mass) + .then_with(|| numeric_words_cmp(left_label, right_label)) + }) + .map(|(label, _)| label.clone()) + .expect("generated syndrome must have at least one explanation") +} + +fn result_mass_map(result: &FrontierResult) -> BTreeMap, f64> { + result + .logical_masses + .iter() + .map(|entry| (entry.logical.words().to_vec(), entry.log_mass)) + .collect() +} + +#[test] +fn unpruned_matches_independent_brute_force_on_seeded_random_dems() { + let mut rng = Xoshiro256PlusPlus::seed_from_u64(0x4652_4f4e_5449_4552); + let mut forced_mechanism_count = 0; + + for case_index in 0..33 { + let mechanism_count = 4 + case_index % 11; + let num_detectors = rng.random_range(1..=10); + let num_observables = rng.random_range(1..=3); + let mut mechanisms = Vec::with_capacity(mechanism_count); + let mut sampled_subset = Vec::with_capacity(mechanism_count); + let mut observed = vec![0_u8; num_detectors]; + + for mechanism_index in 0..mechanism_count { + let forced = rng.random_bool(0.08); + let probability = if forced { + forced_mechanism_count += 1; + 1.0 + } else { + rng.random_range(0.01..0.4) + }; + let max_weight = num_detectors.min(3); + let detector_count = if mechanism_index % 3 == 0 && max_weight >= 2 { + rng.random_range(2..=max_weight) + } else { + 1 + }; + let mut detectors = Vec::with_capacity(detector_count); + while detectors.len() < detector_count { + let detector = u32::try_from(rng.random_range(0..num_detectors)) + .expect("random detector index is at most 9"); + if !detectors.contains(&detector) { + detectors.push(detector); + } + } + detectors.sort_unstable(); + + let observable_count = rng.random_range(0..=num_observables.min(2)); + let mut observables = Vec::with_capacity(observable_count); + while observables.len() < observable_count { + let observable = u32::try_from(rng.random_range(0..num_observables)) + .expect("random observable index is at most 2"); + if !observables.contains(&observable) { + observables.push(observable); + } + } + observables.sort_unstable(); + + let taken = forced || rng.random_bool(0.35); + if taken { + for &detector in &detectors { + observed[detector as usize] ^= 1; + } + } + sampled_subset.push(taken); + mechanisms.push((probability, detectors, observables)); + } + + let dem = sparse_dem(mechanisms, num_detectors, num_observables); + let mut config = exact_config(); + if case_index % 2 == 1 { + config.column_order = Some((0..mechanism_count).rev().collect()); + } + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, config).unwrap(); + let result = decoder.decode(&observed).unwrap(); + let expected = independent_enumeration(&dem, &observed); + let actual = result_mass_map(&result); + + assert_eq!(actual.len(), expected.len(), "case {case_index}"); + for (label, expected_mass) in &expected { + let actual_mass = actual.get(label).expect("logical label must be retained"); + assert!( + (actual_mass - expected_mass).abs() <= 1e-9, + "case {case_index}, label {label:?}: expected {expected_mass}, got {actual_mass}" + ); + } + assert_eq!( + result.predicted.words(), + independent_winner(&expected), + "case {case_index}, sampled subset {sampled_subset:?}" + ); + } + + assert!( + forced_mechanism_count > 0, + "seeded models must exercise forced mechanisms" + ); +} + +#[test] +fn degeneracy_mass_beats_the_single_most_likely_error() { + // For syndrome D0=1, q=0.30's lone L0-flipping fault is the most likely + // single configuration: q(1-p)^2 = 0.192 versus p(1-p)(1-q) = 0.112. + // The no-logical-flip coset nevertheless wins by degeneracy: + // (1-q)2p(1-p) = 0.224 versus q((1-p)^2+p^2) = 0.204. + let dem = sparse_dem( + vec![ + (0.30, vec![0], vec![0]), + (0.20, vec![0], vec![]), + (0.20, vec![0], vec![]), + ], + 1, + 1, + ); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let result = decoder.decode(&[1]).unwrap(); + + assert!(result.predicted.is_zero()); + assert!((result.log_evidence.exp() - (0.224 + 0.204)).abs() < 1e-12); + assert!((result.logical_masses[0].log_mass.exp() - 0.224).abs() < 1e-12); + assert!((result.runner_up_gap.unwrap() - (0.224_f64 / 0.204).ln()).abs() < 1e-12); +} + +#[test] +fn supports_wide_detectors_and_observables_without_truncation() { + let dem = sparse_dem(vec![(0.2, vec![69], vec![69])], 70, 70); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let mut syndrome = vec![0; 70]; + syndrome[69] = 7; + + let mask = decoder.decode_obs(&syndrome).unwrap(); + assert!(mask.get(69)); + assert_eq!(mask.count_ones(), 1); + assert!(matches!( + decoder.decode_to_observables(&syndrome), + Err(DecoderError::InvalidConfiguration(_)) + )); +} + +#[test] +fn wide_logical_ties_use_numeric_label_order() { + let dem = sparse_dem( + vec![(0.5, vec![0], vec![0]), (0.5, vec![0], vec![64])], + 1, + 65, + ); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let result = decoder.decode(&[1]).unwrap(); + + assert_eq!(result.predicted, ObsMask::from_u64(1)); + assert_eq!(result.runner_up_gap, Some(0.0)); + assert_eq!(result.logical_masses[0].logical, ObsMask::from_u64(1)); +} + +#[test] +fn rejects_duplicate_detector_and_observable_indices() { + let duplicate_detector = sparse_dem(vec![(0.1, vec![0, 0], vec![])], 1, 0); + let detector_error = + FrontierDecoder::from_sparse_dem(&duplicate_detector, exact_config()).unwrap_err(); + let detector_message = detector_error.to_string(); + assert!(detector_message.contains("mechanism 0")); + assert!(detector_message.contains("detector index 0")); + + let duplicate_observable = sparse_dem(vec![(0.1, vec![], vec![1, 1])], 0, 2); + let observable_error = + FrontierDecoder::from_sparse_dem(&duplicate_observable, exact_config()).unwrap_err(); + let observable_message = observable_error.to_string(); + assert!(observable_message.contains("mechanism 0")); + assert!(observable_message.contains("observable index 1")); + + let parsed_error = + FrontierDecoder::from_dem_str("error(0.1) D0 D0\n", FrontierConfig::default()).unwrap_err(); + assert!(parsed_error.to_string().contains("detector index 0")); +} + +#[test] +fn forced_only_detector_requires_its_deterministic_syndrome() { + let dem = sparse_dem(vec![(1.0, vec![0], vec![])], 1, 0); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + + let matching = decoder.decode(&[1]).unwrap(); + assert!(matching.predicted.is_zero()); + assert_eq!(matching.log_evidence.to_bits(), 0.0_f64.to_bits()); + assert_eq!(matching.processed_columns, 0); + assert!(decoder.decode(&[0]).is_err()); +} + +#[test] +fn forced_syndrome_shifts_shared_probabilistic_detector() { + let dem = sparse_dem(vec![(1.0, vec![0], vec![0]), (0.2, vec![0], vec![1])], 1, 2); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + + let probabilistic_skip = decoder.decode(&[1]).unwrap(); + assert_eq!(probabilistic_skip.predicted, ObsMask::from_u64(1)); + assert!((probabilistic_skip.log_evidence.exp() - 0.8).abs() < 1e-12); + + let probabilistic_take = decoder.decode(&[0]).unwrap(); + assert_eq!(probabilistic_take.predicted, ObsMask::from_u64(3)); + assert!((probabilistic_take.log_evidence.exp() - 0.2).abs() < 1e-12); +} + +#[test] +fn forced_logical_flip_seeds_the_winning_label() { + let dem = sparse_dem(vec![(1.0, vec![], vec![1]), (0.1, vec![], vec![0])], 0, 2); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let result = decoder.decode(&[]).unwrap(); + + assert_eq!(result.predicted, ObsMask::from_u64(2)); + assert!((result.logical_masses[0].log_mass.exp() - 0.9).abs() < 1e-12); +} + +#[test] +fn truly_empty_dem_decodes_to_empty_label_with_unit_evidence() { + let dem = sparse_dem(Vec::new(), 0, 0); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let result = decoder.decode(&[]).unwrap(); + + assert!(result.predicted.is_zero()); + assert_eq!(result.log_evidence.to_bits(), 0.0_f64.to_bits()); + assert_eq!(result.logical_masses.len(), 1); +} + +#[test] +fn interleaved_syndromes_do_not_retain_decode_state() { + let dem = sparse_dem(vec![(0.2, vec![0], vec![0]), (0.3, vec![1], vec![1])], 2, 2); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + + let first_a = decoder.decode(&[1, 0]).unwrap(); + let b = decoder.decode(&[0, 1]).unwrap(); + let second_a = decoder.decode(&[1, 0]).unwrap(); + + assert_eq!(first_a, second_a); + assert_eq!(first_a.predicted, ObsMask::from_u64(1)); + assert_eq!(b.predicted, ObsMask::from_u64(2)); + assert_ne!(first_a.predicted, b.predicted); +} + +#[test] +fn default_batch_decode_matches_individual_shots() { + let dem = sparse_dem(vec![(0.2, vec![0], vec![0]), (0.3, vec![1], vec![1])], 2, 2); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let expected = vec![ + decoder.decode_to_observables(&[1, 0]).unwrap(), + decoder.decode_to_observables(&[0, 1]).unwrap(), + decoder.decode_to_observables(&[1, 1]).unwrap(), + ]; + + let batched = decoder + .decode_batch_to_observables(&[1, 0, 0, 1, 1, 1], 3, 2) + .unwrap(); + assert_eq!(batched, expected); +} + +#[test] +fn observable_only_mechanism_has_both_terminal_labels() { + let dem = sparse_dem(vec![(0.2, vec![], vec![0])], 0, 1); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let result = decoder.decode(&[]).unwrap(); + + assert!(result.predicted.is_zero()); + assert_eq!(result.logical_masses.len(), 2); + assert!((result.logical_masses[0].log_mass.exp() - 0.8).abs() < 1e-12); + assert!((result.logical_masses[1].log_mass.exp() - 0.2).abs() < 1e-12); + assert!(result.log_evidence.abs() < 1e-12); +} + +#[test] +fn fails_loud_for_untouched_and_unachievable_syndromes() { + let untouched_dem = sparse_dem(vec![(0.2, vec![0], vec![])], 2, 0); + let mut untouched = FrontierDecoder::from_sparse_dem(&untouched_dem, exact_config()).unwrap(); + let untouched_error = untouched.decode(&[0, 1]).unwrap_err(); + assert!(untouched_error.to_string().contains("unexplainable")); + + let parity_locked_dem = sparse_dem(vec![(0.2, vec![0, 1], vec![])], 2, 0); + let mut parity_locked = + FrontierDecoder::from_sparse_dem(&parity_locked_dem, exact_config()).unwrap(); + let impossible_error = parity_locked.decode(&[1, 0]).unwrap_err(); + assert!(impossible_error.to_string().contains("pruning parameters")); +} + +#[test] +fn overpruning_can_remove_the_only_eventually_feasible_prefix() { + let dem = sparse_dem( + vec![(0.4, vec![0], vec![]), (0.1, vec![0, 1], vec![0])], + 2, + 1, + ); + let tight = FrontierConfig { + k: 1, + delta: 0.01, + score_alpha: 0.0, + column_order: None, + }; + let mut overpruned = FrontierDecoder::from_sparse_dem(&dem, tight).unwrap(); + assert!(overpruned.decode(&[0, 1]).is_err()); + + let mut default_decoder = + FrontierDecoder::from_sparse_dem(&dem, FrontierConfig::default()).unwrap(); + assert_eq!( + default_decoder.decode(&[0, 1]).unwrap().predicted, + ObsMask::from_u64(1) + ); +} + +#[test] +fn width_and_delta_pruning_can_change_the_logical_answer() { + let dem = sparse_dem( + vec![ + (0.20, vec![0], vec![]), + (0.20, vec![0], vec![]), + (0.30, vec![0], vec![0]), + ], + 1, + 1, + ); + let mut exact = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + assert!(exact.decode(&[1]).unwrap().predicted.is_zero()); + + let mut greedy = FrontierDecoder::from_sparse_dem( + &dem, + FrontierConfig { + k: 1, + delta: f64::INFINITY, + score_alpha: 0.0, + column_order: None, + }, + ) + .unwrap(); + assert_eq!(greedy.decode(&[1]).unwrap().predicted, ObsMask::from_u64(1)); + + let mut delta_pruned = FrontierDecoder::from_sparse_dem( + &dem, + FrontierConfig { + k: usize::MAX, + delta: 0.1, + score_alpha: 0.0, + column_order: None, + }, + ) + .unwrap(); + assert_eq!( + delta_pruned.decode(&[1]).unwrap().predicted, + ObsMask::from_u64(1) + ); +} + +#[test] +fn suffix_compatibility_changes_the_greedy_survivor() { + // After column 0, prefix-only scoring prefers skip (mass 0.6) over take + // (mass 0.4). For observed D0=1, the future p=0.1 column gives residual + // compatibility rho=0.1 after skip and rho=0.9 after take. At alpha=0.8: + // skip: ln(0.6) + 0.8 ln(0.1) = -2.353 + // take: ln(0.4) + 0.8 ln(0.9) = -1.001 + // Taking column 0 then skipping column 1 is the correct logical class; + // prefix-only K=1 instead keeps skip and must take logical-flipping column 1. + let dem = sparse_dem(vec![(0.4, vec![0], vec![]), (0.1, vec![0], vec![0])], 1, 1); + let mut prefix_only = FrontierDecoder::from_sparse_dem( + &dem, + FrontierConfig { + k: 1, + delta: f64::INFINITY, + score_alpha: 0.0, + column_order: None, + }, + ) + .unwrap(); + let mut suffix_scored = FrontierDecoder::from_sparse_dem( + &dem, + FrontierConfig { + k: 1, + delta: f64::INFINITY, + score_alpha: 0.8, + column_order: None, + }, + ) + .unwrap(); + + assert_eq!( + prefix_only.decode(&[1]).unwrap().predicted, + ObsMask::from_u64(1) + ); + assert!(suffix_scored.decode(&[1]).unwrap().predicted.is_zero()); +} + +#[test] +fn decoding_is_bitwise_deterministic_including_a_tie() { + let dem = sparse_dem(vec![(0.5, vec![0], vec![0]), (0.5, vec![0], vec![1])], 1, 2); + let mut first_decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let first = first_decoder.decode(&[1]).unwrap(); + let second = first_decoder.decode(&[1]).unwrap(); + let mut fresh_decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let fresh = fresh_decoder.decode(&[1]).unwrap(); + + assert_eq!(first, second); + assert_eq!(first, fresh); + assert_eq!(first.log_evidence.to_bits(), second.log_evidence.to_bits()); + assert_eq!(first.log_evidence.to_bits(), fresh.log_evidence.to_bits()); + assert_eq!(first.predicted, ObsMask::from_u64(1)); + assert_eq!(first.runner_up_gap, Some(0.0)); +} + +#[test] +fn works_through_observable_decoder_trait_object() { + let dem = sparse_dem(vec![(0.2, vec![0], vec![0])], 1, 1); + let decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + let mut boxed: Box = Box::new(decoder); + + assert_eq!(boxed.decode_obs(&[1]).unwrap(), ObsMask::from_u64(1)); +} + +#[test] +fn parses_a_stim_dem_string() { + let dem_text = "\ + detector(0, 0, 0) D0\n\ + detector(1, 0, 0) D1\n\ + logical_observable L0\n\ + error(0.1) D0 D1 L0\n\ + error(0.2) D1\n"; + let mut decoder = FrontierDecoder::from_dem_str(dem_text, exact_config()).unwrap(); + + assert_eq!(decoder.decode_obs(&[1, 1]).unwrap(), ObsMask::from_u64(1)); +} + +#[test] +fn validates_probabilities_indices_order_and_pruning_configuration() { + assert_eq!( + FrontierConfig::default().score_alpha.to_bits(), + 0.8_f64.to_bits() + ); + let probability_one_dem = sparse_dem(vec![(1.0, vec![], vec![])], 0, 0); + assert!(FrontierDecoder::from_sparse_dem(&probability_one_dem, exact_config()).is_ok()); + + for probability in [1.000_000_1, 1.1, -0.1, f64::NAN, f64::INFINITY] { + let dem = sparse_dem(vec![(probability, vec![], vec![])], 0, 0); + assert!(matches!( + FrontierDecoder::from_sparse_dem(&dem, FrontierConfig::default()), + Err(DecoderError::InvalidConfiguration(_)) + )); + } + + let zero_dem = sparse_dem(vec![(0.0, vec![], vec![])], 0, 0); + let mut zero_decoder = FrontierDecoder::from_sparse_dem(&zero_dem, exact_config()).unwrap(); + let zero_result = zero_decoder.decode(&[]).unwrap(); + assert_eq!(zero_result.processed_columns, 0); + assert_eq!(zero_result.log_evidence.to_bits(), 0.0_f64.to_bits()); + + let invalid_zero_dem = sparse_dem(vec![(0.0, vec![0], vec![])], 0, 0); + assert!(FrontierDecoder::from_sparse_dem(&invalid_zero_dem, exact_config()).is_err()); + + let bad_detector = sparse_dem(vec![(0.1, vec![1], vec![])], 1, 0); + assert!(FrontierDecoder::from_sparse_dem(&bad_detector, exact_config()).is_err()); + let bad_observable = sparse_dem(vec![(0.1, vec![], vec![1])], 0, 1); + assert!(FrontierDecoder::from_sparse_dem(&bad_observable, exact_config()).is_err()); + + let two_columns = sparse_dem(vec![(0.1, vec![], vec![]); 2], 0, 0); + for bad_order in [vec![0], vec![0, 0], vec![0, 2]] { + let config = FrontierConfig { + column_order: Some(bad_order), + ..FrontierConfig::default() + }; + assert!(FrontierDecoder::from_sparse_dem(&two_columns, config).is_err()); + } + + let k_error = FrontierDecoder::from_sparse_dem( + &zero_dem, + FrontierConfig { + k: 0, + ..FrontierConfig::default() + }, + ) + .unwrap_err(); + assert!(k_error.to_string().contains('k')); + + let negative_delta_error = FrontierDecoder::from_sparse_dem( + &zero_dem, + FrontierConfig { + delta: -0.1, + ..FrontierConfig::default() + }, + ) + .unwrap_err(); + assert!(negative_delta_error.to_string().contains("delta")); + + let nan_delta_error = FrontierDecoder::from_sparse_dem( + &zero_dem, + FrontierConfig { + delta: f64::NAN, + ..FrontierConfig::default() + }, + ) + .unwrap_err(); + assert!(nan_delta_error.to_string().contains("delta")); + + for score_alpha in [-0.1, f64::NAN, f64::INFINITY] { + let alpha_error = FrontierDecoder::from_sparse_dem( + &zero_dem, + FrontierConfig { + score_alpha, + ..FrontierConfig::default() + }, + ) + .unwrap_err(); + assert!(alpha_error.to_string().contains("score_alpha")); + } +} + +#[test] +fn rejects_wrong_syndrome_length() { + let dem = sparse_dem(vec![(0.2, vec![0], vec![])], 1, 0); + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, exact_config()).unwrap(); + assert!(matches!( + decoder.decode(&[]), + Err(DecoderError::InvalidDimensions { + expected: 1, + actual: 0 + }) + )); +} diff --git a/exp/pecos-frontier/tests/order_committee.rs b/exp/pecos-frontier/tests/order_committee.rs new file mode 100644 index 000000000..78dfe2777 --- /dev/null +++ b/exp/pecos-frontier/tests/order_committee.rs @@ -0,0 +1,291 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +use pecos_decoder_core::ObservableDecoder; +use pecos_decoder_core::dem::SparseDem; +use pecos_decoder_core::obs_mask::ObsMask; +use pecos_frontier::{ + CommitteeDirection, CommitteeStatus, FrontierCommittee, FrontierConfig, + backward_deadline_column_order, deadline_column_order, +}; +use serde::Deserialize; +use std::collections::BTreeMap; + +const FIXTURES_JSON: &str = include_str!("fixtures/upstream_order_fixtures.json"); + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct FixtureFile { + generator: String, + fixtures: Vec, +} + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct Fixture { + name: String, + mechanisms: Vec<(f64, Vec, Vec)>, + num_detectors: usize, + num_observables: usize, + syndromes: Vec, + pruned: PruningConfig, + forward_ordering: Vec, + backward_ordering: Vec, + expected_committee: Vec, +} + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct PruningConfig { + k: usize, + delta: f64, +} + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct ExpectedCommitteeResult { + syndrome: u128, + status: String, + logical_hat: Option, + direction: String, + log_evidence: Option, + engine: String, +} + +fn sparse_dem( + mechanisms: Vec<(f64, Vec, Vec)>, + num_detectors: usize, + num_observables: usize, +) -> SparseDem { + SparseDem { + mechanisms, + detector_coords: BTreeMap::new(), + num_detectors, + num_observables, + } +} + +fn dense_syndrome(mask: u128, num_detectors: usize) -> Vec { + (0..num_detectors) + .map(|bit| u8::from(mask & (1_u128 << bit) != 0)) + .collect() +} + +fn mask_as_u128(mask: &ObsMask) -> u128 { + assert!( + mask.words().iter().skip(2).all(|&word| word == 0), + "fixture labels must fit in u128" + ); + u128::from(mask.words().first().copied().unwrap_or(0)) + | (u128::from(mask.words().get(1).copied().unwrap_or(0)) << 64) +} + +fn expected_direction(direction: &str) -> CommitteeDirection { + match direction { + "forward" => CommitteeDirection::Forward, + "backward" => CommitteeDirection::Backward, + unknown => panic!("unknown committee direction {unknown}"), + } +} + +#[test] +fn ordering_and_committee_match_upstream_fixtures() { + let fixture_file: FixtureFile = + serde_json::from_str(FIXTURES_JSON).expect("order fixtures must parse"); + assert_eq!( + fixture_file.generator, + "generate_upstream_order_fixtures.py" + ); + let mut backward_selections = 0; + let mut single_leg_selections = 0; + + for fixture in fixture_file.fixtures { + let dem = sparse_dem( + fixture.mechanisms, + fixture.num_detectors, + fixture.num_observables, + ); + assert_eq!( + deadline_column_order(&dem).unwrap(), + fixture.forward_ordering, + "{}: forward ordering", + fixture.name + ); + assert_eq!( + backward_deadline_column_order(&dem).unwrap(), + fixture.backward_ordering, + "{}: backward ordering", + fixture.name + ); + + // The fixture mechanisms are in time order. Passing forward_ordering + // builds the same forward-ordered model used by the upstream generator; + // FrontierCommittee then makes its second leg by plain reversal. + let mut committee = FrontierCommittee::from_sparse_dem( + &dem, + FrontierConfig { + k: fixture.pruned.k, + delta: fixture.pruned.delta, + score_alpha: 0.8, + column_order: Some(fixture.forward_ordering.clone()), + }, + ) + .unwrap(); + + let mut expected_by_syndrome = BTreeMap::new(); + for expected in fixture.expected_committee { + assert!( + expected_by_syndrome + .insert(expected.syndrome, expected) + .is_none(), + "{}: duplicate expected syndrome", + fixture.name + ); + } + assert_eq!( + expected_by_syndrome.len(), + fixture.syndromes.len(), + "{}: committee result count", + fixture.name + ); + + for syndrome_mask in fixture.syndromes { + let expected = expected_by_syndrome + .remove(&syndrome_mask) + .unwrap_or_else(|| panic!("{}: missing expected syndrome", fixture.name)); + assert_eq!(expected.engine, "native_binary"); + let syndrome = dense_syndrome(syndrome_mask, fixture.num_detectors); + let decoded = committee.decode(&syndrome); + + match expected.status.as_str() { + "ok" => { + let result = decoded.unwrap_or_else(|error| { + panic!( + "{} syndrome {syndrome_mask}: expected success, got {error}", + fixture.name + ) + }); + let direction = expected_direction(&expected.direction); + assert_eq!( + result.direction, direction, + "{} syndrome {syndrome_mask}: direction", + fixture.name + ); + assert_eq!( + mask_as_u128(&result.selected.predicted), + expected.logical_hat.expect("ok result needs logical_hat"), + "{} syndrome {syndrome_mask}: logical label", + fixture.name + ); + let expected_evidence = + expected.log_evidence.expect("ok result needs log_evidence"); + assert!( + (result.selected.log_evidence - expected_evidence).abs() <= 1e-9, + "{} syndrome {syndrome_mask}: expected evidence {expected_evidence}, got {}", + fixture.name, + result.selected.log_evidence + ); + let selected_member = match result.direction { + CommitteeDirection::Forward => result.forward, + CommitteeDirection::Backward => result.backward, + }; + assert_eq!(selected_member.status, CommitteeStatus::Ok); + assert_eq!( + selected_member.log_evidence.to_bits(), + result.selected.log_evidence.to_bits() + ); + if result.direction == CommitteeDirection::Backward { + backward_selections += 1; + } + if result.forward.status != result.backward.status { + single_leg_selections += 1; + let failed_member = if result.forward.status == CommitteeStatus::NoPath { + result.forward + } else { + result.backward + }; + assert_eq!( + failed_member.log_evidence.to_bits(), + f64::NEG_INFINITY.to_bits() + ); + } + } + "no_path" => { + assert!(expected.logical_hat.is_none()); + assert!(expected.log_evidence.is_none()); + assert!( + decoded.is_err(), + "{} syndrome {syndrome_mask}", + fixture.name + ); + } + status => panic!("{}: unknown committee status {status}", fixture.name), + } + } + assert!(expected_by_syndrome.is_empty()); + } + + assert!(backward_selections > 0, "fixtures must select backward"); + assert!( + single_leg_selections > 0, + "fixtures must exercise one-leg no-path selection" + ); +} + +#[test] +fn deadline_order_closes_earlier_rows_first() { + // D0 has first/last touches 0/2; D1 has 1/3. The keys therefore place + // both D0 columns before both D1 columns, preserving original-index ties: + // [0(D0), 2(D0), 1(D1), 3(D1)]. + let dem = sparse_dem( + vec![ + (0.1, vec![0], vec![]), + (0.1, vec![1], vec![]), + (0.1, vec![0], vec![]), + (0.1, vec![1], vec![]), + ], + 2, + 0, + ); + assert_eq!(deadline_column_order(&dem).unwrap(), vec![0, 2, 1, 3]); +} + +#[test] +fn detector_free_mechanisms_sort_last_and_empty_dem_stays_empty() { + let dem = sparse_dem(vec![(0.1, vec![], vec![0]), (0.1, vec![0], vec![])], 1, 1); + assert_eq!(deadline_column_order(&dem).unwrap(), vec![1, 0]); + + let empty = sparse_dem(Vec::new(), 0, 0); + assert!(deadline_column_order(&empty).unwrap().is_empty()); + assert!(backward_deadline_column_order(&empty).unwrap().is_empty()); +} + +#[test] +fn ordering_rejects_invalid_and_duplicate_detector_indices() { + let out_of_range = sparse_dem(vec![(0.1, vec![1], vec![])], 1, 0); + assert!(deadline_column_order(&out_of_range).is_err()); + + let duplicate = sparse_dem(vec![(0.1, vec![0, 0], vec![])], 1, 0); + assert!(backward_deadline_column_order(&duplicate).is_err()); +} + +#[test] +fn committee_ties_select_forward_and_trait_uses_selected_mask() { + let mut committee = FrontierCommittee::from_dem_str("", FrontierConfig::default()).unwrap(); + let result = committee.decode(&[]).unwrap(); + assert_eq!(result.direction, CommitteeDirection::Forward); + assert_eq!(result.forward.status, CommitteeStatus::Ok); + assert_eq!(result.backward.status, CommitteeStatus::Ok); + + let mut boxed: Box = Box::new(committee); + assert!(boxed.decode_obs(&[]).unwrap().is_zero()); +} diff --git a/exp/pecos-frontier/tests/upstream_fixtures.rs b/exp/pecos-frontier/tests/upstream_fixtures.rs new file mode 100644 index 000000000..43ef7e580 --- /dev/null +++ b/exp/pecos-frontier/tests/upstream_fixtures.rs @@ -0,0 +1,253 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +use pecos_decoder_core::dem::SparseDem; +use pecos_decoder_core::obs_mask::ObsMask; +use pecos_frontier::{FrontierConfig, FrontierDecoder, FrontierResult}; +use serde::Deserialize; +use std::collections::BTreeMap; + +const FIXTURES_JSON: &str = include_str!("fixtures/upstream_fixtures.json"); + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct FixtureFile { + generator: String, + fixtures: Vec, +} + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct Fixture { + name: String, + mechanisms: Vec<(f64, Vec, Vec)>, + num_detectors: usize, + num_observables: usize, + syndromes: Vec, + pruned: PruningConfig, + expected_unpruned: Vec, + expected_pruned: Vec, +} + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct PruningConfig { + k: usize, + delta: f64, +} + +#[derive(Debug, Deserialize)] +#[serde(deny_unknown_fields)] +struct ExpectedResult { + syndrome: u128, + status: String, + logical_hat: Option, + log_evidence: Option, + terminal_log_masses: BTreeMap, + engine: String, +} + +fn parse_fixtures() -> FixtureFile { + serde_json::from_str(FIXTURES_JSON).expect("upstream fixture file must parse") +} + +fn dense_syndrome(mask: u128, num_detectors: usize) -> Vec { + (0..num_detectors) + .map(|bit| u8::from(mask & (1_u128 << bit) != 0)) + .collect() +} + +fn mask_as_u128(mask: &ObsMask) -> u128 { + assert!( + mask.words().iter().skip(2).all(|&word| word == 0), + "fixture labels must fit in u128" + ); + u128::from(mask.words().first().copied().unwrap_or(0)) + | (u128::from(mask.words().get(1).copied().unwrap_or(0)) << 64) +} + +fn actual_masses(result: &FrontierResult) -> BTreeMap { + result + .logical_masses + .iter() + .map(|entry| (mask_as_u128(&entry.logical), entry.log_mass)) + .collect() +} + +fn assert_expected_results( + fixture: &Fixture, + config: FrontierConfig, + expected_results: &[ExpectedResult], + regime: &str, +) { + let dem = SparseDem { + mechanisms: fixture.mechanisms.clone(), + detector_coords: BTreeMap::new(), + num_detectors: fixture.num_detectors, + num_observables: fixture.num_observables, + }; + let mut decoder = FrontierDecoder::from_sparse_dem(&dem, config).unwrap_or_else(|error| { + panic!( + "{} {regime}: decoder construction failed: {error}", + fixture.name + ) + }); + + let mut expected_by_syndrome = BTreeMap::new(); + for expected in expected_results { + assert!( + expected_by_syndrome + .insert(expected.syndrome, expected) + .is_none(), + "{} {regime}: duplicate expected syndrome", + fixture.name + ); + } + assert_eq!( + expected_by_syndrome.len(), + fixture.syndromes.len(), + "{} {regime}: expected result count differs from syndrome count", + fixture.name + ); + + for &syndrome_mask in &fixture.syndromes { + let expected = expected_by_syndrome + .remove(&syndrome_mask) + .unwrap_or_else(|| panic!("{} {regime}: missing expected result", fixture.name)); + assert_eq!( + expected.engine, "native_binary", + "{} {regime} syndrome {syndrome_mask}: unexpected engine", + fixture.name + ); + let syndrome = dense_syndrome(syndrome_mask, fixture.num_detectors); + let decoded = decoder.decode(&syndrome); + + match expected.status.as_str() { + "ok" => { + let result = decoded.unwrap_or_else(|error| { + panic!( + "{} {regime} syndrome {syndrome_mask}: expected success, got {error}", + fixture.name + ) + }); + let logical_hat = expected.logical_hat.expect("ok result needs logical_hat"); + assert_eq!( + mask_as_u128(&result.predicted), + logical_hat, + "{} {regime} syndrome {syndrome_mask}: predicted label", + fixture.name + ); + + let expected_masses: BTreeMap = expected + .terminal_log_masses + .iter() + .map(|(label, &mass)| { + (label.parse().expect("logical label must fit in u128"), mass) + }) + .collect(); + let actual_masses = actual_masses(&result); + assert_eq!( + actual_masses.len(), + expected_masses.len(), + "{} {regime} syndrome {syndrome_mask}: terminal-label count", + fixture.name + ); + for (label, expected_mass) in expected_masses { + let actual_mass = actual_masses.get(&label).unwrap_or_else(|| { + panic!( + "{} {regime} syndrome {syndrome_mask}: missing label {label}", + fixture.name + ) + }); + assert!( + (actual_mass - expected_mass).abs() <= 1e-9, + "{} {regime} syndrome {syndrome_mask}, label {label}: expected {expected_mass}, got {actual_mass}", + fixture.name + ); + } + + let expected_evidence = expected + .log_evidence + .expect("ok result needs finite log_evidence"); + assert!( + (result.log_evidence - expected_evidence).abs() <= 1e-9, + "{} {regime} syndrome {syndrome_mask}: expected evidence {expected_evidence}, got {}", + fixture.name, + result.log_evidence + ); + } + "no_path" => { + assert!( + expected.logical_hat.is_none(), + "{} {regime} syndrome {syndrome_mask}: no-path logical_hat must be null", + fixture.name + ); + assert!( + expected.log_evidence.is_none(), + "{} {regime} syndrome {syndrome_mask}: no-path log_evidence must be null", + fixture.name + ); + assert!( + expected.terminal_log_masses.is_empty(), + "{} {regime} syndrome {syndrome_mask}: no-path masses must be empty", + fixture.name + ); + assert!( + decoded.is_err(), + "{} {regime} syndrome {syndrome_mask}: expected no path", + fixture.name + ); + } + status => panic!( + "{} {regime} syndrome {syndrome_mask}: unknown fixture status {status}", + fixture.name + ), + } + } + + assert!( + expected_by_syndrome.is_empty(), + "{} {regime}: expected results contain extra syndromes", + fixture.name + ); +} + +#[test] +fn unpruned_and_pruned_results_match_upstream_golden_fixtures() { + let fixture_file = parse_fixtures(); + assert_eq!(fixture_file.generator, "generate_upstream_fixtures.py"); + + for fixture in fixture_file.fixtures { + assert_expected_results( + &fixture, + FrontierConfig { + k: usize::MAX, + delta: f64::INFINITY, + score_alpha: 0.8, + column_order: None, + }, + &fixture.expected_unpruned, + "unpruned", + ); + assert_expected_results( + &fixture, + FrontierConfig { + k: fixture.pruned.k, + delta: fixture.pruned.delta, + score_alpha: 0.8, + column_order: None, + }, + &fixture.expected_pruned, + "pruned", + ); + } +} diff --git a/python/pecos-rslib-exp/Cargo.toml b/python/pecos-rslib-exp/Cargo.toml index bddcee29a..fe1dd707a 100644 --- a/python/pecos-rslib-exp/Cargo.toml +++ b/python/pecos-rslib-exp/Cargo.toml @@ -27,6 +27,7 @@ extension-module = [ [dependencies] pecos-core.workspace = true pecos-eeg.workspace = true +pecos-frontier.workspace = true pecos-neo.workspace = true pecos-qec.workspace = true pecos-quantum.workspace = true diff --git a/python/pecos-rslib-exp/src/frontier_bindings.rs b/python/pecos-rslib-exp/src/frontier_bindings.rs new file mode 100644 index 000000000..129ebc44f --- /dev/null +++ b/python/pecos-rslib-exp/src/frontier_bindings.rs @@ -0,0 +1,392 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +use pecos_frontier::{ + CommitteeDirection, CommitteeMember, CommitteeStatus, DecoderError, + FrontierCommittee as RustFrontierCommittee, + FrontierCommitteeResult as RustFrontierCommitteeResult, FrontierConfig as RustFrontierConfig, + FrontierDecoder as RustFrontierDecoder, FrontierResult as RustFrontierResult, ObsMask, + SparseDem, backward_deadline_column_order, deadline_column_order, +}; +use pyo3::Borrowed; +use pyo3::exceptions::{PyRuntimeError, PyValueError}; +use pyo3::prelude::*; +use pyo3::types::{PyAny, PyBytes, PyInt}; + +enum ColumnOrderArgument { + Name(String), + Explicit(Vec), +} + +impl<'a, 'py> FromPyObject<'a, 'py> for ColumnOrderArgument { + type Error = PyErr; + + fn extract(object: Borrowed<'a, 'py, PyAny>) -> PyResult { + if let Ok(name) = object.extract::() { + return Ok(Self::Name(name)); + } + if let Ok(order) = object.extract::>() { + return Ok(Self::Explicit(order)); + } + Err(PyValueError::new_err( + "column_order must be 'deadline_reorder', 'time_order', \ + 'backward_deadline_reorder', or a list of mechanism indices", + )) + } +} + +impl Default for ColumnOrderArgument { + fn default() -> Self { + Self::Name("deadline_reorder".to_owned()) + } +} + +fn runtime_error(error: &DecoderError) -> PyErr { + PyRuntimeError::new_err(error.to_string()) +} + +fn sparse_index_error(index: u64, num_detectors: usize) -> PyErr { + PyRuntimeError::new_err(format!( + "Invalid node index {index}: must be < {num_detectors}" + )) +} + +fn sparse_to_dense(indices: &[u64], num_detectors: usize) -> PyResult> { + let mut syndrome = vec![0; num_detectors]; + for &index in indices { + let detector = + usize::try_from(index).map_err(|_| sparse_index_error(index, num_detectors))?; + let bit = syndrome + .get_mut(detector) + .ok_or_else(|| sparse_index_error(index, num_detectors))?; + *bit = 1; + } + Ok(syndrome) +} + +fn resolve_column_order( + dem: &SparseDem, + column_order: ColumnOrderArgument, +) -> PyResult>> { + match column_order { + ColumnOrderArgument::Name(name) => match name.as_str() { + "deadline_reorder" => deadline_column_order(dem) + .map(Some) + .map_err(|e| runtime_error(&e)), + "time_order" => Ok(None), + "backward_deadline_reorder" => backward_deadline_column_order(dem) + .map(Some) + .map_err(|e| runtime_error(&e)), + _ => Err(PyValueError::new_err(format!( + "invalid column_order {name:?}; expected 'deadline_reorder', 'time_order', \ + 'backward_deadline_reorder', or a list of mechanism indices" + ))), + }, + ColumnOrderArgument::Explicit(order) => Ok(Some(order)), + } +} + +fn parse_dem_and_config( + dem_str: &str, + k: usize, + delta: f64, + score_alpha: f64, + column_order: ColumnOrderArgument, +) -> PyResult<(SparseDem, RustFrontierConfig)> { + let dem = SparseDem::from_dem_str(dem_str).map_err(|e| runtime_error(&e))?; + let column_order = resolve_column_order(&dem, column_order)?; + Ok(( + dem, + RustFrontierConfig { + k, + delta, + score_alpha, + column_order, + }, + )) +} + +fn obs_mask_to_py(py: Python<'_>, mask: &ObsMask) -> PyResult> { + if let Some(value) = mask.to_u64() { + return Ok(value.into_pyobject(py)?.into_any().unbind()); + } + + let mut bytes = Vec::with_capacity(std::mem::size_of_val(mask.words())); + for &word in mask.words() { + bytes.extend_from_slice(&word.to_le_bytes()); + } + let py_bytes = PyBytes::new(py, &bytes); + Ok(py + .get_type::() + .call_method1("from_bytes", (py_bytes, "little"))? + .unbind()) +} + +fn observable_bits(mask: &ObsMask, num_observables: usize) -> Vec { + (0..num_observables) + .map(|observable| i32::from(mask.get(observable))) + .collect() +} + +fn logical_masses(py: Python<'_>, result: &RustFrontierResult) -> PyResult, f64)>> { + result + .logical_masses + .iter() + .map(|mass| Ok((obs_mask_to_py(py, &mass.logical)?, mass.log_mass))) + .collect() +} + +fn member_log_evidence(member: CommitteeMember) -> Option { + match member.status { + CommitteeStatus::Ok => Some(member.log_evidence), + CommitteeStatus::NoPath => None, + } +} + +/// Result returned by the native experimental Frontier decoder. +#[pyclass(name = "FrontierResult", module = "pecos_rslib_exp")] +pub struct PyFrontierResult { + inner: RustFrontierResult, +} + +#[pymethods] +impl PyFrontierResult { + /// Predicted observable flips as an arbitrary-precision Python integer. + #[getter] + fn observables_mask(&self, py: Python<'_>) -> PyResult> { + obs_mask_to_py(py, &self.inner.predicted) + } + + /// Total retained log evidence. + #[getter] + fn log_evidence(&self) -> f64 { + self.inner.log_evidence + } + + /// Winner-to-runner-up terminal log-mass gap, if a runner-up exists. + #[getter] + fn runner_up_gap(&self) -> Option { + self.inner.runner_up_gap + } + + /// Peak number of retained dynamic-programming states. + #[getter] + fn peak_retained_states(&self) -> usize { + self.inner.peak_retained_states + } + + /// Number of probabilistic columns processed. + #[getter] + fn processed_columns(&self) -> usize { + self.inner.processed_columns + } + + /// Terminal logical labels and their unnormalized joint log masses. + #[getter] + fn logical_masses(&self, py: Python<'_>) -> PyResult, f64)>> { + logical_masses(py, &self.inner) + } + + /// Return the predicted observable flips as a dense bit vector. + fn observable_bits(&self, num_observables: usize) -> Vec { + observable_bits(&self.inner.predicted, num_observables) + } +} + +/// Result returned by the native experimental forward/backward committee. +#[pyclass(name = "FrontierCommitteeResult", module = "pecos_rslib_exp")] +pub struct PyFrontierCommitteeResult { + inner: RustFrontierCommitteeResult, +} + +#[pymethods] +impl PyFrontierCommitteeResult { + /// Predicted observable flips from the selected leg. + #[getter] + fn observables_mask(&self, py: Python<'_>) -> PyResult> { + obs_mask_to_py(py, &self.inner.selected.predicted) + } + + /// Total retained log evidence from the selected leg. + #[getter] + fn log_evidence(&self) -> f64 { + self.inner.selected.log_evidence + } + + /// Winner-to-runner-up terminal log-mass gap from the selected leg. + #[getter] + fn runner_up_gap(&self) -> Option { + self.inner.selected.runner_up_gap + } + + /// Peak number of states retained by the selected leg. + #[getter] + fn peak_retained_states(&self) -> usize { + self.inner.selected.peak_retained_states + } + + /// Number of probabilistic columns processed by the selected leg. + #[getter] + fn processed_columns(&self) -> usize { + self.inner.selected.processed_columns + } + + /// Terminal logical labels and log masses from the selected leg. + #[getter] + fn logical_masses(&self, py: Python<'_>) -> PyResult, f64)>> { + logical_masses(py, &self.inner.selected) + } + + /// Selected committee direction: `"forward"` or `"backward"`. + #[getter] + fn direction(&self) -> &'static str { + match self.inner.direction { + CommitteeDirection::Forward => "forward", + CommitteeDirection::Backward => "backward", + } + } + + /// Forward-leg evidence, or `None` when that leg found no path. + #[getter] + fn forward_log_evidence(&self) -> Option { + member_log_evidence(self.inner.forward) + } + + /// Backward-leg evidence, or `None` when that leg found no path. + #[getter] + fn backward_log_evidence(&self) -> Option { + member_log_evidence(self.inner.backward) + } + + /// Return the selected observable flips as a dense bit vector. + fn observable_bits(&self, num_observables: usize) -> Vec { + observable_bits(&self.inner.selected.predicted, num_observables) + } +} + +/// Native implementation of the Frontier decoder (Leverrier & Urbanke, +/// arXiv:2606.20513). This experimental decoder is implemented in Rust and does +/// not wrap the upstream `frontier` package. +#[pyclass(name = "FrontierDecoder", module = "pecos_rslib_exp", unsendable)] +pub struct PyFrontierDecoder { + inner: RustFrontierDecoder, + num_detectors: usize, +} + +#[pymethods] +impl PyFrontierDecoder { + /// Construct a native Frontier decoder from a Stim-format DEM string. + #[staticmethod] + #[pyo3( + signature = (dem, *, k=64, delta=50.0, score_alpha=0.8, column_order=ColumnOrderArgument::default()), + text_signature = "(dem, *, k=64, delta=50.0, score_alpha=0.8, column_order='deadline_reorder')" + )] + fn from_dem( + dem: &str, + k: usize, + delta: f64, + score_alpha: f64, + column_order: ColumnOrderArgument, + ) -> PyResult { + let (dem, config) = parse_dem_and_config(dem, k, delta, score_alpha, column_order)?; + let num_detectors = dem.num_detectors; + let inner = + RustFrontierDecoder::from_sparse_dem(&dem, config).map_err(|e| runtime_error(&e))?; + Ok(Self { + inner, + num_detectors, + }) + } + + /// Decode sparse fired-detector indices. + fn decode(&mut self, detection_events: Vec) -> PyResult { + let syndrome = sparse_to_dense(&detection_events, self.num_detectors)?; + self.decode_syndrome(syndrome) + } + + /// Decode one dense detector syndrome. + fn decode_syndrome(&mut self, syndrome: Vec) -> PyResult { + self.inner + .decode(&syndrome) + .map(|inner| PyFrontierResult { inner }) + .map_err(|e| runtime_error(&e)) + } + + /// Decode a batch of dense detector syndromes in input order. + fn decode_batch(&mut self, shots: Vec>) -> PyResult> { + shots + .into_iter() + .map(|syndrome| self.decode_syndrome(syndrome)) + .collect() + } +} + +/// Native experimental forward/backward committee of two Frontier decoders. +/// The implementation is Rust-native and does not wrap the upstream +/// `frontier` package (Leverrier & Urbanke, arXiv:2606.20513). +#[pyclass( + name = "FrontierCommitteeDecoder", + module = "pecos_rslib_exp", + unsendable +)] +pub struct PyFrontierCommitteeDecoder { + inner: RustFrontierCommittee, + num_detectors: usize, +} + +#[pymethods] +impl PyFrontierCommitteeDecoder { + /// Construct a native forward/backward committee from a Stim-format DEM. + #[staticmethod] + #[pyo3( + signature = (dem, *, k=64, delta=50.0, score_alpha=0.8, column_order=ColumnOrderArgument::default()), + text_signature = "(dem, *, k=64, delta=50.0, score_alpha=0.8, column_order='deadline_reorder')" + )] + fn from_dem( + dem: &str, + k: usize, + delta: f64, + score_alpha: f64, + column_order: ColumnOrderArgument, + ) -> PyResult { + let (dem, config) = parse_dem_and_config(dem, k, delta, score_alpha, column_order)?; + let num_detectors = dem.num_detectors; + let inner = + RustFrontierCommittee::from_sparse_dem(&dem, config).map_err(|e| runtime_error(&e))?; + Ok(Self { + inner, + num_detectors, + }) + } + + /// Decode sparse fired-detector indices with both committee legs. + fn decode(&mut self, detection_events: Vec) -> PyResult { + let syndrome = sparse_to_dense(&detection_events, self.num_detectors)?; + self.decode_syndrome(syndrome) + } + + /// Decode one dense detector syndrome with both committee legs. + fn decode_syndrome(&mut self, syndrome: Vec) -> PyResult { + self.inner + .decode(&syndrome) + .map(|inner| PyFrontierCommitteeResult { inner }) + .map_err(|e| runtime_error(&e)) + } + + /// Decode a batch of dense detector syndromes in input order. + fn decode_batch(&mut self, shots: Vec>) -> PyResult> { + shots + .into_iter() + .map(|syndrome| self.decode_syndrome(syndrome)) + .collect() + } +} diff --git a/python/pecos-rslib-exp/src/lib.rs b/python/pecos-rslib-exp/src/lib.rs index e3ea3b0e5..2e78b6c3c 100644 --- a/python/pecos-rslib-exp/src/lib.rs +++ b/python/pecos-rslib-exp/src/lib.rs @@ -24,13 +24,14 @@ clippy::unused_self )] -//! Python bindings for experimental PECOS simulators. +//! Python bindings for experimental PECOS components. //! //! Exposes `StabMps` (stabilizer + MPS hybrid) and `Mast` (magic state -//! injection) from `pecos-stab-tn` via `PyO3`. +//! injection) from `pecos-stab-tn`, plus native Frontier decoders, via `PyO3`. mod coherent_idle_channel; mod eeg_bindings; +mod frontier_bindings; mod mast_bindings; mod sim_neo_bindings; mod stab_mps_bindings; @@ -64,6 +65,10 @@ pub(crate) fn extract_angle( #[pymodule] fn pecos_rslib_exp(m: &Bound<'_, PyModule>) -> PyResult<()> { + m.add_class::()?; + m.add_class::()?; + m.add_class::()?; + m.add_class::()?; m.add_class::()?; m.add_class::()?; m.add_class::()?; diff --git a/python/pecos-rslib/Cargo.toml b/python/pecos-rslib/Cargo.toml index a252957ed..b52021117 100644 --- a/python/pecos-rslib/Cargo.toml +++ b/python/pecos-rslib/Cargo.toml @@ -105,6 +105,7 @@ nalgebra.workspace = true num-complex.workspace = true parking_lot.workspace = true serde_json.workspace = true +sha2.workspace = true tempfile.workspace = true log.workspace = true libc.workspace = true diff --git a/python/pecos-rslib/src/fault_tolerance_bindings.rs b/python/pecos-rslib/src/fault_tolerance_bindings.rs index fbb855d8b..a0475b6bf 100644 --- a/python/pecos-rslib/src/fault_tolerance_bindings.rs +++ b/python/pecos-rslib/src/fault_tolerance_bindings.rs @@ -74,11 +74,27 @@ use pyo3::prelude::*; use std::collections::BTreeMap; use std::str::FromStr; +mod decoder_comparison; +mod decoder_scoring; +mod sample_corpus; + +use decoder_comparison::{ + PyDecoderComparisonResult, compare_decoder_outcomes, validate_comparison_arguments, +}; +use decoder_scoring::{ + MaskedObservableDecoder, ShotDecodeError, TimedObservableDecoder, count_decoder_mismatches, +}; +use sample_corpus::{CorpusError, CorpusToSave, LoadedCorpus}; + type PyDemMechanismTuple = (f64, Vec, Vec); type PyDemFitResult = (Vec, Vec); /// Per-shot detector rows paired with per-shot observable/DEM-output rows. type PyDetectorObservableRows = (Vec>, Vec>); +fn map_shot_decode_error(error: ShotDecodeError) -> PyErr { + pyo3::exceptions::PyRuntimeError::new_err(error.to_string()) +} + fn parse_p1_weights(weights: BTreeMap) -> PyResult { use pecos_core::pauli::{X, Y, Z}; @@ -3368,6 +3384,11 @@ pub struct PySampleBatch { obs_columns: Vec>, num_detectors: usize, num_shots: usize, + seed: Option, + dem: Option, + metadata_json: Option, + generator: Option, + format_version: Option, } impl PySampleBatch { @@ -3440,6 +3461,7 @@ impl PySampleBatch { det_columns: Vec>, obs_columns: Vec>, num_shots: usize, + seed: Option, ) -> Self { let num_detectors = det_columns.len(); Self { @@ -3447,6 +3469,11 @@ impl PySampleBatch { obs_columns, num_detectors, num_shots, + seed, + dem: None, + metadata_json: None, + generator: None, + format_version: None, } } @@ -3493,6 +3520,54 @@ impl PySampleBatch { obs_columns, num_detectors, num_shots, + seed: None, + dem: None, + metadata_json: None, + generator: None, + format_version: None, + } + } + + fn from_corpus(corpus: LoadedCorpus) -> Self { + Self { + num_detectors: corpus.det_columns.len(), + det_columns: corpus.det_columns, + obs_columns: corpus.obs_columns, + num_shots: corpus.num_shots, + seed: corpus.seed, + dem: Some(corpus.dem), + metadata_json: corpus.metadata_json, + generator: Some(corpus.generator), + format_version: Some(corpus.format_version), + } + } + + fn ensure_dem_matches(&self, dem: &str, allow_dem_mismatch: bool) -> PyResult<()> { + if allow_dem_mismatch { + return Ok(()); + } + if let Some(embedded_dem) = &self.dem + && embedded_dem != dem + { + return Err(pyo3::exceptions::PyValueError::new_err( + "supplied DEM differs from the DEM embedded in this loaded SampleBatch; pass \ + allow_dem_mismatch=True to use a different model deliberately", + )); + } + Ok(()) + } + + fn map_corpus_error(error: CorpusError, path: &std::path::Path) -> PyErr { + match error { + CorpusError::Io(error) => match error.raw_os_error() { + Some(errno) => pyo3::exceptions::PyOSError::new_err(( + errno, + error.to_string(), + path.as_os_str().to_os_string(), + )), + None => error.into(), + }, + CorpusError::Invalid(message) => pyo3::exceptions::PyValueError::new_err(message), } } } @@ -3542,6 +3617,99 @@ impl PySampleBatch { self.num_shots } + /// Resolved random seed used to generate this batch, if known. + #[getter] + const fn seed(&self) -> Option { + self.seed + } + + /// Exact detector error model stored with a loaded corpus, if any. + #[getter] + fn dem(&self) -> Option<&str> { + self.dem.as_deref() + } + + /// Opaque caller metadata JSON stored with a loaded corpus, if any. + #[getter] + fn metadata_json(&self) -> Option<&str> { + self.metadata_json.as_deref() + } + + /// PECOS writer identity stored with a loaded corpus, if any. + #[getter] + fn generator(&self) -> Option<&str> { + self.generator.as_deref() + } + + /// Corpus format version for a loaded batch, if any. + #[getter] + const fn format_version(&self) -> Option { + self.format_version + } + + /// Save this serially captured shot batch as a self-describing corpus. + /// + /// Args: + /// path: Destination file path. + /// dem: DEM text associated with the samples. For a generated or + /// Python-constructed batch, only detector and observable dimensions + /// can be checked. This catches gross mismatches, but cannot prove DEM + /// identity or detect a different model with the same dimensions. For + /// a loaded corpus, the text must exactly match its embedded DEM unless + /// `allow_dem_mismatch` is true. + /// `metadata_json`: Optional syntactically valid JSON string. ``None`` + /// preserves metadata already carried by a loaded batch. A supplied + /// value replaces it. + /// `clear_metadata`: Explicitly omit metadata when true. Cannot be combined + /// with a supplied `metadata_json` value. + /// `allow_dem_mismatch`: Permit deliberately saving a loaded batch with a + /// DEM different from its embedded model. + /// + /// Corpora contain shots captured by the serial ``generate_samples`` path. + /// Parallel sample-and-decode paths discard individual shots and cannot be + /// captured by this API. + #[pyo3(signature = (path, *, dem, metadata_json=None, clear_metadata=false, allow_dem_mismatch=false))] + fn save( + &self, + path: std::path::PathBuf, + dem: &str, + metadata_json: Option<&str>, + clear_metadata: bool, + allow_dem_mismatch: bool, + ) -> PyResult<()> { + self.ensure_dem_matches(dem, allow_dem_mismatch)?; + if clear_metadata && metadata_json.is_some() { + return Err(pyo3::exceptions::PyValueError::new_err( + "metadata_json and clear_metadata=True are mutually exclusive", + )); + } + let metadata_json = if clear_metadata { + None + } else { + metadata_json.or(self.metadata_json.as_deref()) + }; + sample_corpus::save( + &path, + CorpusToSave { + det_columns: &self.det_columns, + obs_columns: &self.obs_columns, + num_shots: self.num_shots, + seed: self.seed, + dem, + metadata_json, + }, + ) + .map_err(|error| Self::map_corpus_error(error, &path)) + } + + /// Load and validate a self-describing shot corpus. + #[staticmethod] + fn load(path: std::path::PathBuf) -> PyResult { + sample_corpus::load(&path) + .map(Self::from_corpus) + .map_err(|error| Self::map_corpus_error(error, &path)) + } + /// Get the syndrome for shot `i` as a list of u8 values. fn get_syndrome(&self, i: usize) -> PyResult> { if i >= self.num_shots { @@ -3588,27 +3756,34 @@ impl PySampleBatch { /// `decoder_type`: "pymatching", "`pymatching_correlated`", /// "`pymatching_uncorrelated`", "tesseract", "`bp_osd`", /// "`bp_lsd`", "`union_find`", "`relay_bp`", or "`min_sum_bp`". + /// `allow_dem_mismatch`: Permit a DEM different from the one embedded in + /// a loaded corpus. /// /// Returns: /// Number of logical errors. - #[pyo3(signature = (dem, decoder_type="pymatching"))] - fn decode_count(&self, dem: &str, decoder_type: &str) -> PyResult { + /// + /// Decoder failures abort with a `RuntimeError`. Callers that need errors + /// reported per shot as a distinct outcome should use `compare_decoders`. + #[pyo3(signature = (dem, decoder_type="pymatching", *, allow_dem_mismatch=false))] + fn decode_count( + &self, + dem: &str, + decoder_type: &str, + allow_dem_mismatch: bool, + ) -> PyResult { + self.ensure_dem_matches(dem, allow_dem_mismatch)?; let mut decoder = create_observable_decoder(dem, decoder_type)?; - let mut errors = 0usize; let mut syndrome = vec![0u8; self.num_detectors]; - for i in 0..self.num_shots { - self.extract_syndrome(i, &mut syndrome); - // Wide ObsMask comparison: inline (one stack word) for the typical - // <=64 observables, correct without truncation beyond. A decode - // failure counts as a logical error (matching the prior sentinel). - let is_error = decoder - .decode_obs(&syndrome) - .map_or(true, |p| p != self.extract_obs_mask_wide(i)); - if is_error { - errors += 1; - } - } - Ok(errors) + count_decoder_mismatches( + 0..self.num_shots, + &mut syndrome, + |shot, buffer| { + self.extract_syndrome(shot, buffer); + self.extract_obs_mask_wide(shot) + }, + decoder.as_mut(), + ) + .map_err(map_shot_decode_error) } /// Decode every shot and return the predicted observable mask per shot. @@ -3620,32 +3795,86 @@ impl PySampleBatch { /// Args: /// dem: DEM string for the decoder. /// `decoder_type`: Decoder type string. + /// `allow_dem_mismatch`: Permit a DEM different from the one embedded in + /// a loaded corpus. /// /// Returns: /// List of predicted observable masks (Python ints; arbitrary precision, /// so more than 64 observables are not truncated), one per shot. - #[pyo3(signature = (dem, decoder_type="pymatching"))] + #[pyo3(signature = (dem, decoder_type="pymatching", *, allow_dem_mismatch=false))] fn decode_each( &self, py: Python<'_>, dem: &str, decoder_type: &str, + allow_dem_mismatch: bool, ) -> PyResult>> { + self.ensure_dem_matches(dem, allow_dem_mismatch)?; let mut decoder = create_observable_decoder(dem, decoder_type)?; let mut predictions = Vec::with_capacity(self.num_shots); let mut syndrome = vec![0u8; self.num_detectors]; - for i in 0..self.num_shots { - self.extract_syndrome(i, &mut syndrome); + for shot in 0..self.num_shots { + self.extract_syndrome(shot, &mut syndrome); // Propagate a decode failure rather than masking it as a sentinel // observable value (which would read as a spurious disagreement). - let predicted = decoder - .decode_obs(&syndrome) - .map_err(|e| PyErr::new::(e.to_string()))?; + let predicted = decoder.decode_obs(&syndrome).map_err(|error| { + pyo3::exceptions::PyRuntimeError::new_err(format!( + "decoder failed on shot {shot}: {error}" + )) + })?; predictions.push(obsmask_to_py(py, &predicted)?); } Ok(predictions) } + /// Decode every shot with a decoder under test (DUT) and a reference decoder. + /// + /// Both decoders receive the same shots in the same order. Each result is + /// independently classified as correct, mismatch, or decode error, and a + /// decode error is counted for that shot without aborting the comparison. + /// Predictions and truth are compared as wide observable masks, with no + /// 64-observable limit. + /// + /// Args: + /// dem: DEM string shared by both decoders. + /// `dut_decoder_type`: Decoder type string for the decoder under test. + /// `reference_decoder_type`: Decoder type string for the reference. + /// alpha: Tail probability for equal-tailed Jeffreys intervals. + /// `allow_dem_mismatch`: Permit deliberate cross-model comparison of a + /// loaded corpus. + /// + /// Returns: + /// A `DecoderComparisonResult` containing the raw 3x3 counts and + /// headline DUT-only-failure and both-failed proportions. + #[pyo3(signature = (dem, dut_decoder_type, reference_decoder_type, alpha=0.05, *, allow_dem_mismatch=false))] + fn compare_decoders( + &self, + dem: &str, + dut_decoder_type: &str, + reference_decoder_type: &str, + alpha: f64, + allow_dem_mismatch: bool, + ) -> PyResult { + validate_comparison_arguments(self.num_shots, alpha) + .map_err(|error| pyo3::exceptions::PyValueError::new_err(error.to_string()))?; + self.ensure_dem_matches(dem, allow_dem_mismatch)?; + let mut dut = create_observable_decoder(dem, dut_decoder_type)?; + let mut reference = create_observable_decoder(dem, reference_decoder_type)?; + let mut syndrome = vec![0u8; self.num_detectors]; + let counts = compare_decoder_outcomes( + self.num_shots, + &mut syndrome, + |shot, buffer| { + self.extract_syndrome(shot, buffer); + self.extract_obs_mask_wide(shot) + }, + dut.as_mut(), + reference.as_mut(), + ); + PyDecoderComparisonResult::new(counts, alpha) + .map_err(|error| pyo3::exceptions::PyRuntimeError::new_err(error.to_string())) + } + /// Parallel decode: distributes samples across rayon workers. /// /// Each worker creates its own decoder instance. Faster for slow decoders. @@ -3654,18 +3883,29 @@ impl PySampleBatch { /// dem: DEM string for the decoder. /// `decoder_type`: Decoder type string. /// `num_workers`: Number of parallel workers (default: number of CPUs). + /// `allow_dem_mismatch`: Permit a DEM different from the one embedded in + /// a loaded corpus. /// /// Returns: /// Number of logical errors. - #[pyo3(signature = (dem, decoder_type="pymatching", num_workers=None))] + /// + /// Decoder failures abort with a `RuntimeError`. Callers that need errors + /// reported per shot as a distinct outcome should use `compare_decoders`. + /// + /// Set `allow_dem_mismatch` to true to use a DEM different from the one + /// embedded in a loaded corpus. + #[pyo3(signature = (dem, decoder_type="pymatching", num_workers=None, *, allow_dem_mismatch=false))] fn decode_count_parallel( &self, dem: &str, decoder_type: &str, num_workers: Option, + allow_dem_mismatch: bool, ) -> PyResult { use rayon::prelude::*; + self.ensure_dem_matches(dem, allow_dem_mismatch)?; + drop(create_observable_decoder(dem, decoder_type)?); let n_workers = num_workers.unwrap_or_else(rayon::current_num_threads); let pool = rayon::ThreadPoolBuilder::new() .num_threads(n_workers) @@ -3688,23 +3928,40 @@ impl PySampleBatch { let observable_masks: Vec = (0..n).map(|i| self.extract_obs_mask_wide(i)).collect(); - let total_errors: usize = pool.install(|| { - (0..n) + let worker_results: Vec> = pool.install(|| { + let chunk_size = n.div_ceil(n_workers); + (0..n_workers) .into_par_iter() - .map_init( - || create_observable_decoder(&dem_str, &dt).unwrap(), - |decoder, i| { - usize::from( - decoder - .decode_obs(&detection_events[i]) - .map_or(true, |p| p != observable_masks[i]), - ) - }, - ) - .sum() + .map(|worker_id| { + let start = worker_id * chunk_size; + let end = (start + chunk_size).min(n); + if start >= end { + return Ok(0); + } + + // Safe after the identical factory call was validated above. + let mut decoder = create_observable_decoder(&dem_str, &dt).unwrap(); + let mut syndrome = vec![0u8; num_dets]; + count_decoder_mismatches( + start..end, + &mut syndrome, + |shot, buffer| { + buffer.copy_from_slice(&detection_events[shot]); + observable_masks[shot].clone() + }, + decoder.as_mut(), + ) + }) + .collect() }); - Ok(total_errors) + worker_results + .into_iter() + .try_fold(0usize, |total, result| { + result + .map(|count| total + count) + .map_err(map_shot_decode_error) + }) } /// Batch decode all samples at once using `PyMatching`'s batch API. @@ -3715,10 +3972,15 @@ impl PySampleBatch { /// /// Returns: /// Number of logical errors. - #[pyo3(signature = (dem))] - fn decode_count_batch(&self, dem: &str) -> PyResult { + /// + /// Batch decoder failures abort with a `RuntimeError`. Callers that need + /// errors reported per shot as a distinct outcome should use + /// `compare_decoders`. + #[pyo3(signature = (dem, *, allow_dem_mismatch=false))] + fn decode_count_batch(&self, dem: &str, allow_dem_mismatch: bool) -> PyResult { use pecos_decoders::{BatchConfig, PyMatchingDecoder}; + self.ensure_dem_matches(dem, allow_dem_mismatch)?; let mut decoder = PyMatchingDecoder::from_dem(dem) .map_err(|e| PyErr::new::(e.to_string()))?; @@ -3774,28 +4036,36 @@ impl PySampleBatch { /// Args: /// dem: DEM string for the decoder. /// `decoder_type`: Decoder type string. + /// `allow_dem_mismatch`: Permit a DEM different from the one embedded in + /// a loaded corpus. /// /// Returns: /// `DecodeStats` with timing breakdown. - #[pyo3(signature = (dem, decoder_type="pymatching"))] - fn decode_stats(&self, dem: &str, decoder_type: &str) -> PyResult { - use std::time::Instant; - - let mut decoder = create_observable_decoder(dem, decoder_type)?; - let mut num_errors = 0usize; - let mut per_shot_seconds: Vec = Vec::with_capacity(self.num_shots); + /// + /// Decoder failures abort with a `RuntimeError`. Callers that need errors + /// reported per shot as a distinct outcome should use `compare_decoders`. + #[pyo3(signature = (dem, decoder_type="pymatching", *, allow_dem_mismatch=false))] + fn decode_stats( + &self, + dem: &str, + decoder_type: &str, + allow_dem_mismatch: bool, + ) -> PyResult { + self.ensure_dem_matches(dem, allow_dem_mismatch)?; + let decoder = create_observable_decoder(dem, decoder_type)?; + let mut decoder = TimedObservableDecoder::new(decoder, self.num_shots); let mut syndrome = vec![0u8; self.num_detectors]; - - for i in 0..self.num_shots { - self.extract_syndrome(i, &mut syndrome); - let t0 = Instant::now(); - let predicted = decoder.decode_obs(&syndrome); - let elapsed = t0.elapsed().as_secs_f64(); - per_shot_seconds.push(elapsed); - if predicted.map_or(true, |p| p != self.extract_obs_mask_wide(i)) { - num_errors += 1; - } - } + let num_errors = count_decoder_mismatches( + 0..self.num_shots, + &mut syndrome, + |shot, buffer| { + self.extract_syndrome(shot, buffer); + self.extract_obs_mask_wide(shot) + }, + &mut decoder, + ) + .map_err(map_shot_decode_error)?; + let per_shot_seconds = decoder.into_times(); Ok(PyDecodeStats::from_times( self.num_shots, @@ -3817,15 +4087,22 @@ impl PySampleBatch { /// dem: DEM string for the decoder. /// `decoder_type`: Decoder type string. /// `num_workers`: Number of parallel workers (default: number of CPUs). - #[pyo3(signature = (dem, decoder_type="mwpf", num_workers=None))] + /// `allow_dem_mismatch`: Permit a DEM different from the one embedded in + /// a loaded corpus. + /// + /// Decoder failures abort with a `RuntimeError`. Callers that need errors + /// reported per shot as a distinct outcome should use `compare_decoders`. + #[pyo3(signature = (dem, decoder_type="mwpf", num_workers=None, *, allow_dem_mismatch=false))] fn decode_stats_parallel( &self, dem: &str, decoder_type: &str, num_workers: Option, + allow_dem_mismatch: bool, ) -> PyResult { use rayon::prelude::*; + self.ensure_dem_matches(dem, allow_dem_mismatch)?; let n_workers = num_workers.unwrap_or_else(rayon::current_num_threads); // Validate decoder type early. @@ -3852,8 +4129,8 @@ impl PySampleBatch { .map(|i| self.extract_obs_mask_wide(i)) .collect(); - // Each worker decodes a slice of shots and returns (errors, per_shot_times). - let results: Vec<(usize, Vec)> = pool.install(|| { + // Each worker decodes a contiguous slice and returns its count and timings. + let results: Vec), ShotDecodeError>> = pool.install(|| { let chunk_size = self.num_shots.div_ceil(n_workers); (0..n_workers) .into_par_iter() @@ -3861,29 +4138,31 @@ impl PySampleBatch { let start = worker_id * chunk_size; let end = (start + chunk_size).min(self.num_shots); if start >= end { - return (0, Vec::new()); + return Ok((0, Vec::new())); } - let mut decoder = create_observable_decoder(&dem_str, &dt).unwrap(); - let mut errors = 0usize; - let mut times = Vec::with_capacity(end - start); - - for i in start..end { - let t0 = std::time::Instant::now(); - let predicted = decoder.decode_obs(&detection_events[i]); - times.push(t0.elapsed().as_secs_f64()); - if predicted.map_or(true, |p| p != observable_masks[i]) { - errors += 1; - } - } - (errors, times) + // Safe after the identical factory call was validated above. + let decoder = create_observable_decoder(&dem_str, &dt).unwrap(); + let mut decoder = TimedObservableDecoder::new(decoder, end - start); + let mut syndrome = vec![0u8; num_dets]; + let errors = count_decoder_mismatches( + start..end, + &mut syndrome, + |shot, buffer| { + buffer.copy_from_slice(&detection_events[shot]); + observable_masks[shot].clone() + }, + &mut decoder, + )?; + Ok((errors, decoder.into_times())) }) .collect() }); let mut total_errors = 0usize; let mut all_times = Vec::with_capacity(self.num_shots); - for (errs, times) in results { + for result in results { + let (errs, times) = result.map_err(map_shot_decode_error)?; total_errors += errs; all_times.extend(times); } @@ -4488,15 +4767,13 @@ impl PyDemSampler { use pecos_random::PecosRng; use rand::RngExt; - let mut rng = match seed { - Some(s) => PecosRng::seed_from_u64(s), - None => PecosRng::seed_from_u64(rand::rng().random()), - }; + let actual_seed = seed.unwrap_or_else(|| rand::rng().random()); + let mut rng = PecosRng::seed_from_u64(actual_seed); // Use geometric columnar sampler via DemSampler. let (det_columns, obs_columns) = self.inner.sample_batch_geometric(num_shots, &mut rng); - PySampleBatch::from_columnar(det_columns, obs_columns, num_shots) + PySampleBatch::from_columnar(det_columns, obs_columns, num_shots, Some(actual_seed)) } /// Compute statistics without storing individual shots. @@ -4609,6 +4886,10 @@ impl PyDemSampler { /// /// Returns: /// Number of logical errors (mismatches between decoder prediction and true flip). + /// + /// Decoder failures abort with a `RuntimeError`. Callers that need errors + /// reported per shot as a distinct outcome should generate a `SampleBatch` + /// and use `SampleBatch.compare_decoders`. #[pyo3(signature = (dem, num_shots, decoder_type="pymatching", seed=None))] fn sample_decode_count( &self, @@ -4623,27 +4904,28 @@ impl PyDemSampler { let actual_seed = seed.unwrap_or_else(|| rand::rng().random()); let mut rng = PecosRng::seed_from_u64(actual_seed); - let mut decoder = create_observable_decoder(dem, decoder_type)?; + let decoder = create_observable_decoder(dem, decoder_type)?; let observable_mask = self.inner.observable_dem_output_mask(); + let mut decoder = MaskedObservableDecoder::new(decoder, observable_mask.clone()); // Tight sample+decode loop -- no Python involvement. // Single-threaded: sample and decode sequentially. - let mut errors = 0usize; - for _ in 0..num_shots { - let (det_events, obs_flips) = self.inner.sample(&mut rng); - let syndrome: Vec = det_events.iter().map(|&b| u8::from(b)).collect(); - let mut predicted = decoder - .decode_obs(&syndrome) - .map_err(|e| PyErr::new::(e.to_string()))?; - predicted &= &observable_mask; - let true_mask = self - .inner - .observable_mask_from_dem_output_flips(&obs_flips, &observable_mask); - if predicted != true_mask { - errors += 1; - } - } - Ok(errors) + let mut syndrome = vec![0u8; self.inner.num_detectors()]; + count_decoder_mismatches( + 0..num_shots, + &mut syndrome, + |_, buffer| { + let (det_events, obs_flips) = self.inner.sample(&mut rng); + debug_assert_eq!(det_events.len(), buffer.len()); + for (value, event) in buffer.iter_mut().zip(det_events) { + *value = u8::from(event); + } + self.inner + .observable_mask_from_dem_output_flips(&obs_flips, &observable_mask) + }, + &mut decoder, + ) + .map_err(map_shot_decode_error) } /// Parallel sample+decode: distributes shots across threads. @@ -4662,6 +4944,10 @@ impl PyDemSampler { /// /// Returns: /// Number of logical errors. + /// + /// Decoder failures abort with a `RuntimeError`. Callers that need errors + /// reported per shot as a distinct outcome should generate a `SampleBatch` + /// and use `SampleBatch.compare_decoders`. #[pyo3(signature = (dem, num_shots, decoder_type="pymatching", seed=None, num_workers=None))] fn sample_decode_count_parallel( &self, @@ -4692,7 +4978,7 @@ impl PyDemSampler { let dem_str = dem.to_string(); let dt = decoder_type.to_string(); - let total_errors: usize = pool.install(|| { + let worker_results: Vec> = pool.install(|| { (0..n_workers) .into_par_iter() .map(|worker_id| { @@ -4700,35 +4986,44 @@ impl PyDemSampler { let my_shots = shots_per_worker + usize::from(worker_id < remainder); if my_shots == 0 { - return 0; + return Ok(0); } + let start = worker_id * shots_per_worker + worker_id.min(remainder); + let end = start + my_shots; let my_sampler = sampler.clone(); let mut my_rng = PecosRng::seed_from_u64(actual_seed.wrapping_add(worker_id as u64)); // unwrap is safe: we validated above - let mut decoder = create_observable_decoder(&dem_str, &dt).unwrap(); - - let mut errors = 0usize; - for _ in 0..my_shots { - let (det_events, obs_flips) = my_sampler.sample(&mut my_rng); - let syndrome: Vec = det_events.iter().map(|&b| u8::from(b)).collect(); - let mut predicted = decoder - .decode_obs(&syndrome) - .unwrap_or_else(|_| observable_mask.clone()); - predicted &= &observable_mask; - let truth = my_sampler - .observable_mask_from_dem_output_flips(&obs_flips, &observable_mask); - if predicted != truth { - errors += 1; - } - } - errors + let decoder = create_observable_decoder(&dem_str, &dt).unwrap(); + let mut decoder = + MaskedObservableDecoder::new(decoder, observable_mask.clone()); + let mut syndrome = vec![0u8; my_sampler.num_detectors()]; + count_decoder_mismatches( + start..end, + &mut syndrome, + |_, buffer| { + let (det_events, obs_flips) = my_sampler.sample(&mut my_rng); + debug_assert_eq!(det_events.len(), buffer.len()); + for (value, event) in buffer.iter_mut().zip(det_events) { + *value = u8::from(event); + } + my_sampler + .observable_mask_from_dem_output_flips(&obs_flips, &observable_mask) + }, + &mut decoder, + ) }) - .sum() + .collect() }); - Ok(total_errors) + worker_results + .into_iter() + .try_fold(0usize, |total, result| { + result + .map(|count| total + count) + .map_err(map_shot_decode_error) + }) } fn __repr__(&self) -> String { @@ -6940,6 +7235,7 @@ pub fn register_qec_module(m: &Bound<'_, PyModule>) -> PyResult<()> { qec.add_class::()?; qec.add_class::()?; qec.add_class::()?; + qec.add_class::()?; qec.add_class::()?; qec.add_class::()?; qec.add_class::()?; diff --git a/python/pecos-rslib/src/fault_tolerance_bindings/decoder_comparison.rs b/python/pecos-rslib/src/fault_tolerance_bindings/decoder_comparison.rs new file mode 100644 index 000000000..7f04031ca --- /dev/null +++ b/python/pecos-rslib/src/fault_tolerance_bindings/decoder_comparison.rs @@ -0,0 +1,486 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +//! Paired DUT/reference decoder comparison over a shared sequence of shots. + +use pecos_decoder_core::obs_mask::ObsMask; +use pecos_decoder_core::{DecoderError, ObservableDecoder}; +use pecos_num::stats::{ + JeffreysError, JeffreysEstimator, JeffreysInterval, jeffreys_interval, jeffreys_point, +}; +use pyo3::prelude::*; + +#[derive(Clone, Copy, Debug, Eq, PartialEq)] +enum DecoderOutcome { + Correct, + Mismatch, + Error, +} + +impl DecoderOutcome { + const fn index(self) -> usize { + match self { + Self::Correct => 0, + Self::Mismatch => 1, + Self::Error => 2, + } + } +} + +fn classify(result: Result, truth: &ObsMask) -> DecoderOutcome { + match result { + Ok(prediction) if prediction == *truth => DecoderOutcome::Correct, + Ok(_) => DecoderOutcome::Mismatch, + Err(_) => DecoderOutcome::Error, + } +} + +/// Validate caller-controlled comparison arguments without running an interval +/// calculation or constructing either decoder. +pub(super) fn validate_comparison_arguments( + num_shots: usize, + alpha: f64, +) -> Result<(), JeffreysError> { + let num_shots = u64::try_from(num_shots).expect("usize always fits in u64"); + // The mean estimator performs the shared zero/maximum-trial validation but + // no special-function solve, keeping this preflight constant-time. + jeffreys_point(0, num_shots, JeffreysEstimator::Mean)?; + if !alpha.is_finite() || alpha <= 0.0 || alpha >= 1.0 { + return Err(JeffreysError::InvalidAlpha { alpha }); + } + Ok(()) +} + +/// Counts indexed by DUT outcome first, then reference outcome. +/// +/// In each dimension the order is correct, mismatch, decode error. +#[derive(Clone, Debug, Default, Eq, PartialEq)] +pub(super) struct DecoderComparisonCounts { + cells: [[u64; 3]; 3], +} + +impl DecoderComparisonCounts { + fn record(&mut self, dut: DecoderOutcome, reference: DecoderOutcome) { + self.cells[dut.index()][reference.index()] += 1; + } + + pub(super) const fn cells(&self) -> &[[u64; 3]; 3] { + &self.cells + } + + fn total_shots(&self) -> u64 { + self.cells.iter().flatten().sum() + } + + const fn dut_only_failures(&self) -> u64 { + self.cells[DecoderOutcome::Mismatch.index()][DecoderOutcome::Correct.index()] + } + + const fn both_failed(&self) -> u64 { + self.cells[DecoderOutcome::Mismatch.index()][DecoderOutcome::Mismatch.index()] + } +} + +/// Compare two decoders on the same shots in the same order. +/// +/// `prepare_shot` writes the selected syndrome into the reusable buffer and +/// returns that shot's wide true-observable mask. +pub(super) fn compare_decoder_outcomes( + num_shots: usize, + syndrome: &mut [u8], + mut prepare_shot: impl FnMut(usize, &mut [u8]) -> ObsMask, + dut: &mut dyn ObservableDecoder, + reference: &mut dyn ObservableDecoder, +) -> DecoderComparisonCounts { + let mut counts = DecoderComparisonCounts::default(); + for shot in 0..num_shots { + let truth = prepare_shot(shot, syndrome); + // Run both decoders before classifying either result. In particular, a + // DUT error must not prevent the reference from seeing this shot. + let dut_result = dut.decode_obs(syndrome); + let reference_result = reference.decode_obs(syndrome); + counts.record( + classify(dut_result, &truth), + classify(reference_result, &truth), + ); + } + counts +} + +#[derive(Clone, Copy, Debug)] +struct HeadlineProportion { + point: f64, + interval: JeffreysInterval, +} + +impl HeadlineProportion { + fn new(count: u64, total_shots: u64, alpha: f64) -> Result { + let interval = jeffreys_interval(count, total_shots, alpha)?; + Ok(Self { + point: interval.point, + interval, + }) + } +} + +/// Python-facing paired decoder contingency counts and headline proportions. +#[pyclass( + name = "DecoderComparisonResult", + module = "pecos_rslib.qec", + skip_from_py_object +)] +#[derive(Clone, Debug)] +pub(super) struct PyDecoderComparisonResult { + counts: DecoderComparisonCounts, + total_shots: u64, + alpha: f64, + dut_only_failure: HeadlineProportion, + both_failed: HeadlineProportion, +} + +impl PyDecoderComparisonResult { + pub(super) fn new(counts: DecoderComparisonCounts, alpha: f64) -> Result { + let total_shots = counts.total_shots(); + let dut_only_failure = + HeadlineProportion::new(counts.dut_only_failures(), total_shots, alpha)?; + let both_failed = HeadlineProportion::new(counts.both_failed(), total_shots, alpha)?; + Ok(Self { + counts, + total_shots, + alpha, + dut_only_failure, + both_failed, + }) + } +} + +#[pymethods] +impl PyDecoderComparisonResult { + /// Raw 3x3 counts in correct, mismatch, error order on both axes. + #[getter] + fn counts(&self) -> Vec> { + self.counts.cells().iter().map(|row| row.to_vec()).collect() + } + + /// Number of shots compared. + #[getter] + const fn total_shots(&self) -> u64 { + self.total_shots + } + + /// Tail probability used for the equal-tailed Jeffreys intervals. + #[getter] + const fn alpha(&self) -> f64 { + self.alpha + } + + #[getter] + const fn dut_correct_reference_correct(&self) -> u64 { + self.counts.cells[0][0] + } + + #[getter] + const fn dut_correct_reference_mismatch(&self) -> u64 { + self.counts.cells[0][1] + } + + #[getter] + const fn dut_correct_reference_error(&self) -> u64 { + self.counts.cells[0][2] + } + + #[getter] + const fn dut_mismatch_reference_correct(&self) -> u64 { + self.counts.cells[1][0] + } + + #[getter] + const fn dut_mismatch_reference_mismatch(&self) -> u64 { + self.counts.cells[1][1] + } + + #[getter] + const fn dut_mismatch_reference_error(&self) -> u64 { + self.counts.cells[1][2] + } + + #[getter] + const fn dut_error_reference_correct(&self) -> u64 { + self.counts.cells[2][0] + } + + #[getter] + const fn dut_error_reference_mismatch(&self) -> u64 { + self.counts.cells[2][1] + } + + #[getter] + const fn dut_error_reference_error(&self) -> u64 { + self.counts.cells[2][2] + } + + /// DUT mismatches on shots where the reference was correct. + #[getter] + const fn dut_only_failures(&self) -> u64 { + self.counts.dut_only_failures() + } + + /// Jeffreys posterior-mean proportion for DUT-only failures. + #[getter] + const fn dut_only_failure_proportion(&self) -> f64 { + self.dut_only_failure.point + } + + /// Equal-tailed Jeffreys interval for the DUT-only-failure proportion. + #[getter] + const fn dut_only_failure_interval(&self) -> (f64, f64) { + ( + self.dut_only_failure.interval.lo, + self.dut_only_failure.interval.hi, + ) + } + + /// Shots on which both decoders returned mismatching predictions. + #[getter] + const fn both_failed(&self) -> u64 { + self.counts.both_failed() + } + + /// Jeffreys posterior-mean proportion for shots where both decoders failed. + #[getter] + const fn both_failed_proportion(&self) -> f64 { + self.both_failed.point + } + + /// Equal-tailed Jeffreys interval for the both-failed proportion. + #[getter] + const fn both_failed_interval(&self) -> (f64, f64) { + (self.both_failed.interval.lo, self.both_failed.interval.hi) + } + + fn __repr__(&self) -> String { + format!( + "DecoderComparisonResult(shots={}, dut_only_failures={}, both_failed={})", + self.total_shots, + self.counts.dut_only_failures(), + self.counts.both_failed(), + ) + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[derive(Clone, Debug)] + enum StubResult { + Prediction(ObsMask), + Error, + } + + struct StubDecoder { + expected_syndromes: Vec>, + results: Vec, + next: usize, + } + + impl StubDecoder { + fn new(expected_syndromes: &[Vec], results: Vec) -> Self { + assert_eq!(expected_syndromes.len(), results.len()); + Self { + expected_syndromes: expected_syndromes.to_vec(), + results, + next: 0, + } + } + } + + impl ObservableDecoder for StubDecoder { + fn decode_obs(&mut self, syndrome: &[u8]) -> Result { + assert_eq!(syndrome, self.expected_syndromes[self.next]); + let result = match &self.results[self.next] { + StubResult::Prediction(mask) => Ok(mask.clone()), + StubResult::Error => Err(DecoderError::DecodingFailed("stub error".into())), + }; + self.next += 1; + result + } + } + + fn mask(bits: &[usize]) -> ObsMask { + let mut mask = ObsMask::new(); + for &bit in bits { + mask.set(bit); + } + mask + } + + fn predictions(masks: &[ObsMask]) -> Vec { + masks.iter().cloned().map(StubResult::Prediction).collect() + } + + fn compare( + shots: &[(Vec, ObsMask)], + dut_results: Vec, + reference_results: Vec, + ) -> DecoderComparisonCounts { + let syndromes: Vec> = shots.iter().map(|(s, _)| s.clone()).collect(); + let mut dut = StubDecoder::new(&syndromes, dut_results); + let mut reference = StubDecoder::new(&syndromes, reference_results); + let mut syndrome = vec![0; syndromes.first().map_or(0, Vec::len)]; + compare_decoder_outcomes( + shots.len(), + &mut syndrome, + |shot, buffer| { + buffer.copy_from_slice(&shots[shot].0); + shots[shot].1.clone() + }, + &mut dut, + &mut reference, + ) + } + + fn sample_shots() -> Vec<(Vec, ObsMask)> { + vec![ + (vec![0, 0], mask(&[])), + (vec![1, 0], mask(&[0])), + (vec![0, 1], mask(&[1])), + (vec![1, 1], mask(&[0, 1])), + ] + } + + #[test] + fn both_decoders_correct_puts_all_mass_in_correct_correct() { + let shots = sample_shots(); + let truths: Vec = shots.iter().map(|(_, truth)| truth.clone()).collect(); + let counts = compare(&shots, predictions(&truths), predictions(&truths)); + + assert_eq!(counts.cells(), &[[4, 0, 0], [0, 0, 0], [0, 0, 0]]); + assert_eq!(counts.dut_only_failures(), 0); + } + + #[test] + fn dut_only_failures_count_a_known_wrong_subset() { + let shots = sample_shots(); + let truths: Vec = shots.iter().map(|(_, truth)| truth.clone()).collect(); + let mut dut = truths.clone(); + dut[1] = mask(&[]); + dut[3] = mask(&[]); + + let counts = compare(&shots, predictions(&dut), predictions(&truths)); + + // Shots 1 and 3 are deliberately wrong for the DUT: 2 DUT-only failures. + assert_eq!(counts.dut_only_failures(), 2); + assert_eq!(counts.cells(), &[[2, 0, 0], [2, 0, 0], [0, 0, 0]]); + } + + #[test] + fn dut_errors_are_not_mismatches_and_do_not_abort() { + let shots = sample_shots(); + let truths: Vec = shots.iter().map(|(_, truth)| truth.clone()).collect(); + let mut dut = predictions(&truths); + dut[1] = StubResult::Error; + dut[3] = StubResult::Error; + + let counts = compare(&shots, dut, predictions(&truths)); + + assert_eq!(counts.cells(), &[[2, 0, 0], [0, 0, 0], [2, 0, 0]]); + assert_eq!(counts.cells()[DecoderOutcome::Mismatch.index()][0], 0); + } + + #[test] + fn reference_errors_are_counted_and_do_not_abort() { + let shots = sample_shots(); + let truths: Vec = shots.iter().map(|(_, truth)| truth.clone()).collect(); + let mut reference = predictions(&truths); + reference[0] = StubResult::Error; + reference[2] = StubResult::Error; + + let counts = compare(&shots, predictions(&truths), reference); + + assert_eq!(counts.cells(), &[[2, 0, 2], [0, 0, 0], [0, 0, 0]]); + } + + #[test] + fn wide_observable_difference_above_bit_63_is_preserved() { + let wide_truth = mask(&[70]); + let shots = vec![(vec![1], wide_truth.clone())]; + let counts = compare( + &shots, + predictions(&[ObsMask::new()]), + predictions(&[wide_truth]), + ); + + assert_eq!(counts.cells(), &[[0, 0, 0], [1, 0, 0], [0, 0, 0]]); + assert_eq!(counts.dut_only_failures(), 1); + } + + #[test] + fn headline_interval_matches_pecos_num_helper() { + let shots = sample_shots(); + let truths: Vec = shots.iter().map(|(_, truth)| truth.clone()).collect(); + let mut dut = truths.clone(); + dut[1] = mask(&[]); + let summary = PyDecoderComparisonResult::new( + compare(&shots, predictions(&dut), predictions(&truths)), + 0.05, + ) + .expect("valid Jeffreys inputs"); + let expected = jeffreys_interval(1, 4, 0.05).expect("valid direct helper inputs"); + + assert_eq!(summary.dut_only_failure.interval, expected); + // Both sides come from the same helper call, so the point estimate must be + // bit-identical; compare bit patterns rather than floats. + assert_eq!( + summary.dut_only_failure.point.to_bits(), + expected.point.to_bits() + ); + } + + #[test] + fn comparison_arguments_reject_invalid_shot_counts() { + assert_eq!( + validate_comparison_arguments(0, 0.05), + Err(JeffreysError::ZeroTrials) + ); + assert!(matches!( + validate_comparison_arguments(100_000_001, 0.05), + Err(JeffreysError::TrialsExceedSupported { + n: 100_000_001, + max: 100_000_000 + }) + )); + } + + #[test] + fn comparison_arguments_reject_alpha_outside_open_unit_interval() { + for alpha in [f64::NAN, f64::NEG_INFINITY, 0.0, 1.0, f64::INFINITY] { + assert!(matches!( + validate_comparison_arguments(1, alpha), + Err(JeffreysError::InvalidAlpha { .. }) + )); + } + } + + #[test] + fn comparison_is_deterministic_for_the_same_batch() { + let shots = sample_shots(); + let truths: Vec = shots.iter().map(|(_, truth)| truth.clone()).collect(); + let mut dut = truths.clone(); + dut[2] = mask(&[]); + + let first = compare(&shots, predictions(&dut), predictions(&truths)); + let second = compare(&shots, predictions(&dut), predictions(&truths)); + + assert_eq!(first, second); + } +} diff --git a/python/pecos-rslib/src/fault_tolerance_bindings/decoder_scoring.rs b/python/pecos-rslib/src/fault_tolerance_bindings/decoder_scoring.rs new file mode 100644 index 000000000..9a4f88fed --- /dev/null +++ b/python/pecos-rslib/src/fault_tolerance_bindings/decoder_scoring.rs @@ -0,0 +1,234 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +//! Shared fail-loud scoring for observable decoders. + +use pecos_decoder_core::obs_mask::ObsMask; +use pecos_decoder_core::{DecoderError, ObservableDecoder}; +use std::fmt; +use std::ops::Range; + +/// A decoder failure annotated with the shot that caused it. +#[derive(Debug)] +pub(super) struct ShotDecodeError { + shot_index: usize, + source: DecoderError, +} + +impl fmt::Display for ShotDecodeError { + fn fmt(&self, formatter: &mut fmt::Formatter<'_>) -> fmt::Result { + write!( + formatter, + "decoder failed on shot {}: {}", + self.shot_index, self.source + ) + } +} + +impl std::error::Error for ShotDecodeError { + fn source(&self) -> Option<&(dyn std::error::Error + 'static)> { + Some(&self.source) + } +} + +/// Decode and score a contiguous range of shots, aborting on the first failure. +/// +/// `access_shot` fills the reusable syndrome buffer and returns that shot's +/// true observable mask. Shot indices are kept absolute so callers can combine +/// independently scored worker ranges without losing error context. +pub(super) fn count_decoder_mismatches( + shots: Range, + syndrome: &mut [u8], + mut access_shot: impl FnMut(usize, &mut [u8]) -> ObsMask, + decoder: &mut dyn ObservableDecoder, +) -> Result { + let mut mismatches = 0; + for shot_index in shots { + let truth = access_shot(shot_index, syndrome); + let prediction = decoder + .decode_obs(syndrome) + .map_err(|source| ShotDecodeError { shot_index, source })?; + mismatches += usize::from(prediction != truth); + } + Ok(mismatches) +} + +/// Observable-decoder adapter that keeps only caller-selected observables. +pub(super) struct MaskedObservableDecoder { + inner: Box, + mask: ObsMask, +} + +impl MaskedObservableDecoder { + pub(super) fn new(inner: Box, mask: ObsMask) -> Self { + Self { inner, mask } + } +} + +impl ObservableDecoder for MaskedObservableDecoder { + fn decode_obs(&mut self, syndrome: &[u8]) -> Result { + let mut prediction = self.inner.decode_obs(syndrome)?; + prediction &= &self.mask; + Ok(prediction) + } +} + +/// Observable-decoder adapter that records one elapsed time per decode call. +pub(super) struct TimedObservableDecoder { + inner: Box, + per_shot_seconds: Vec, +} + +impl TimedObservableDecoder { + pub(super) fn new(inner: Box, capacity: usize) -> Self { + Self { + inner, + per_shot_seconds: Vec::with_capacity(capacity), + } + } + + pub(super) fn into_times(self) -> Vec { + self.per_shot_seconds + } +} + +impl ObservableDecoder for TimedObservableDecoder { + fn decode_obs(&mut self, syndrome: &[u8]) -> Result { + let start = std::time::Instant::now(); + let result = self.inner.decode_obs(syndrome); + self.per_shot_seconds.push(start.elapsed().as_secs_f64()); + result + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[derive(Clone)] + enum StubResult { + Prediction(ObsMask), + Error, + } + + struct StubDecoder { + expected_syndromes: Vec>, + results: Vec, + next: usize, + } + + impl StubDecoder { + fn new(expected_syndromes: &[Vec], results: Vec) -> Self { + assert_eq!(expected_syndromes.len(), results.len()); + Self { + expected_syndromes: expected_syndromes.to_vec(), + results, + next: 0, + } + } + } + + impl ObservableDecoder for StubDecoder { + fn decode_obs(&mut self, syndrome: &[u8]) -> Result { + assert_eq!(syndrome, self.expected_syndromes[self.next]); + let result = match &self.results[self.next] { + StubResult::Prediction(mask) => Ok(mask.clone()), + StubResult::Error => Err(DecoderError::DecodingFailed("stub error".into())), + }; + self.next += 1; + result + } + } + + fn mask(bits: &[usize]) -> ObsMask { + let mut mask = ObsMask::new(); + for &bit in bits { + mask.set(bit); + } + mask + } + + fn score( + shots: &[(Vec, ObsMask)], + results: Vec, + ) -> Result { + let syndromes: Vec> = shots.iter().map(|(syndrome, _)| syndrome.clone()).collect(); + let mut decoder = StubDecoder::new(&syndromes, results); + let mut syndrome = vec![0; syndromes.first().map_or(0, Vec::len)]; + count_decoder_mismatches( + 0..shots.len(), + &mut syndrome, + |shot_index, buffer| { + buffer.copy_from_slice(&shots[shot_index].0); + shots[shot_index].1.clone() + }, + &mut decoder, + ) + } + + #[test] + fn decoder_error_aborts_with_shot_index_and_source() { + let shots = vec![ + (vec![0], mask(&[])), + (vec![1], mask(&[0])), + (vec![0], mask(&[])), + ]; + let error = score( + &shots, + vec![ + StubResult::Prediction(mask(&[])), + StubResult::Error, + StubResult::Prediction(mask(&[])), + ], + ) + .unwrap_err(); + + assert_eq!(error.shot_index, 1); + assert!(error.to_string().contains("shot 1")); + assert!(error.to_string().contains("stub error")); + } + + #[test] + fn healthy_decoder_preserves_exact_mismatch_count() { + let shots = vec![ + (vec![0, 0], mask(&[])), + (vec![1, 0], mask(&[0])), + (vec![0, 1], mask(&[1])), + (vec![1, 1], mask(&[0, 1])), + ]; + let count = score( + &shots, + vec![ + StubResult::Prediction(mask(&[])), + StubResult::Prediction(mask(&[])), + StubResult::Prediction(mask(&[1])), + StubResult::Prediction(mask(&[1])), + ], + ) + .unwrap(); + + assert_eq!(count, 2); + } + + #[test] + fn always_erroring_decoder_does_not_match_all_observables_flipped_truth() { + let all_observables = mask(&[0, 1, 2]); + let shots = vec![(vec![1, 1], all_observables)]; + + // The old parallel sampler substituted the full observable-selection + // mask on error, which incorrectly scored this exact truth as correct. + let error = score(&shots, vec![StubResult::Error]).unwrap_err(); + + assert_eq!(error.shot_index, 0); + assert!(error.to_string().contains("stub error")); + } +} diff --git a/python/pecos-rslib/src/fault_tolerance_bindings/sample_corpus.rs b/python/pecos-rslib/src/fault_tolerance_bindings/sample_corpus.rs new file mode 100644 index 000000000..08c7986aa --- /dev/null +++ b/python/pecos-rslib/src/fault_tolerance_bindings/sample_corpus.rs @@ -0,0 +1,848 @@ +// Copyright 2026 The PECOS Developers +// +// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except +// in compliance with the License. You may obtain a copy of the License at +// +// https://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software distributed under the License +// is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express +// or implied. See the License for the specific language governing permissions and limitations under +// the License. + +//! Versioned, self-describing shot-corpus serialization. + +use pecos_decoder_core::dem::SparseDem; +use serde_json::{Map, Value}; +use sha2::{Digest, Sha256}; +use std::path::Path; + +const MAGIC: &[u8; 12] = b"PECOSCORPUS\0"; +pub(super) const FORMAT_VERSION: u32 = 1; +const SHA256_LEN: usize = 32; +const HEADER_LEN_END: usize = MAGIC.len() + size_of::(); +const PREFIX_LEN: usize = HEADER_LEN_END + SHA256_LEN; + +// These limits cover the full verified Jeffreys comparison regime and corpora +// far larger than typical research workloads. Capping each column family at one +// million also bounds degenerate zero-width column vectors to tens of MiB of +// descriptors, rather than allowing attacker-selected multi-gigabyte allocations. +pub(super) const MAX_SHOTS: usize = 100_000_000; +pub(super) const MAX_DETECTORS: usize = 1_000_000; +pub(super) const MAX_OBSERVABLES: usize = 1_000_000; + +#[derive(Debug)] +pub(super) enum CorpusError { + Io(std::io::Error), + Invalid(String), +} + +impl From for CorpusError { + fn from(error: std::io::Error) -> Self { + Self::Io(error) + } +} + +pub(super) struct CorpusToSave<'a> { + pub det_columns: &'a [Vec], + pub obs_columns: &'a [Vec], + pub num_shots: usize, + pub seed: Option, + pub dem: &'a str, + pub metadata_json: Option<&'a str>, +} + +#[derive(Debug, Eq, PartialEq)] +pub(super) struct LoadedCorpus { + pub det_columns: Vec>, + pub obs_columns: Vec>, + pub num_shots: usize, + pub seed: Option, + pub dem: String, + pub metadata_json: Option, + pub generator: String, + pub format_version: u32, +} + +fn invalid(message: impl Into) -> CorpusError { + CorpusError::Invalid(message.into()) +} + +fn sha256(bytes: &[u8]) -> [u8; SHA256_LEN] { + Sha256::digest(bytes).into() +} + +fn hex_encode(bytes: &[u8]) -> String { + const HEX: &[u8; 16] = b"0123456789abcdef"; + let mut output = String::with_capacity(bytes.len() * 2); + for &byte in bytes { + output.push(char::from(HEX[usize::from(byte >> 4)])); + output.push(char::from(HEX[usize::from(byte & 0x0f)])); + } + output +} + +fn sha256_hex(bytes: &[u8]) -> String { + hex_encode(&sha256(bytes)) +} + +fn validate_dimensions( + num_shots: usize, + num_detectors: usize, + num_observables: usize, +) -> Result<(), CorpusError> { + if num_shots > MAX_SHOTS { + return Err(invalid(format!( + "corpus num_shots={num_shots} exceeds the format limit MAX_SHOTS={MAX_SHOTS}" + ))); + } + if num_detectors > MAX_DETECTORS { + return Err(invalid(format!( + "corpus num_detectors={num_detectors} exceeds the format limit MAX_DETECTORS={MAX_DETECTORS}" + ))); + } + if num_observables > MAX_OBSERVABLES { + return Err(invalid(format!( + "corpus num_observables={num_observables} exceeds the format limit MAX_OBSERVABLES={MAX_OBSERVABLES}" + ))); + } + Ok(()) +} + +fn checked_payload_len( + num_detectors: usize, + num_observables: usize, + words_per_column: usize, +) -> Result { + num_detectors + .checked_add(num_observables) + .and_then(|columns| columns.checked_mul(words_per_column)) + .and_then(|words| words.checked_mul(size_of::())) + .ok_or_else(|| invalid("corpus dimensions overflow the supported payload size")) +} + +fn validate_columns(columns: &[Vec], words_per_column: usize) -> Result<(), CorpusError> { + if let Some((index, column)) = columns + .iter() + .enumerate() + .find(|(_, column)| column.len() != words_per_column) + { + return Err(invalid(format!( + "sample column {index} has {} word(s), expected {words_per_column}", + column.len() + ))); + } + Ok(()) +} + +/// Mask selecting the meaningful low bits in a column's final word. +/// +/// The format requires every unused high padding bit to be zero. +fn final_word_mask(num_shots: usize) -> u64 { + let used_bits = num_shots % 64; + if used_bits == 0 { + u64::MAX + } else { + (1_u64 << used_bits) - 1 + } +} + +pub(super) fn save(path: &Path, corpus: CorpusToSave<'_>) -> Result<(), CorpusError> { + validate_dimensions( + corpus.num_shots, + corpus.det_columns.len(), + corpus.obs_columns.len(), + )?; + let parsed_dem = SparseDem::from_dem_str(corpus.dem) + .map_err(|error| invalid(format!("invalid DEM supplied to SampleBatch.save: {error}")))?; + if parsed_dem.num_detectors != corpus.det_columns.len() + || parsed_dem.num_observables != corpus.obs_columns.len() + { + return Err(invalid(format!( + "DEM dimensions do not match SampleBatch: DEM has {} detector(s) and {} observable(s), batch has {} detector(s) and {} observable(s)", + parsed_dem.num_detectors, + parsed_dem.num_observables, + corpus.det_columns.len(), + corpus.obs_columns.len() + ))); + } + + if let Some(metadata) = corpus.metadata_json { + serde_json::from_str::(metadata) + .map_err(|error| invalid(format!("metadata_json is not valid JSON: {error}")))?; + } + + let words_per_column = corpus.num_shots.div_ceil(64); + validate_columns(corpus.det_columns, words_per_column)?; + validate_columns(corpus.obs_columns, words_per_column)?; + let payload_len = checked_payload_len( + corpus.det_columns.len(), + corpus.obs_columns.len(), + words_per_column, + )?; + let mut payload = Vec::with_capacity(payload_len); + for column in corpus.det_columns.iter().chain(corpus.obs_columns) { + for (word_index, &word) in column.iter().enumerate() { + let word = if word_index + 1 == words_per_column { + word & final_word_mask(corpus.num_shots) + } else { + word + }; + payload.extend_from_slice(&word.to_le_bytes()); + } + } + + let header = serde_json::json!({ + "format_version": FORMAT_VERSION, + "num_shots": corpus.num_shots, + "num_detectors": corpus.det_columns.len(), + "num_observables": corpus.obs_columns.len(), + "words_per_column": words_per_column, + "seed": corpus.seed, + "dem": corpus.dem, + "dem_sha256": sha256_hex(corpus.dem.as_bytes()), + "metadata_json": corpus.metadata_json, + "generator": concat!("pecos-rslib ", env!("CARGO_PKG_VERSION")), + }); + let header_bytes = serde_json::to_vec(&header) + .map_err(|error| invalid(format!("could not serialize corpus header: {error}")))?; + let header_len = u32::try_from(header_bytes.len()) + .map_err(|_| invalid("corpus JSON header is too large to encode"))?; + let mut content_hasher = Sha256::new(); + content_hasher.update(&header_bytes); + content_hasher.update(&payload); + let content_sha256: [u8; SHA256_LEN] = content_hasher.finalize().into(); + let file_len = PREFIX_LEN + .checked_add(header_bytes.len()) + .and_then(|len| len.checked_add(payload.len())) + .ok_or_else(|| invalid("corpus file size overflows this platform"))?; + let mut bytes = Vec::with_capacity(file_len); + bytes.extend_from_slice(MAGIC); + bytes.extend_from_slice(&header_len.to_le_bytes()); + bytes.extend_from_slice(&content_sha256); + bytes.extend_from_slice(&header_bytes); + bytes.extend_from_slice(&payload); + std::fs::write(path, bytes)?; + Ok(()) +} + +fn required_u64(header: &Map, field: &str) -> Result { + header.get(field).and_then(Value::as_u64).ok_or_else(|| { + invalid(format!( + "corpus header field {field:?} must be an unsigned integer" + )) + }) +} + +fn required_usize(header: &Map, field: &str) -> Result { + usize::try_from(required_u64(header, field)?).map_err(|_| { + invalid(format!( + "corpus header field {field:?} is too large for this platform" + )) + }) +} + +fn required_string<'a>( + header: &'a Map, + field: &str, +) -> Result<&'a str, CorpusError> { + header + .get(field) + .and_then(Value::as_str) + .ok_or_else(|| invalid(format!("corpus header field {field:?} must be a string"))) +} + +fn nullable_u64(header: &Map, field: &str) -> Result, CorpusError> { + match header.get(field) { + Some(Value::Null) => Ok(None), + Some(value) => value.as_u64().map(Some).ok_or_else(|| { + invalid(format!( + "corpus header field {field:?} must be an unsigned integer or null" + )) + }), + None => Err(invalid(format!( + "corpus header is missing required field {field:?}" + ))), + } +} + +fn nullable_string( + header: &Map, + field: &str, +) -> Result, CorpusError> { + match header.get(field) { + Some(Value::Null) => Ok(None), + Some(Value::String(value)) => Ok(Some(value.clone())), + Some(_) => Err(invalid(format!( + "corpus header field {field:?} must be a string or null" + ))), + None => Err(invalid(format!( + "corpus header is missing required field {field:?}" + ))), + } +} + +pub(super) fn load(path: &Path) -> Result { + let bytes = std::fs::read(path)?; + if bytes.get(..MAGIC.len()) != Some(MAGIC.as_slice()) { + return Err(invalid( + "bad shot-corpus magic: expected PECOSCORPUS followed by a NUL byte", + )); + } + if bytes.len() < HEADER_LEN_END { + return Err(invalid("shot corpus is missing its 4-byte header length")); + } + if bytes.len() < PREFIX_LEN { + return Err(invalid( + "shot corpus is missing its 32-byte content SHA-256", + )); + } + let mut header_len_bytes = [0_u8; size_of::()]; + header_len_bytes.copy_from_slice(&bytes[MAGIC.len()..HEADER_LEN_END]); + let header_len = usize::try_from(u32::from_le_bytes(header_len_bytes)) + .map_err(|_| invalid("corpus header length is too large for this platform"))?; + let header_end = PREFIX_LEN + .checked_add(header_len) + .ok_or_else(|| invalid("corpus header length overflows this platform"))?; + if header_end > bytes.len() { + return Err(invalid(format!( + "corpus header length declares {header_len} byte(s), but the file contains only {} after the prefix", + bytes.len() - PREFIX_LEN + ))); + } + + let expected_content_sha = &bytes[HEADER_LEN_END..PREFIX_LEN]; + let actual_content_sha = sha256(&bytes[PREFIX_LEN..]); + if expected_content_sha != actual_content_sha { + return Err(invalid(format!( + "corpus content SHA-256 mismatch: expected {}, computed {}", + hex_encode(expected_content_sha), + hex_encode(&actual_content_sha) + ))); + } + + // The digest covers the exact header bytes, so a duplicate JSON key cannot + // be injected into an existing corpus without breaking content integrity. + let header_value: Value = serde_json::from_slice(&bytes[PREFIX_LEN..header_end]) + .map_err(|error| invalid(format!("invalid corpus header JSON: {error}")))?; + let header = header_value + .as_object() + .ok_or_else(|| invalid("invalid corpus header JSON: top-level value must be an object"))?; + + let version = required_u64(header, "format_version")?; + if version != u64::from(FORMAT_VERSION) { + return Err(invalid(format!( + "unsupported corpus format_version {version}; this PECOS build supports version {FORMAT_VERSION}" + ))); + } + + let num_shots = required_usize(header, "num_shots")?; + let num_detectors = required_usize(header, "num_detectors")?; + let num_observables = required_usize(header, "num_observables")?; + validate_dimensions(num_shots, num_detectors, num_observables)?; + let words_per_column = required_usize(header, "words_per_column")?; + let expected_words = num_shots.div_ceil(64); + if words_per_column != expected_words { + return Err(invalid(format!( + "corpus words_per_column is {words_per_column}, but num_shots={num_shots} requires {expected_words}" + ))); + } + let seed = nullable_u64(header, "seed")?; + let dem = required_string(header, "dem")?.to_owned(); + let expected_dem_sha = required_string(header, "dem_sha256")?; + let metadata_json = nullable_string(header, "metadata_json")?; + let generator = required_string(header, "generator")?.to_owned(); + + let payload = &bytes[header_end..]; + let expected_payload_len = + checked_payload_len(num_detectors, num_observables, words_per_column)?; + if payload.len() != expected_payload_len { + return Err(invalid(format!( + "corpus payload length is {} byte(s), but declared dimensions require {expected_payload_len} byte(s)", + payload.len() + ))); + } + let actual_dem_sha = sha256_hex(dem.as_bytes()); + if expected_dem_sha != actual_dem_sha { + return Err(invalid(format!( + "corpus DEM SHA-256 mismatch: expected {expected_dem_sha}, computed {actual_dem_sha}" + ))); + } + if let Some(metadata) = &metadata_json { + serde_json::from_str::(metadata) + .map_err(|error| invalid(format!("corpus metadata_json is not valid JSON: {error}")))?; + } + let parsed_dem = SparseDem::from_dem_str(&dem) + .map_err(|error| invalid(format!("corpus contains an invalid DEM: {error}")))?; + if parsed_dem.num_detectors != num_detectors || parsed_dem.num_observables != num_observables { + return Err(invalid(format!( + "corpus DEM dimensions disagree with its header: DEM has {} detector(s) and {} observable(s), header declares {num_detectors} detector(s) and {num_observables} observable(s)", + parsed_dem.num_detectors, parsed_dem.num_observables + ))); + } + + if words_per_column != 0 { + let padding_mask = !final_word_mask(num_shots); + for (column_index, final_word) in payload + .chunks_exact(size_of::()) + .skip(words_per_column - 1) + .step_by(words_per_column) + .enumerate() + { + let mut word_bytes = [0_u8; size_of::()]; + word_bytes.copy_from_slice(final_word); + if u64::from_le_bytes(word_bytes) & padding_mask != 0 { + return Err(invalid(format!( + "corpus payload column {column_index} has nonzero padding bits above num_shots={num_shots}; unused high bits in the final word must be zero" + ))); + } + } + } + + let mut words = Vec::with_capacity(payload.len() / size_of::()); + for chunk in payload.chunks_exact(size_of::()) { + let mut bytes = [0_u8; size_of::()]; + bytes.copy_from_slice(chunk); + words.push(u64::from_le_bytes(bytes)); + } + let mut offset = 0; + let mut read_columns = |count: usize| { + let mut columns = Vec::with_capacity(count); + for _ in 0..count { + let end = offset + words_per_column; + columns.push(words[offset..end].to_vec()); + offset = end; + } + columns + }; + let det_columns = read_columns(num_detectors); + let obs_columns = read_columns(num_observables); + + Ok(LoadedCorpus { + det_columns, + obs_columns, + num_shots, + seed, + dem, + metadata_json, + generator, + format_version: FORMAT_VERSION, + }) +} + +#[cfg(test)] +mod tests { + use super::*; + use tempfile::TempDir; + + const DEM: &str = "error(0.125) D0 L0\n"; + + fn corpus_path() -> (TempDir, std::path::PathBuf) { + let directory = tempfile::tempdir().unwrap(); + let path = directory.path().join("shots.pecos"); + (directory, path) + } + + fn save_test_corpus(path: &Path) { + save( + path, + CorpusToSave { + det_columns: &[vec![0b10]], + obs_columns: &[vec![0b10]], + num_shots: 2, + seed: Some(42), + dem: DEM, + metadata_json: Some(r#"{ "decoder": "pymatching" }"#), + }, + ) + .unwrap(); + } + + fn invalid_message(result: Result) -> String { + match result.unwrap_err() { + CorpusError::Invalid(message) => message, + CorpusError::Io(error) => panic!("unexpected I/O error: {error}"), + } + } + + fn header_end(bytes: &[u8]) -> usize { + let mut length = [0_u8; 4]; + length.copy_from_slice(&bytes[MAGIC.len()..HEADER_LEN_END]); + PREFIX_LEN + usize::try_from(u32::from_le_bytes(length)).unwrap() + } + + fn with_valid_content_sha(mut bytes: Vec) -> Vec { + let digest = sha256(&bytes[PREFIX_LEN..]); + bytes[HEADER_LEN_END..PREFIX_LEN].copy_from_slice(&digest); + bytes + } + + fn replace_header(bytes: &[u8], update: impl FnOnce(&mut Map)) -> Vec { + let old_header_end = header_end(bytes); + let mut header: Value = serde_json::from_slice(&bytes[PREFIX_LEN..old_header_end]).unwrap(); + update(header.as_object_mut().unwrap()); + let new_header = serde_json::to_vec(&header).unwrap(); + let new_header_len = u32::try_from(new_header.len()).unwrap(); + let mut updated = Vec::new(); + updated.extend_from_slice(MAGIC); + updated.extend_from_slice(&new_header_len.to_le_bytes()); + updated.extend_from_slice(&bytes[HEADER_LEN_END..PREFIX_LEN]); + updated.extend_from_slice(&new_header); + updated.extend_from_slice(&bytes[old_header_end..]); + updated + } + + #[test] + fn round_trip_preserves_columns_dimensions_and_provenance() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + + let loaded = load(&path).unwrap(); + assert_eq!(loaded.det_columns, vec![vec![0b10]]); + assert_eq!(loaded.obs_columns, vec![vec![0b10]]); + assert_eq!(loaded.num_shots, 2); + assert_eq!(loaded.seed, Some(42)); + assert_eq!(loaded.dem, DEM); + assert_eq!( + loaded.metadata_json.as_deref(), + Some(r#"{ "decoder": "pymatching" }"#) + ); + assert!(loaded.generator.starts_with("pecos-rslib ")); + assert_eq!(loaded.format_version, FORMAT_VERSION); + + let bytes = std::fs::read(path).unwrap(); + let header: Value = serde_json::from_slice(&bytes[PREFIX_LEN..header_end(&bytes)]).unwrap(); + assert!(header.get("payload_sha256").is_none()); + } + + #[test] + fn wide_observable_column_round_trips_without_narrowing() { + let (_directory, path) = corpus_path(); + let mut observables = vec![vec![0]; 65]; + observables[64][0] = 1; + save( + &path, + CorpusToSave { + det_columns: &[vec![1]], + obs_columns: &observables, + num_shots: 1, + seed: None, + dem: "error(0.125) D0 L64\n", + metadata_json: None, + }, + ) + .unwrap(); + + let loaded = load(&path).unwrap(); + assert_eq!(loaded.obs_columns.len(), 65); + assert_eq!(loaded.obs_columns[64], vec![1]); + assert_eq!(loaded.seed, None); + } + + #[test] + fn corrupted_payload_fails_content_checksum() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let mut bytes = std::fs::read(&path).unwrap(); + let last = bytes.last_mut().unwrap(); + *last ^= 0x80; + std::fs::write(&path, bytes).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("content SHA-256 mismatch"), "{message}"); + } + + #[test] + fn truncated_payload_fails_content_checksum() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let mut bytes = std::fs::read(&path).unwrap(); + bytes.pop(); + std::fs::write(&path, bytes).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("content SHA-256 mismatch"), "{message}"); + } + + #[test] + fn authenticated_truncated_payload_fails_length_check() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let mut bytes = std::fs::read(&path).unwrap(); + bytes.pop(); + let bytes = with_valid_content_sha(bytes); + std::fs::write(&path, bytes).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("payload length"), "{message}"); + } + + #[test] + fn changing_num_shots_and_header_len_breaks_content_checksum() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = replace_header(&bytes, |header| { + header.insert("num_shots".to_owned(), Value::from(3)); + }); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("content SHA-256 mismatch"), "{message}"); + } + + #[test] + fn changing_seed_and_header_len_breaks_content_checksum() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = replace_header(&bytes, |header| { + header.insert("seed".to_owned(), Value::from(43)); + }); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("content SHA-256 mismatch"), "{message}"); + } + + #[test] + fn bad_magic_is_rejected_first() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let mut bytes = std::fs::read(&path).unwrap(); + bytes[0] ^= 1; + std::fs::write(&path, bytes).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("bad shot-corpus magic"), "{message}"); + } + + #[test] + fn future_format_version_is_actionable() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("format_version".to_owned(), Value::from(999)); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!( + message.contains("unsupported corpus format_version 999"), + "{message}" + ); + } + + #[test] + fn invalid_header_json_is_rejected_without_panicking() { + let (_directory, path) = corpus_path(); + let mut bytes = Vec::from(MAGIC.as_slice()); + bytes.extend_from_slice(&1_u32.to_le_bytes()); + bytes.extend_from_slice(&[0_u8; SHA256_LEN]); + bytes.push(b'{'); + let bytes = with_valid_content_sha(bytes); + std::fs::write(&path, bytes).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("invalid corpus header JSON"), "{message}"); + } + + #[test] + fn dem_checksum_is_verified_after_content_checksum() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("dem".to_owned(), Value::from("error(0.25) D0 L0\n")); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("DEM SHA-256 mismatch"), "{message}"); + } + + #[test] + fn loaded_dem_dimensions_must_match_header() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let replacement_dem = "error(0.25) D1 L0\n"; + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("dem".to_owned(), Value::from(replacement_dem)); + header.insert( + "dem_sha256".to_owned(), + Value::from(sha256_hex(replacement_dem.as_bytes())), + ); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!( + message.contains("corpus DEM dimensions disagree with its header"), + "{message}" + ); + } + + #[test] + fn loaded_metadata_json_must_be_valid() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("metadata_json".to_owned(), Value::from("{")); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!( + message.contains("corpus metadata_json is not valid JSON"), + "{message}" + ); + } + + #[test] + fn declared_shots_above_limit_are_rejected_with_zero_columns() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("num_shots".to_owned(), Value::from(MAX_SHOTS + 1)); + header.insert("num_detectors".to_owned(), Value::from(0)); + header.insert("num_observables".to_owned(), Value::from(0)); + header.insert( + "words_per_column".to_owned(), + Value::from((MAX_SHOTS + 1).div_ceil(64)), + ); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!( + message.contains(&format!("MAX_SHOTS={MAX_SHOTS}")), + "{message}" + ); + } + + #[test] + fn declared_detectors_above_limit_are_rejected_with_zero_shots() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("num_shots".to_owned(), Value::from(0)); + header.insert("num_detectors".to_owned(), Value::from(MAX_DETECTORS + 1)); + header.insert("num_observables".to_owned(), Value::from(0)); + header.insert("words_per_column".to_owned(), Value::from(0)); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!( + message.contains(&format!("MAX_DETECTORS={MAX_DETECTORS}")), + "{message}" + ); + } + + #[test] + fn declared_observables_above_limit_are_rejected_with_zero_shots() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let bytes = std::fs::read(&path).unwrap(); + let updated = with_valid_content_sha(replace_header(&bytes, |header| { + header.insert("num_shots".to_owned(), Value::from(0)); + header.insert("num_detectors".to_owned(), Value::from(0)); + header.insert( + "num_observables".to_owned(), + Value::from(MAX_OBSERVABLES + 1), + ); + header.insert("words_per_column".to_owned(), Value::from(0)); + })); + std::fs::write(&path, updated).unwrap(); + + let message = invalid_message(load(&path)); + assert!( + message.contains(&format!("MAX_OBSERVABLES={MAX_OBSERVABLES}")), + "{message}" + ); + } + + #[test] + fn save_masks_unused_high_padding_bits() { + let (_directory, path) = corpus_path(); + save( + &path, + CorpusToSave { + det_columns: &[vec![u64::MAX]], + obs_columns: &[vec![u64::MAX]], + num_shots: 1, + seed: None, + dem: DEM, + metadata_json: None, + }, + ) + .unwrap(); + + let loaded = load(&path).unwrap(); + assert_eq!(loaded.det_columns, vec![vec![1]]); + assert_eq!(loaded.obs_columns, vec![vec![1]]); + } + + #[test] + fn nonzero_payload_padding_is_rejected() { + let (_directory, path) = corpus_path(); + save_test_corpus(&path); + let mut bytes = std::fs::read(&path).unwrap(); + let payload_start = header_end(&bytes); + bytes[payload_start + size_of::() - 1] |= 0x80; + let bytes = with_valid_content_sha(bytes); + std::fs::write(&path, bytes).unwrap(); + + let message = invalid_message(load(&path)); + assert!(message.contains("nonzero padding bits"), "{message}"); + } + + #[test] + fn mismatched_dem_dimensions_are_rejected_before_writing() { + let (_directory, path) = corpus_path(); + let result = save( + &path, + CorpusToSave { + det_columns: &[vec![0]], + obs_columns: &[vec![0]], + num_shots: 1, + seed: None, + dem: "error(0.125) D1 L0\n", + metadata_json: None, + }, + ); + + let CorpusError::Invalid(message) = result.unwrap_err() else { + panic!("expected malformed-input error"); + }; + assert!(message.contains("DEM dimensions do not match SampleBatch")); + assert!(!path.exists()); + } + + #[test] + fn invalid_metadata_json_is_rejected_before_writing() { + let (_directory, path) = corpus_path(); + let result = save( + &path, + CorpusToSave { + det_columns: &[vec![0]], + obs_columns: &[vec![0]], + num_shots: 1, + seed: None, + dem: DEM, + metadata_json: Some("{"), + }, + ); + + let CorpusError::Invalid(message) = result.unwrap_err() else { + panic!("expected malformed-input error"); + }; + assert!(message.contains("metadata_json is not valid JSON")); + assert!(!path.exists()); + } +} diff --git a/python/quantum-pecos/tests/qec/test_decoder_comparison.py b/python/quantum-pecos/tests/qec/test_decoder_comparison.py new file mode 100644 index 000000000..82b180f8c --- /dev/null +++ b/python/quantum-pecos/tests/qec/test_decoder_comparison.py @@ -0,0 +1,52 @@ +# Copyright 2026 The PECOS Developers +# Licensed under the Apache License, Version 2.0 + +"""Python coverage for paired DUT/reference decoder comparison.""" + +from __future__ import annotations + +import math + +import pytest + +pytest.importorskip("pecos_rslib") + +from pecos_rslib.qec import SampleBatch + + +def test_sample_batch_compare_decoders_exposes_joint_counts() -> None: + dem = "error(0.1) D0 L0\n" + batch = SampleBatch([[0], [1], [0], [1]], [0, 1, 0, 1]) + + first = batch.compare_decoders(dem, "pymatching", "pymatching") + second = batch.compare_decoders(dem, "pymatching", "pymatching") + + assert first.total_shots == 4 + assert first.counts == [[4, 0, 0], [0, 0, 0], [0, 0, 0]] + assert first.dut_correct_reference_correct == 4 + assert first.dut_only_failures == 0 + assert first.both_failed == 0 + assert first.dut_only_failure_interval[0] >= 0.0 + assert first.dut_only_failure_interval[1] <= 1.0 + assert second.counts == first.counts + + +def test_compare_decoders_rejects_empty_batch_before_decoder_construction() -> None: + batch = SampleBatch([], []) + + with pytest.raises(ValueError, match="n must be greater than zero"): + batch.compare_decoders("not a DEM", "not a decoder", "not a decoder") + + +@pytest.mark.parametrize("alpha", [math.nan, -1.0, 0.0, 1.0, 2.0, math.inf]) +def test_compare_decoders_rejects_alpha_outside_open_unit_interval(alpha: float) -> None: + dem = "error(0.1) D0 L0\n" + batch = SampleBatch([[0]], [0]) + + with pytest.raises(ValueError, match=r"alpha must be finite and in \(0, 1\)"): + batch.compare_decoders( + dem, + "pymatching", + "pymatching", + alpha=alpha, + ) diff --git a/python/quantum-pecos/tests/qec/test_frontier_decoder.py b/python/quantum-pecos/tests/qec/test_frontier_decoder.py new file mode 100644 index 000000000..a60947d07 --- /dev/null +++ b/python/quantum-pecos/tests/qec/test_frontier_decoder.py @@ -0,0 +1,108 @@ +# Copyright 2026 The PECOS Developers +# +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use +# this file except in compliance with the License. You may obtain a copy of the +# License at +# +# https://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. See the License for the +# specific language governing permissions and limitations under the License. + +"""Tests for the experimental native Frontier decoder bindings.""" + +from __future__ import annotations + +import pytest + +pecos_rslib_exp = pytest.importorskip("pecos_rslib_exp") + +from pecos_rslib_exp import ( # noqa: E402 + FrontierCommitteeDecoder, + FrontierDecoder, +) + +SMALL_DEM = """\ +error(0.1) D0 L0 +error(0.2) D1 +""" + + +def test_sparse_and_dense_decode_agree() -> None: + decoder = FrontierDecoder.from_dem(SMALL_DEM) + + # Syndrome D0=1, D1=0 forces the first mechanism on and the second off, + # so logical observable L0 flips and the expected mask is 1. + dense = decoder.decode_syndrome([1, 0]) + sparse = decoder.decode([0]) + + assert dense.observables_mask == 1 + assert sparse.observables_mask == dense.observables_mask + assert sparse.log_evidence == dense.log_evidence + assert sparse.logical_masses == dense.logical_masses + assert dense.observable_bits(2) == [1, 0] + + +def test_wide_observable_mask_is_a_python_int() -> None: + decoder = FrontierDecoder.from_dem("error(0.9) D0 L64\n") + + result = decoder.decode_syndrome([1]) + + assert isinstance(result.observables_mask, int) + assert result.observables_mask == 1 << 64 + assert result.observable_bits(65)[64] == 1 + assert result.logical_masses[0][0] == 1 << 64 + + +def test_column_order_variants_and_validation() -> None: + expected = 1 + for column_order in ( + "deadline_reorder", + "time_order", + "backward_deadline_reorder", + [1, 0], + ): + decoder = FrontierDecoder.from_dem(SMALL_DEM, column_order=column_order) + assert decoder.decode_syndrome([1, 0]).observables_mask == expected + + with pytest.raises(ValueError, match="invalid column_order"): + FrontierDecoder.from_dem(SMALL_DEM, column_order="not_an_order") + with pytest.raises(ValueError, match="column_order must be"): + FrontierDecoder.from_dem(SMALL_DEM, column_order=object()) + with pytest.raises(RuntimeError, match="permutation"): + FrontierDecoder.from_dem(SMALL_DEM, column_order=[0, 0]) + + +def test_committee_easy_tie_selects_forward() -> None: + committee = FrontierCommitteeDecoder.from_dem("error(0.1) D0 L0\n", column_order="time_order") + + result = committee.decode([0]) + + assert result.observables_mask == 1 + assert result.direction == "forward" + assert result.forward_log_evidence == result.log_evidence + assert result.backward_log_evidence == result.log_evidence + + +def test_unexplainable_and_out_of_range_syndromes_raise_runtime_error() -> None: + decoder = FrontierDecoder.from_dem("error(0.1) D0 D1\n") + + with pytest.raises(RuntimeError, match="unexplainable"): + decoder.decode_syndrome([1, 0]) + with pytest.raises(RuntimeError, match="Invalid node index"): + decoder.decode([2]) + + +def test_decode_batch_matches_individual_decodes() -> None: + shots = [[0, 0], [1, 0], [0, 1], [1, 1]] + batch_decoder = FrontierDecoder.from_dem(SMALL_DEM) + individual_decoder = FrontierDecoder.from_dem(SMALL_DEM) + + batch = batch_decoder.decode_batch(shots) + individual = [individual_decoder.decode_syndrome(shot) for shot in shots] + + assert [result.observables_mask for result in batch] == [result.observables_mask for result in individual] + assert [result.log_evidence for result in batch] == [result.log_evidence for result in individual] + assert [result.logical_masses for result in batch] == [result.logical_masses for result in individual] diff --git a/python/quantum-pecos/tests/qec/test_sample_batch.py b/python/quantum-pecos/tests/qec/test_sample_batch.py index 4170da0ac..ac896e425 100644 --- a/python/quantum-pecos/tests/qec/test_sample_batch.py +++ b/python/quantum-pecos/tests/qec/test_sample_batch.py @@ -93,3 +93,18 @@ def test_decode_count(self, d3_setup): errors = batch.decode_count(dem_str, "pymatching") assert isinstance(errors, int) assert 0 <= errors <= 1000 + + +def test_seeded_healthy_decoder_counts_and_stats_are_unchanged() -> None: + dem = "error(0.1) D0 L0\nerror(0.2) D0\n" + batch = DemSampler.from_dem_string(dem).generate_samples(257, seed=314159) + + assert batch.decode_count(dem, "pymatching") == 48 + assert batch.decode_count_batch(dem) == 48 + assert batch.decode_count_parallel(dem, "pymatching", num_workers=3) == 48 + + stats = batch.decode_stats(dem, "pymatching") + parallel_stats = batch.decode_stats_parallel(dem, "pymatching", num_workers=3) + assert stats.num_shots == parallel_stats.num_shots == 257 + assert stats.num_errors == parallel_stats.num_errors == 48 + assert stats.logical_error_rate == parallel_stats.logical_error_rate == 48 / 257 diff --git a/python/quantum-pecos/tests/qec/test_sample_corpus.py b/python/quantum-pecos/tests/qec/test_sample_corpus.py new file mode 100644 index 000000000..32af4c138 --- /dev/null +++ b/python/quantum-pecos/tests/qec/test_sample_corpus.py @@ -0,0 +1,161 @@ +# Copyright 2026 The PECOS Developers +# Licensed under the Apache License, Version 2.0 + +"""Python coverage for self-describing SampleBatch shot corpora.""" + +from __future__ import annotations + +import errno + +import pytest + +pytest.importorskip("pecos_rslib") + +from pecos_rslib.qec import DemSampler, SampleBatch + + +def test_generated_batch_round_trip_preserves_shots_and_provenance(tmp_path) -> None: + dem = "error(0.125) D0 L0\n" + metadata = '{ "decoder": "pymatching", "decoder_seed": 17 }' + batch = DemSampler.from_dem_string(dem).generate_samples(130, seed=42) + path = tmp_path / "round-trip.pecos" + + batch.save(path, dem=dem, metadata_json=metadata) + loaded = SampleBatch.load(path) + + assert loaded.num_shots == batch.num_shots + assert loaded.seed == batch.seed == 42 + assert loaded.dem == dem + assert loaded.metadata_json == metadata + assert loaded.generator.startswith("pecos-rslib ") + assert loaded.format_version == 1 + assert [loaded.get_syndrome(i) for i in range(loaded.num_shots)] == [ + batch.get_syndrome(i) for i in range(batch.num_shots) + ] + assert [loaded.get_observable_mask_wide(i) for i in range(loaded.num_shots)] == [ + batch.get_observable_mask_wide(i) for i in range(batch.num_shots) + ] + + +def test_wide_observable_above_bit_63_round_trips(tmp_path) -> None: + dem = "error(0.125) D0 L64\n" + batch = SampleBatch([[1], [0]], [1 << 64, 0]) + path = tmp_path / "wide.pecos" + + batch.save(path, dem=dem) + loaded = SampleBatch.load(path) + + assert loaded.get_observable_mask_wide(0) == 1 << 64 + assert loaded.get_observable_mask_wide(1) == 0 + + +def test_save_rejects_mismatched_dem_dimensions(tmp_path) -> None: + batch = SampleBatch([[0]], [0]) + + with pytest.raises(ValueError, match="DEM dimensions do not match SampleBatch"): + batch.save(tmp_path / "wrong-dem.pecos", dem="error(0.1) D1\n") + + +def test_save_rejects_invalid_metadata_json(tmp_path) -> None: + batch = SampleBatch([[0]], [0]) + + with pytest.raises(ValueError, match="metadata_json is not valid JSON"): + batch.save( + tmp_path / "bad-metadata.pecos", + dem="error(0.1) D0\n", + metadata_json="{", + ) + + +def test_load_maps_malformed_files_to_value_error(tmp_path) -> None: + path = tmp_path / "bad-magic.pecos" + path.write_bytes(b"not a PECOS corpus") + + with pytest.raises(ValueError, match="bad shot-corpus magic"): + SampleBatch.load(path) + + +def test_load_maps_filesystem_failures_to_io_error(tmp_path) -> None: + path = tmp_path / "missing.pecos" + + with pytest.raises(FileNotFoundError, match="No such file or directory") as exc_info: + SampleBatch.load(path) + + assert exc_info.value.errno == errno.ENOENT + assert exc_info.value.filename == str(path) + + +def test_resave_preserves_metadata_unless_explicitly_cleared(tmp_path) -> None: + dem = "error(0.125) D0 L0\n" + metadata = '{"source": "original"}' + original_path = tmp_path / "original.pecos" + preserved_path = tmp_path / "preserved.pecos" + cleared_path = tmp_path / "cleared.pecos" + batch = SampleBatch([[1]], [1]) + batch.save(original_path, dem=dem, metadata_json=metadata) + loaded = SampleBatch.load(original_path) + + loaded.save(preserved_path, dem=dem) + loaded.save(cleared_path, dem=dem, clear_metadata=True) + + assert SampleBatch.load(preserved_path).metadata_json == metadata + assert SampleBatch.load(cleared_path).metadata_json is None + + +def test_loaded_batch_requires_its_embedded_dem_unless_opted_out(tmp_path) -> None: + embedded_dem = "error(0.125) D0 L0\n" + different_dem = "error(0.25) D0 L0\n" + path = tmp_path / "dem-bound.pecos" + batch = SampleBatch([[1]], [1]) + batch.save(path, dem=embedded_dem) + loaded = SampleBatch.load(path) + + with pytest.raises(ValueError, match="differs from the DEM embedded"): + loaded.compare_decoders( + different_dem, + "pymatching", + "pymatching", + ) + + result = loaded.compare_decoders( + different_dem, + "pymatching", + "pymatching", + allow_dem_mismatch=True, + ) + assert result.total_shots == 1 + + +def test_generate_samples_records_resolved_and_explicit_seeds() -> None: + sampler = DemSampler.from_dem_string("error(0.125) D0 L0\n") + + resolved = sampler.generate_samples(130) + explicit = sampler.generate_samples(1, seed=0xDEADBEEF) + replayed = sampler.generate_samples(resolved.num_shots, seed=resolved.seed) + + assert isinstance(resolved.seed, int) + assert 0 <= resolved.seed <= (1 << 64) - 1 + assert explicit.seed == 0xDEADBEEF + assert [replayed.get_syndrome(i) for i in range(replayed.num_shots)] == [ + resolved.get_syndrome(i) for i in range(resolved.num_shots) + ] + assert [replayed.get_observable_mask_wide(i) for i in range(replayed.num_shots)] == [ + resolved.get_observable_mask_wide(i) for i in range(resolved.num_shots) + ] + + +def test_compare_decoders_counts_survive_corpus_round_trip(tmp_path) -> None: + dem = "error(0.1) D0 L0\nerror(0.1) D0\n" + batch = DemSampler.from_dem_string(dem).generate_samples(257, seed=314159) + before = batch.compare_decoders(dem, "pymatching", "pymatching") + path = tmp_path / "comparison.pecos" + + batch.save(path, dem=dem) + loaded = SampleBatch.load(path) + after = loaded.compare_decoders( + loaded.dem, + "pymatching", + "pymatching", + ) + + assert after.counts == before.counts diff --git a/ruff.toml b/ruff.toml index 909707189..18d8ca2e5 100644 --- a/ruff.toml +++ b/ruff.toml @@ -108,6 +108,10 @@ ignore = [ "INP001", # Implicit namespace package - standalone scripts, not a package "S311", # Standard pseudo-random - fixed-seed case sampling, not crypto ] +"exp/*/tests/fixtures/*.py" = [ + "INP001", # Implicit namespace package - standalone scripts, not a package + "S311", # Standard pseudo-random - fixed-seed case sampling, not crypto +] # Main pecos __init__.py - special case for module initialization "python/quantum-pecos/src/pecos/__init__.py" = [