diff --git a/libs/qec/include/cudaq/qec/code.h b/libs/qec/include/cudaq/qec/code.h index f7b9ea1d..e41177a7 100644 --- a/libs/qec/include/cudaq/qec/code.h +++ b/libs/qec/include/cudaq/qec/code.h @@ -201,6 +201,47 @@ class code : public cudaqx::extension_point { /// @return Tensor representing Lz cudaqx::tensor get_observables_z() const; + /// @brief Get the inlined feedback matrix for detector construction + /// + /// Shape is [numCols x numCols], where numCols = + /// get_num_ancilla_qubits() (one measurement record per ancilla per + /// round, in [Z][X] order). Entry (j, k) = 1 means the cross-round + /// detector comparing record j between consecutive rounds additionally + /// XORs record k of the earlier round, and the final boundary detector + /// for record j additionally XORs record k of the last round. + /// @return Tensor representing the inlined feedback matrix. An empty + /// tensor (the default) means no feedback is applied and detectors are + /// formed from same-record comparisons only. + virtual cudaqx::tensor get_inlined_feedback() const; + + /// @brief Get the inlined feedback matrix for logical observables measured + /// in the Z basis + /// + /// Shape is [num_observables x numCols], where numCols = + /// get_num_ancilla_qubits() (one measurement record per ancilla per + /// round, in [Z][X] order). Entry (m, k) = 1 means logical observable m + /// additionally XORs record k of every round. This getter is consumed when + /// the memory experiment measures its observables in the Z basis (prep0 / + /// prep1 state preparations). + /// @return Tensor representing the Z-basis observable inlined feedback + /// matrix. An empty tensor (the default) means no observable feedback is + /// applied in the Z basis. + virtual cudaqx::tensor get_observable_inlined_feedback_z() const; + + /// @brief Get the inlined feedback matrix for logical observables measured + /// in the X basis + /// + /// Shape is [num_observables x numCols], where numCols = + /// get_num_ancilla_qubits() (one measurement record per ancilla per + /// round, in [Z][X] order). Entry (m, k) = 1 means logical observable m + /// additionally XORs record k of every round. This getter is consumed when + /// the memory experiment measures its observables in the X basis (prepp / + /// prepm state preparations). + /// @return Tensor representing the X-basis observable inlined feedback + /// matrix. An empty tensor (the default) means no observable feedback is + /// applied in the X basis. + virtual cudaqx::tensor get_observable_inlined_feedback_x() const; + /// @brief Get the stabilizer generators /// @return Reference to stabilizers const std::vector &get_stabilizers() const { diff --git a/libs/qec/include/cudaq/qec/detector_error_model.h b/libs/qec/include/cudaq/qec/detector_error_model.h index 64729f3d..88b05330 100644 --- a/libs/qec/include/cudaq/qec/detector_error_model.h +++ b/libs/qec/include/cudaq/qec/detector_error_model.h @@ -82,6 +82,49 @@ struct detector_error_model { bool remove_zero_syndrome_errors = false); }; +/// @brief Record-to-detector composition layout derived from a code's +/// declared inlined-feedback matrices (CSR-style). +/// +/// Entry range [detector_offsets[j], detector_offsets[j+1]) of +/// detector_indices lists the herald record columns XOR-ed into record j's +/// cross-round detector (records of the earlier round) and into its final +/// boundary detector (records of the last round). Entry range +/// [observable_offsets[m], observable_offsets[m+1]) of observable_indices +/// lists the record columns XOR-ed into logical observable m on every round. +/// Empty offsets vectors mean no feedback is declared for that target +/// (identity detector/observable structure). +struct inlined_feedback_layout { + std::vector detector_indices; + std::vector detector_offsets; + std::vector observable_indices; + std::vector observable_offsets; +}; + +/// @brief Build the record-to-detector composition layout from a code's +/// declared inlined-feedback matrices. +/// +/// @param feedback Detector feedback matrix, shape +/// [num_syndromes_per_round x num_syndromes_per_round] with 0/1 +/// entries, or an empty tensor for none. Entry (j, k) = 1 means the +/// cross-round detector for record j additionally XORs record k of +/// the earlier round, and the final boundary detector for record j +/// additionally XORs record k of the last round. +/// @param obs_feedback Observable feedback matrix, shape +/// [num_observables x num_syndromes_per_round] with 0/1 entries, or +/// an empty tensor for none. Entry (m, k) = 1 means logical +/// observable m additionally XORs record k of every round. +/// @param num_syndromes_per_round Records per round (number of ancilla +/// qubits). +/// @param num_observables Number of logical observables. +/// @return The CSR layout; see inlined_feedback_layout. +/// @throws std::runtime_error if a non-empty tensor's shape does not match +/// the expected dimensions. +inlined_feedback_layout +build_inlined_feedback_layout(const cudaqx::tensor &feedback, + const cudaqx::tensor &obs_feedback, + std::size_t num_syndromes_per_round, + std::size_t num_observables); + /// Parse the Stim DEM string @p dem_text into detector/observable flip /// matrices and error rates. DEM-native decoders should consume raw DEM text /// instead. By default (@p use_decomp_suggestions = false) the '^' separators diff --git a/libs/qec/lib/code.cpp b/libs/qec/lib/code.cpp index 19693b8c..39be5efe 100644 --- a/libs/qec/lib/code.cpp +++ b/libs/qec/lib/code.cpp @@ -61,6 +61,18 @@ cudaqx::tensor code::get_observables_z() const { return to_parity_matrix(m_pauli_observables, stabilizer_type::Z); } +cudaqx::tensor code::get_inlined_feedback() const { + return cudaqx::tensor(); +} + +cudaqx::tensor code::get_observable_inlined_feedback_z() const { + return cudaqx::tensor(); +} + +cudaqx::tensor code::get_observable_inlined_feedback_x() const { + return cudaqx::tensor(); +} + std::unique_ptr get_code(const std::string &name, const std::vector &stab, const heterogeneous_map options) { diff --git a/libs/qec/lib/detector_error_model.cpp b/libs/qec/lib/detector_error_model.cpp index 79a89b79..338f455c 100644 --- a/libs/qec/lib/detector_error_model.cpp +++ b/libs/qec/lib/detector_error_model.cpp @@ -296,4 +296,49 @@ void detector_error_model::canonicalize_for_rounds( this->observables_flips_matrix, final_column_order); } +/// @brief Convert one 0/1 feedback matrix into a CSR row layout. An empty +/// tensor (the identity default) leaves @p indices and @p offsets empty; a +/// non-empty tensor must have shape [num_rows, num_cols]. +static void feedback_rows_to_csr(const cudaqx::tensor &matrix, + std::size_t num_rows, std::size_t num_cols, + const char *name, + std::vector &indices, + std::vector &offsets) { + if (matrix.rank() == 0 || matrix.size() == 0) + return; + if (matrix.rank() != 2 || matrix.shape()[0] != num_rows || + matrix.shape()[1] != num_cols) { + std::string actual; + for (std::size_t d = 0; d < matrix.rank(); ++d) + actual += (d ? ", " : "") + std::to_string(matrix.shape()[d]); + throw std::runtime_error( + std::string(name) + " has invalid shape [" + actual + "] - expected [" + + std::to_string(num_rows) + ", " + std::to_string(num_cols) + + "] or an empty tensor."); + } + offsets.resize(num_rows + 1); + offsets[0] = 0; + for (std::size_t r = 0; r < num_rows; ++r) { + for (std::size_t c = 0; c < num_cols; ++c) + if (matrix.at({r, c}) != 0) + indices.push_back(c); + offsets[r + 1] = indices.size(); + } +} + +inlined_feedback_layout +build_inlined_feedback_layout(const cudaqx::tensor &feedback, + const cudaqx::tensor &obs_feedback, + std::size_t num_syndromes_per_round, + std::size_t num_observables) { + inlined_feedback_layout layout; + feedback_rows_to_csr(feedback, num_syndromes_per_round, + num_syndromes_per_round, "inlined feedback", + layout.detector_indices, layout.detector_offsets); + feedback_rows_to_csr(obs_feedback, num_observables, num_syndromes_per_round, + "observable inlined feedback", layout.observable_indices, + layout.observable_offsets); + return layout; +} + } // namespace cudaq::qec diff --git a/libs/qec/lib/device/inlined_feedback.h b/libs/qec/lib/device/inlined_feedback.h new file mode 100644 index 00000000..9884ba45 --- /dev/null +++ b/libs/qec/lib/device/inlined_feedback.h @@ -0,0 +1,161 @@ +/****************************************************************-*- C++ -*-**** + * Copyright (c) 2024 - 2025 NVIDIA Corporation & Affiliates. * + * All rights reserved. * + * * + * This source code and the accompanying materials are made available under * + * the terms of the Apache License 2.0 which accompanies this distribution. * + ******************************************************************************/ + +#pragma once + +#include "cudaq.h" + +namespace cudaq::qec { + +/// @brief Number of nonzero entries in row @p row of a row-major 0/1 matrix +/// that has @p num_cols columns (i.e. the weight of one feedback/observable/ +/// stabilizer row). +inline __qpu__ std::size_t +row_support_weight(const std::vector &matrix_flat, std::size_t row, + std::size_t num_cols) { + std::size_t weight = 0; + for (std::size_t k = 0; k < num_cols; ++k) + if (matrix_flat[row * num_cols + k] != 0) + weight++; + return weight; +} + +/// @brief Records for the cross-round detector of syndrome record @p j: the +/// earlier-vs-current comparison `{prev[j], curr[j]}` followed by the +/// earlier-round herald records `prev[k]` for each column k (ascending) with +/// `feedback_flat(j, k) != 0`. +inline __qpu__ std::vector cross_round_detector_records( + const std::vector &prev, + const std::vector &curr, std::size_t j, + const std::vector &feedback_flat, std::size_t num_cols) { + std::size_t weight = row_support_weight(feedback_flat, j, num_cols); + std::vector det(2 + weight); + det[0] = prev[j]; + det[1] = curr[j]; + std::size_t idx = 2; + for (std::size_t k = 0; k < num_cols; ++k) + if (feedback_flat[j * num_cols + k] != 0) + det[idx++] = prev[k]; + return det; +} + +/// @brief Records for the boundary detector of syndrome record @p record_row: +/// `{last_syndrome[record_row]}`, then the data-qubit readouts +/// `data_results[q]` for each data qubit q (ascending) in the stabilizer +/// support `stabilizers[row_base + q]`, then the last-round herald records +/// `last_syndrome[k]` for each column k (ascending) with +/// `feedback_flat(record_row, k) != 0`. When @p feedback_flat is empty the +/// herald part contributes nothing, so one helper serves both the +/// feedback and legacy branches. +inline __qpu__ std::vector boundary_detector_records( + const std::vector &last_syndrome, + std::size_t record_row, + const std::vector &data_results, + const std::vector &stabilizers, std::size_t row_base, + std::size_t num_data, const std::vector &feedback_flat, + std::size_t num_cols) { + std::size_t support_weight = 0; + for (std::size_t q = 0; q < num_data; ++q) + if (stabilizers[row_base + q] != 0) + support_weight++; + std::size_t fb_weight = + feedback_flat.size() > 0 + ? row_support_weight(feedback_flat, record_row, num_cols) + : 0; + std::vector support(1 + support_weight + fb_weight); + support[0] = last_syndrome[record_row]; + std::size_t idx = 1; + for (std::size_t q = 0; q < num_data; ++q) + if (stabilizers[row_base + q] != 0) + support[idx++] = data_results[q]; + if (feedback_flat.size() > 0) + for (std::size_t k = 0; k < num_cols; ++k) + if (feedback_flat[record_row * num_cols + k] != 0) + support[idx++] = last_syndrome[k]; + return support; +} + +/// @brief Round-major slice offsets into the observable-feedback buffer. Entry +/// m is the start index of observable m's slice; each observable occupies +/// `num_rounds * (weight of obs_feedback_flat row m)` records. The returned +/// vector has `num_observables + 1` entries; the last entry is the total +/// buffer size. Returns an empty vector when @p obs_feedback_flat is empty +/// (no observable feedback declared). +inline __qpu__ std::vector +observable_feedback_offsets(const std::vector &obs_feedback_flat, + std::size_t num_observables, std::size_t num_cols, + std::size_t num_rounds) { + // Sized (not default) construction: the __qpu__ dialect forbids the default + // std::vector constructor. An empty obs_feedback_flat yields a size-0 vector. + std::size_t num_offsets = + obs_feedback_flat.size() > 0 ? num_observables + 1 : 0; + std::vector offsets(num_offsets); + if (obs_feedback_flat.size() == 0) + return offsets; + std::size_t total = 0; + for (std::size_t m = 0; m < num_observables; ++m) { + offsets[m] = total; + total += num_rounds * row_support_weight(obs_feedback_flat, m, num_cols); + } + offsets[num_observables] = total; + return offsets; +} + +/// @brief Write the records that observable feedback collects during @p round +/// into the round-major @p obs_fb_records buffer. For each observable m, the +/// round-r block starts at `offsets[m] + round * (weight of obs_feedback_flat +/// row m)` and holds `syndrome[k]` for each column k (ascending) with +/// `obs_feedback_flat(m, k) != 0`. Call with round 0 for the first round's +/// syndrome and with the running round index for each subsequent round. +inline __qpu__ void collect_observable_feedback_round( + std::vector &obs_fb_records, + const std::vector &offsets, + const std::vector &syndrome, std::size_t round, + const std::vector &obs_feedback_flat, + std::size_t num_observables, std::size_t num_cols) { + for (std::size_t m = 0; m < num_observables; ++m) { + std::size_t weight = row_support_weight(obs_feedback_flat, m, num_cols); + std::size_t idx = offsets[m] + round * weight; + for (std::size_t k = 0; k < num_cols; ++k) + if (obs_feedback_flat[m * num_cols + k] != 0) + obs_fb_records[idx++] = syndrome[k]; + } +} + +/// @brief Support records for logical observable @p obs: the feedback records +/// collected in @p obs_fb_records (round-major, occupying +/// `[offsets[obs], offsets[obs + 1])`), followed by the data-qubit readouts +/// `data_results[q]` for each data qubit q (ascending) with +/// `obs_matrix_flat(obs, q) != 0`. When @p offsets is empty (no observable +/// feedback declared) the feedback prefix is empty and only the data support +/// is returned. +inline __qpu__ std::vector observable_support_records( + std::size_t obs, const std::vector &obs_matrix_flat, + std::size_t num_data, + const std::vector &data_results, + const std::vector &obs_fb_records, + const std::vector &offsets) { + std::size_t support_weight = + row_support_weight(obs_matrix_flat, obs, num_data); + std::size_t fb_count = 0; + std::size_t base = 0; + if (offsets.size() > 0) { + base = offsets[obs]; + fb_count = offsets[obs + 1] - offsets[obs]; + } + std::vector obs_support(fb_count + support_weight); + for (std::size_t i = 0; i < fb_count; ++i) + obs_support[i] = obs_fb_records[base + i]; + std::size_t idx = fb_count; + for (std::size_t q = 0; q < num_data; ++q) + if (obs_matrix_flat[obs * num_data + q] != 0) + obs_support[idx++] = data_results[q]; + return obs_support; +} + +} // namespace cudaq::qec diff --git a/libs/qec/lib/device/memory_circuit.cpp b/libs/qec/lib/device/memory_circuit.cpp index b00754a2..ad1a8811 100644 --- a/libs/qec/lib/device/memory_circuit.cpp +++ b/libs/qec/lib/device/memory_circuit.cpp @@ -6,7 +6,7 @@ * the terms of the Apache License 2.0 which accompanies this distribution. * ******************************************************************************/ #include "memory_circuit.h" -#include +#include "inlined_feedback.h" namespace cudaq::qec { @@ -18,7 +18,9 @@ __qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round, const std::vector &z_stabilizers, const std::vector &obs_matrix_flat, std::size_t num_observables, - bool measure_in_x_basis) { + bool measure_in_x_basis, + const std::vector &feedback_flat, + const std::vector &obs_feedback_flat) { // Allocate the data and ancilla qubits cudaq::qvector data(num_data), xstab_anc(numAncx), zstab_anc(numAncz); @@ -39,10 +41,42 @@ __qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round, cudaq::detector(final_syndrome[fixed_offset + i]); } + std::size_t numCols = numAncx + numAncz; + + // Observable inlined feedback accumulation. Nested std::vector is not + // supported in __qpu__ code, so instead of one record vector per observable + // we use a single flat buffer with per-observable round-major slices; see + // observable_feedback_offsets / collect_observable_feedback_round in + // inlined_feedback.h. + std::vector obs_fb_offsets = observable_feedback_offsets( + obs_feedback_flat, num_observables, numCols, num_rounds); + std::size_t obs_fb_total = + obs_feedback_flat.size() > 0 ? obs_fb_offsets[num_observables] : 0; + std::vector obs_fb_records(obs_fb_total); + + // Collect the first-round feedback records for each observable. + if (obs_feedback_flat.size() > 0) + collect_observable_feedback_round(obs_fb_records, obs_fb_offsets, + final_syndrome, 0, obs_feedback_flat, + num_observables, numCols); + // Generate syndrome data for (std::size_t round = 1; round < num_rounds; ++round) { auto syndrome = stabilizer_round(logical, x_stabilizers, z_stabilizers); - cudaq::detectors(final_syndrome, syndrome); + if (obs_feedback_flat.size() > 0) + collect_observable_feedback_round(obs_fb_records, obs_fb_offsets, + syndrome, round, obs_feedback_flat, + num_observables, numCols); + if (feedback_flat.size() == 0) { + cudaq::detectors(final_syndrome, syndrome); + } else { + // Cross-round detector for record j: earlier vs current round record, + // augmented with the earlier-round records declared in row j of the + // feedback matrix. + for (std::size_t j = 0; j < numCols; ++j) + cudaq::detector(cross_round_detector_records( + final_syndrome, syndrome, j, feedback_flat, numCols)); + } final_syndrome = syndrome; } @@ -51,48 +85,23 @@ __qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round, } auto data_results = mz(data); - // Emit one logical_observable per row of the observable matrix. - for (std::size_t obs = 0; obs < num_observables; ++obs) { - std::size_t support_weight = 0; - for (std::size_t q = 0; q < num_data; ++q) { - if (obs_matrix_flat[obs * num_data + q] != 0) - support_weight++; - } - std::vector obs_support(support_weight); - std::size_t idx = 0; - for (std::size_t q = 0; q < num_data; ++q) { - if (obs_matrix_flat[obs * num_data + q] != 0) - obs_support[idx++] = data_results[q]; - } - cudaq::logical_observable(obs_support); - } + // Emit one logical_observable per row of the observable matrix: the per-round + // feedback records followed by the data-qubit support. + for (std::size_t obs = 0; obs < num_observables; ++obs) + cudaq::logical_observable( + observable_support_records(obs, obs_matrix_flat, num_data, data_results, + obs_fb_records, obs_fb_offsets)); // For each stabilizer, form detectors from data qubit readout connected with - // final stabilizer round. + // final stabilizer round. With inlined feedback, the boundary detector for + // record (fixed_offset + x) is extended with the declared last-round records. const std::vector &stabilizers = measure_in_x_basis ? x_stabilizers : z_stabilizers; - for (std::size_t x = 0; x < num_fixed_measurements; ++x) { - std::size_t row_base = x * num_data; - - std::size_t support_weight = 0; - for (std::size_t q = 0; q < num_data; ++q) { - if (stabilizers[row_base + q] != 0) { - support_weight++; - } - } - - std::vector support(support_weight + 1); - support[0] = final_syndrome[fixed_offset + x]; - std::size_t support_idx = 1; - for (std::size_t q = 0; q < num_data; ++q) { - if (stabilizers[row_base + q] != 0) { - support[support_idx++] = data_results[q]; - } - } - - cudaq::detector(support); - } + for (std::size_t x = 0; x < num_fixed_measurements; ++x) + cudaq::detector(boundary_detector_records( + final_syndrome, fixed_offset + x, data_results, stabilizers, + x * num_data, num_data, feedback_flat, numCols)); } } // namespace cudaq::qec diff --git a/libs/qec/lib/device/memory_circuit.h b/libs/qec/lib/device/memory_circuit.h index 1f0280d9..a15b2584 100644 --- a/libs/qec/lib/device/memory_circuit.h +++ b/libs/qec/lib/device/memory_circuit.h @@ -29,6 +29,16 @@ namespace cudaq::qec { /// (num_observables × numData entries, values 0/1). /// @param num_observables Number of rows in the observable matrix (k). /// @param measure_in_x_basis Performing X- or Z-memory circuit +/// @param feedback_flat Row-major flattened inlined feedback matrix for +/// detectors. Size is 0 (no feedback; legacy detector structure) or +/// numCols*numCols with numCols = numAncx + numAncz. Entry (j, k) = 1 +/// means the cross-round detector for record j additionally XORs +/// record k of the earlier round, and the final boundary detector for +/// record j additionally XORs record k of the last round. +/// @param obs_feedback_flat Row-major flattened inlined feedback matrix for +/// logical observables. Size is 0 (no feedback) or +/// num_observables*numCols. Entry (m, k) = 1 means logical observable +/// m additionally XORs record k of every round. __qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round, const code::one_qubit_encoding &statePrep, std::size_t numData, std::size_t numAncx, @@ -37,5 +47,7 @@ __qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round, const std::vector &z_stabilizers, const std::vector &obs_matrix_flat, std::size_t num_observables, - bool measure_in_x_basis); + bool measure_in_x_basis, + const std::vector &feedback_flat, + const std::vector &obs_feedback_flat); } // namespace cudaq::qec diff --git a/libs/qec/lib/experiments.cpp b/libs/qec/lib/experiments.cpp index 9dbb9d0d..d4a9f178 100644 --- a/libs/qec/lib/experiments.cpp +++ b/libs/qec/lib/experiments.cpp @@ -254,6 +254,29 @@ stabilizer_detector_ranges(std::size_t numRounds, std::size_t numZStabs, return ranges; } +/// @brief Validate and flatten (row-major) an inlined feedback tensor +/// declared by a code. An empty tensor (the identity default) flattens to an +/// empty vector; a non-empty tensor must have shape +/// [expected_rows, expected_cols]. +static std::vector +flatten_feedback_tensor(const cudaqx::tensor &t, + std::size_t expected_rows, std::size_t expected_cols, + const std::string &name) { + if (t.size() == 0) + return {}; + if (t.rank() != 2 || t.shape()[0] != expected_rows || + t.shape()[1] != expected_cols) { + std::string actual; + for (std::size_t d = 0; d < t.rank(); ++d) + actual += (d ? ", " : "") + std::to_string(t.shape()[d]); + throw std::runtime_error(name + " has invalid shape [" + actual + + "] - expected [" + std::to_string(expected_rows) + + ", " + std::to_string(expected_cols) + + "] or an empty tensor."); + } + return std::vector(t.data(), t.data() + t.size()); +} + /// @brief Given a memory circuit setup, sample syndrome and data qubit /// measurements. This is the main driver function that all of the /// sample_memory_circuit overloads invoke. @@ -309,6 +332,15 @@ sample_memory_circuit(const code &code, operation statePrep, const std::size_t numCols = numAncx + numAncz; + auto feedback_flat = flatten_feedback_tensor( + code.get_inlined_feedback(), numCols, numCols, "get_inlined_feedback()"); + auto obs_feedback_flat = flatten_feedback_tensor( + is_z_prep ? code.get_observable_inlined_feedback_z() + : code.get_observable_inlined_feedback_x(), + num_obs, numCols, + is_z_prep ? "get_observable_inlined_feedback_z()" + : "get_observable_inlined_feedback_x()"); + // Obtain the Measurement-to-Detector (M2D) sparse matrix. // m2d.rows[d] = set of chronological measurement indices whose XOR = detector // d. @@ -316,14 +348,16 @@ sample_memory_circuit(const code &code, operation statePrep, cudaq::M2OSparseMatrix m2o; cudaq::dem_from_kernel(memory_circuit, &noise, /*options=*/{}, m2d, m2o, stabRound, prep, numData, numAncx, numAncz, numRounds, - xVec, zVec, obs_flat, num_obs, !is_z_prep); + xVec, zVec, obs_flat, num_obs, !is_z_prep, + feedback_flat, obs_feedback_flat); // Sample the memory circuit and collect all raw measurements. cudaq::sample_options opts{ .shots = numShots, .noise = noise, .explicit_measurements = true}; - auto result = cudaq::sample(opts, memory_circuit, stabRound, prep, numData, - numAncx, numAncz, numRounds, xVec, zVec, obs_flat, - num_obs, !is_z_prep); + auto result = + cudaq::sample(opts, memory_circuit, stabRound, prep, numData, numAncx, + numAncz, numRounds, xVec, zVec, obs_flat, num_obs, + !is_z_prep, feedback_flat, obs_feedback_flat); // mzTable[shot, meas_idx] = raw 0/1 outcome; shape (numShots, // numMeasPerShot). Measurement layout per shot: numRounds*numCols ancilla, @@ -479,11 +513,22 @@ dem_from_memory_circuit(const code &code, operation statePrep, std::vector obs_flat(logical_obs.data(), logical_obs.data() + logical_obs.size()); + const std::size_t numAnc = numAncx + numAncz; + auto feedback_flat = flatten_feedback_tensor( + code.get_inlined_feedback(), numAnc, numAnc, "get_inlined_feedback()"); + auto obs_feedback_flat = flatten_feedback_tensor( + is_z_prep ? code.get_observable_inlined_feedback_z() + : code.get_observable_inlined_feedback_x(), + num_obs, numAnc, + is_z_prep ? "get_observable_inlined_feedback_z()" + : "get_observable_inlined_feedback_x()"); + cudaq::dem_options dem_opts; dem_opts.decompose_errors = decompose_errors; auto dem_text = cudaq::dem_from_kernel( memory_circuit, &noise, dem_opts, stabRound, prep, numData, numAncx, - numAncz, numRounds, xVec, zVec, obs_flat, num_obs, !is_z_prep); + numAncz, numRounds, xVec, zVec, obs_flat, num_obs, !is_z_prep, + feedback_flat, obs_feedback_flat); auto dem = cudaq::qec::dem_from_stim_text(dem_text, decompose_errors); const auto numXStabs = code.get_num_x_stabilizers(); diff --git a/libs/qec/python/bindings/py_code.cpp b/libs/qec/python/bindings/py_code.cpp index 8e3c493d..52a59b43 100644 --- a/libs/qec/python/bindings/py_code.cpp +++ b/libs/qec/python/bindings/py_code.cpp @@ -29,11 +29,61 @@ using namespace cudaqx; namespace cudaq::qec { +/// @brief Convert the array returned by a Python inlined-feedback getter to +/// an owning cudaqx::tensor. +/// @details The result is coerced to a C-contiguous uint8 numpy array first, +/// so any integer-convertible input (other dtypes, nested lists) is accepted. +/// toTensor only borrows the numpy buffer, so the data is copied into the +/// returned tensor to outlive the temporary Python array. +static cudaqx::tensor feedbackArrayToTensor(nb::object array, + const char *methodName) { + auto np = nb::module_::import_("numpy"); + nb::object coerced = np.attr("ascontiguousarray")(array, np.attr("uint8")); + auto ndarray = nb::cast>(coerced); + if (ndarray.ndim() != 2) + throw std::runtime_error(std::string(methodName) + + " must return a 2-D array."); + auto borrowed = cudaqx::toTensor(ndarray); + cudaqx::tensor owned; + owned.copy(borrowed.data(), borrowed.shape()); + return owned; +} + class PyCode : public qec::code { public: - NB_TRAMPOLINE(qec::code, 6); + NB_TRAMPOLINE(qec::code, 9); protected: + // Trampoline methods for the optional inlined-feedback virtuals: forward to + // a Python override when one exists (converting its numpy result), else + // fall back to the base-class empty default. This is NB_OVERRIDE with a + // custom return conversion, since cudaqx::tensor has no type caster. + cudaqx::tensor get_inlined_feedback() const override { + nb::detail::ticket nb_ticket(nb_trampoline, "get_inlined_feedback", false); + if (nb_ticket.key.is_valid()) + return feedbackArrayToTensor(nb_trampoline.base().attr(nb_ticket.key)(), + "get_inlined_feedback"); + return NBBase::get_inlined_feedback(); + } + + cudaqx::tensor get_observable_inlined_feedback_z() const override { + nb::detail::ticket nb_ticket(nb_trampoline, + "get_observable_inlined_feedback_z", false); + if (nb_ticket.key.is_valid()) + return feedbackArrayToTensor(nb_trampoline.base().attr(nb_ticket.key)(), + "get_observable_inlined_feedback_z"); + return NBBase::get_observable_inlined_feedback_z(); + } + + cudaqx::tensor get_observable_inlined_feedback_x() const override { + nb::detail::ticket nb_ticket(nb_trampoline, + "get_observable_inlined_feedback_x", false); + if (nb_ticket.key.is_valid()) + return feedbackArrayToTensor(nb_trampoline.base().attr(nb_ticket.key)(), + "get_observable_inlined_feedback_x"); + return NBBase::get_observable_inlined_feedback_x(); + } + // Trampoline methods for pure virtual functions std::size_t get_num_data_qubits() const override { NB_OVERRIDE_PURE(get_num_data_qubits); @@ -186,6 +236,32 @@ class PyCodeHandle : public qec::code { return m_py_operation_encodings; } + /// @brief Forward the optional inlined-feedback declarations to the + /// registered Python code when it defines them; otherwise return the + /// base-class empty default (identity detector layout). + cudaqx::tensor get_inlined_feedback() const override { + if (!nb::hasattr(pyCode, "get_inlined_feedback")) + return code::get_inlined_feedback(); + return feedbackArrayToTensor(pyCode.attr("get_inlined_feedback")(), + "get_inlined_feedback"); + } + + cudaqx::tensor get_observable_inlined_feedback_z() const override { + if (!nb::hasattr(pyCode, "get_observable_inlined_feedback_z")) + return code::get_observable_inlined_feedback_z(); + return feedbackArrayToTensor( + pyCode.attr("get_observable_inlined_feedback_z")(), + "get_observable_inlined_feedback_z"); + } + + cudaqx::tensor get_observable_inlined_feedback_x() const override { + if (!nb::hasattr(pyCode, "get_observable_inlined_feedback_x")) + return code::get_observable_inlined_feedback_x(); + return feedbackArrayToTensor( + pyCode.attr("get_observable_inlined_feedback_x")(), + "get_observable_inlined_feedback_x"); + } + protected: // Trampoline methods for pure virtual functions std::size_t get_num_data_qubits() const override { diff --git a/libs/qec/python/tests/test_code.py b/libs/qec/python/tests/test_code.py index 0280d037..0a08c3bc 100644 --- a/libs/qec/python/tests/test_code.py +++ b/libs/qec/python/tests/test_code.py @@ -315,5 +315,228 @@ def test_version(): assert "CUDA-Q QEC" in qec.__version__ +# Kernels for the inlined-feedback toy code: two data qubits with stabilizers +# XX and ZZ, both extracted through a single superdense (Bell-pair) ancilla +# pair (ancx[0] carries the XX record, ancz[0] the ZZ record). The round ends +# with an *uncorrected* record-conditioned byproduct: after the Bell decode, +# ancx[0] sits in the computational basis holding its future measurement +# outcome r_X, so the trailing CX(ancx[0], data[1]) is the unitary equivalent +# of the classically controlled byproduct X^{r_X} on data[1] that a hardware +# frame update would otherwise track. Same gadget and feedback matrices as the +# C++ toy in unittests/feedback_toy_device.cpp / test_qec.cpp (see the full +# derivation there): feedback = [[0, 1], [0, 0]], observable feedback = +# [[0, 1]], with records per round in [Z][X] order. +@cudaq.kernel +def _feedback_toy_prep0(q: patch): + reset(q.data[0]) + reset(q.data[1]) + + +@cudaq.kernel +def _feedback_toy_round(q: patch, x_stabilizers: list[int], + z_stabilizers: list[int]) -> list[cudaq.measure_handle]: + # Bell-prepare the superdense ancilla pair. + h(q.ancx[0]) + x.ctrl(q.ancx[0], q.ancz[0]) + # Couple XX through the X ancilla and ZZ through the Z ancilla. + x.ctrl(q.ancx[0], q.data[0]) + x.ctrl(q.ancx[0], q.data[1]) + x.ctrl(q.data[0], q.ancz[0]) + x.ctrl(q.data[1], q.ancz[0]) + # Decode the Bell pair. + x.ctrl(q.ancx[0], q.ancz[0]) + h(q.ancx[0]) + # Uncorrected record-conditioned byproduct: X^{r_X} on data[1]. + x.ctrl(q.ancx[0], q.data[1]) + # Records in [Z][X] order. + results = mz([*q.ancz, *q.ancx]) + reset(q.ancz[0]) + reset(q.ancx[0]) + return results + + +def test_python_inlined_feedback_toy(): + + @qec.code('py-feedback-toy') + class FeedbackToy: + + def __init__(self): + qec.Code.__init__(self) + self.stabilizers = [ + cudaq.SpinOperator.from_word(w) for w in ["XX", "ZZ"] + ] + self.pauli_observables = [cudaq.SpinOperator.from_word("ZZ")] + self.operation_encodings = { + qec.operation.prep0: _feedback_toy_prep0, + qec.operation.stabilizer_round: _feedback_toy_round + } + + def get_num_data_qubits(self): + return 2 + + def get_num_ancilla_qubits(self): + return 2 + + def get_num_ancilla_x_qubits(self): + return 1 + + def get_num_ancilla_z_qubits(self): + return 1 + + def get_num_x_stabilizers(self): + return 1 + + def get_num_z_stabilizers(self): + return 1 + + def get_inlined_feedback(self): + return np.array([[0, 1], [0, 0]], dtype=np.uint8) + + def get_observable_inlined_feedback_z(self): + # Exercise the bridge's dtype coercion: int64 instead of uint8. + return np.array([[0, 1]]) + + # Sample on stim, mirroring the C++ unit tests (test_qec links the stim + # target). The toy's X record is genuinely random in round 1, and the + # default (qpp) target does not preserve cross-round mid-circuit + # measurement correlations when sampling Python kernels - independent of + # the inlined-feedback machinery, but this toy (unlike the steane + # example, whose records are all deterministically 0) exposes it. + cudaq.set_target('stim') + try: + # numRounds >= 2 exercises the final boundary detector (same-round + # herald) distinctly from the cross-round detectors (earlier-round + # herald). + numShots, numRounds = 20, 4 + syndromes, data = qec.sample_memory_circuit( + qec.get_code('py-feedback-toy'), + numShots=numShots, + numRounds=numRounds) + + # 1 first-round boundary + 2 * (numRounds - 1) cross-round + 1 final + # boundary detectors, all deterministic (zero) in the noiseless + # circuit. + assert syndromes.shape == (numShots, 2 * numRounds) + assert not np.any(syndromes) + + # The individual data qubits are randomized by the XX measurement, + # but the ZZ parity must close deterministically every shot. + assert data.shape == (numShots, 2) + assert not np.any(data[:, 0] ^ data[:, 1]) + finally: + cudaq.reset_target() + + +def test_python_inlined_feedback_toy_negative_control(): + # The identical toy without the feedback declarations: the uncorrected + # byproduct makes the record-0 detectors non-deterministic and stim must + # reject the circuit. This proves the Python-declared feedback is doing + # the work in the positive test above. + @qec.code('py-feedback-toy-nofb') + class FeedbackToyNoFeedback: + + def __init__(self): + qec.Code.__init__(self) + self.stabilizers = [ + cudaq.SpinOperator.from_word(w) for w in ["XX", "ZZ"] + ] + self.pauli_observables = [cudaq.SpinOperator.from_word("ZZ")] + self.operation_encodings = { + qec.operation.prep0: _feedback_toy_prep0, + qec.operation.stabilizer_round: _feedback_toy_round + } + + def get_num_data_qubits(self): + return 2 + + def get_num_ancilla_qubits(self): + return 2 + + def get_num_ancilla_x_qubits(self): + return 1 + + def get_num_ancilla_z_qubits(self): + return 1 + + def get_num_x_stabilizers(self): + return 1 + + def get_num_z_stabilizers(self): + return 1 + + # stim's determinism analysis (a std::invalid_argument surfaced as + # ValueError) rejects the circuit inside dem_from_kernel, before any + # sampling happens. + with pytest.raises(ValueError, match="non-deterministic detectors"): + qec.sample_memory_circuit(qec.get_code('py-feedback-toy-nofb'), + numShots=20, + numRounds=4) + + +def test_python_inlined_feedback_observable_basis_selection(): + # The bridge routes get_observable_inlined_feedback_x to the X-basis getter, + # so it is never queried on the Z (prep0) path. Plant a deliberately + # wrong-shape _x matrix alongside a valid _z matrix: were the bridge to feed + # _x into the Z path, flatten_feedback_tensor's shape check would reject the + # circuit. The toy registers only prep0, so this proves selection from the + # Z path (the wrong-shape _x is inert there). + @qec.code('py-feedback-toy-basis-select') + class FeedbackToyBasisSelect: + + def __init__(self): + qec.Code.__init__(self) + self.stabilizers = [ + cudaq.SpinOperator.from_word(w) for w in ["XX", "ZZ"] + ] + self.pauli_observables = [cudaq.SpinOperator.from_word("ZZ")] + self.operation_encodings = { + qec.operation.prep0: _feedback_toy_prep0, + qec.operation.stabilizer_round: _feedback_toy_round + } + + def get_num_data_qubits(self): + return 2 + + def get_num_ancilla_qubits(self): + return 2 + + def get_num_ancilla_x_qubits(self): + return 1 + + def get_num_ancilla_z_qubits(self): + return 1 + + def get_num_x_stabilizers(self): + return 1 + + def get_num_z_stabilizers(self): + return 1 + + def get_inlined_feedback(self): + return np.array([[0, 1], [0, 0]], dtype=np.uint8) + + def get_observable_inlined_feedback_z(self): + return np.array([[0, 1]], dtype=np.uint8) + + def get_observable_inlined_feedback_x(self): + # Wrong shape for [num_obs=1 x numCols=2]; must never reach the + # Z-basis flatten on the prep0 path. + return np.zeros((3, 5), dtype=np.uint8) + + cudaq.set_target('stim') + try: + numShots, numRounds = 20, 4 + syndromes, data = qec.sample_memory_circuit( + qec.get_code('py-feedback-toy-basis-select'), + numShots=numShots, + numRounds=numRounds) + assert syndromes.shape == (numShots, 2 * numRounds) + assert not np.any(syndromes) + assert data.shape == (numShots, 2) + assert not np.any(data[:, 0] ^ data[:, 1]) + finally: + cudaq.reset_target() + + if __name__ == "__main__": pytest.main() diff --git a/libs/qec/unittests/CMakeLists.txt b/libs/qec/unittests/CMakeLists.txt index 344d0544..58883482 100644 --- a/libs/qec/unittests/CMakeLists.txt +++ b/libs/qec/unittests/CMakeLists.txt @@ -50,7 +50,16 @@ add_dependencies(CUDAQXQECUnitTests test_decoders_yaml) gtest_discover_tests(test_decoders_yaml) add_executable(test_qec test_qec.cpp) -target_link_libraries(test_qec PRIVATE GTest::gtest_main cudaq-qec cudaq::cudaq-stim-target) +# The inlined-feedback toy code's kernels are nvq++-compiled (like the +# registered codes' *_device.cpp files) so the memory-circuit kernel can +# resolve them through the kernel registry; cudaq::cudaq and nvqir provide +# the runtime symbols that the device object references. +cudaqx_add_device_code(test_qec + SOURCES + feedback_toy_device.cpp +) +target_link_directories(test_qec PRIVATE ${CUDAQ_INSTALL_DIR}/lib) +target_link_libraries(test_qec PRIVATE GTest::gtest_main cudaq-qec cudaq::cudaq nvqir cudaq::cudaq-stim-target) add_dependencies(CUDAQXQECUnitTests test_qec) gtest_discover_tests(test_qec) diff --git a/libs/qec/unittests/feedback_toy_device.cpp b/libs/qec/unittests/feedback_toy_device.cpp new file mode 100644 index 00000000..50b5e2a0 --- /dev/null +++ b/libs/qec/unittests/feedback_toy_device.cpp @@ -0,0 +1,56 @@ +/****************************************************************-*- C++ -*-**** + * Copyright (c) 2026 NVIDIA Corporation & Affiliates. * + * All rights reserved. * + * * + * This source code and the accompanying materials are made available under * + * the terms of the Apache License 2.0 which accompanies this distribution. * + ******************************************************************************/ + +#include "cudaq.h" +#include "cudaq/qec/patch.h" + +// Device kernels for the inlined-feedback toy code used by test_qec.cpp +// (see the QECCodeTester.checkInlinedFeedbackToy* tests). This file is +// nvq++-compiled, like the registered codes' *_device.cpp files, so the +// kernels are visible to the kernel registry that the memory-circuit +// machinery uses to resolve qkernel arguments. + +namespace cudaq::qec::feedback_toy { + +__qpu__ void prep0(patch p) { + for (std::size_t i = 0; i < p.data.size(); i++) + reset(p.data[i]); +} + +// A superdense (Bell-pair) round for a 2-data-qubit code with stabilizers XX +// and ZZ that deliberately ends with an *uncorrected* record-conditioned +// byproduct: after the Bell decode, ancx[0] sits in the computational basis +// holding its future measurement outcome r_X, so the trailing +// CX(ancx[0], data[1]) is the unitary equivalent of the classically +// controlled byproduct X^{r_X} on data[1] that a hardware frame update would +// otherwise track. The stabilizer supports are hardcoded (a single XX and ZZ +// plaquette), so the x/z_stabilizers arguments are intentionally unused. +__qpu__ std::vector +stabilizer_round(patch p, const std::vector &x_stabilizers, + const std::vector &z_stabilizers) { + // Bell-prepare the superdense ancilla pair. + h(p.ancx[0]); + cudaq::x(p.ancx[0], p.ancz[0]); + // Couple XX through the X ancilla and ZZ through the Z ancilla. + cudaq::x(p.ancx[0], p.data[0]); + cudaq::x(p.ancx[0], p.data[1]); + cudaq::x(p.data[0], p.ancz[0]); + cudaq::x(p.data[1], p.ancz[0]); + // Decode the Bell pair. + cudaq::x(p.ancx[0], p.ancz[0]); + h(p.ancx[0]); + // Uncorrected record-conditioned byproduct: X^{r_X} on data[1]. + cudaq::x(p.ancx[0], p.data[1]); + // Records in [Z][X] order. + auto results = mz(p.ancz, p.ancx); + reset(p.ancz[0]); + reset(p.ancx[0]); + return results; +} + +} // namespace cudaq::qec::feedback_toy diff --git a/libs/qec/unittests/test_qec.cpp b/libs/qec/unittests/test_qec.cpp index 47326f08..4b44344d 100644 --- a/libs/qec/unittests/test_qec.cpp +++ b/libs/qec/unittests/test_qec.cpp @@ -11,9 +11,11 @@ #include #include #include +#include #include #include "cuda-qx/core/library_utils.h" +#include "cudaq/algorithms/dem.h" #include "cudaq/qec/codes/surface_code.h" #include "cudaq/qec/experiments.h" #include "cudaq/qec/pcm_utils.h" @@ -2192,3 +2194,423 @@ TEST(PluginLoaderTester, checkCleanupPluginsEdgeCases) { // not do anything. cudaq::qec::cleanup_plugins(cudaq::qec::PluginType::CODE); } + +TEST(QECCodeTester, checkInlinedFeedbackDefaults) { + auto steane = cudaq::qec::get_code("steane"); + auto feedback = steane->get_inlined_feedback(); + EXPECT_EQ(feedback.rank(), 0); + EXPECT_EQ(feedback.size(), 0); + auto obs_feedback_z = steane->get_observable_inlined_feedback_z(); + EXPECT_EQ(obs_feedback_z.rank(), 0); + EXPECT_EQ(obs_feedback_z.size(), 0); + auto obs_feedback_x = steane->get_observable_inlined_feedback_x(); + EXPECT_EQ(obs_feedback_x.rank(), 0); + EXPECT_EQ(obs_feedback_x.size(), 0); +} + +// Device kernels for the inlined-feedback toy code, nvq++-compiled in +// feedback_toy_device.cpp (host-compiled kernels cannot be resolved through +// the kernel registry that the memory-circuit kernel uses for its qkernel +// arguments). +namespace cudaq::qec::feedback_toy { +__qpu__ void prep0(patch p); +__qpu__ std::vector +stabilizer_round(patch p, const std::vector &x_stabilizers, + const std::vector &z_stabilizers); +} // namespace cudaq::qec::feedback_toy + +namespace { + +// --------------------------------------------------------------------------- +// Toy code for the inlined-feedback tests below. +// +// Two data qubits with stabilizers XX and ZZ and logical observable ZZ. Both +// stabilizers are extracted through a single superdense (Bell-pair) ancilla +// pair: ancx[0] carries the XX record and ancz[0] the ZZ record. The round +// (defined in feedback_toy_device.cpp) deliberately ends with an +// *uncorrected* record-conditioned byproduct: after the Bell decode, ancx[0] +// sits in the computational basis holding its future measurement outcome +// r_X, so the trailing CX(ancx[0], data[1]) is the unitary equivalent of the +// classically-controlled byproduct X^{r_X} on data[1] that a hardware frame +// update would otherwise track. +// +// Derivation of the feedback matrices (records per round in [Z][X] order: +// record 0 = ancz outcome r_Z, record 1 = ancx outcome r_X): +// - r_Z(t) reads the ZZ parity of the data entering round t, i.e. the +// accumulated byproduct parity b_t, with b_1 = 0 and +// b_{t+1} = b_t XOR r_X(t). +// - r_X(t) reads XX: uniformly random in round 1 of a Z-basis memory and +// constant afterwards (the X byproduct commutes with XX). +// - Cross-round detector for record 0 without feedback: +// r_Z(t) XOR r_Z(t+1) = r_X(t) -> non-deterministic (the negative-control +// test relies on stim rejecting exactly this). XOR-ing in the earlier +// round's record 1 fixes it: feedback(0, 1) = 1. +// - Record 1 compares deterministically on its own: feedback row 1 = 0. +// - Final boundary detector for record 0: +// r_Z(T) XOR d0 XOR d1 = b_T XOR b_{T+1} = r_X(T); the same +// feedback(0, 1) entry extends the boundary detector with the last +// round's record 1. +// - The ZZ observable closes as d0 XOR d1 = b_{T+1} = XOR of every round's +// record 1: observable_feedback(0, 1) = 1. +// Note that the byproduct must be heralded by the *other* (X) record: a +// byproduct conditioned on the Z record - which is identically 0 in the +// noiseless circuit - would never fire and could not break determinism. +// Both matrices were validated against stim's determinism analysis for +// 1 through 6 rounds. +// --------------------------------------------------------------------------- + +class feedback_toy_code : public cudaq::qec::code { +protected: + std::size_t get_num_data_qubits() const override { return 2; } + std::size_t get_num_ancilla_qubits() const override { return 2; } + std::size_t get_num_ancilla_x_qubits() const override { return 1; } + std::size_t get_num_ancilla_z_qubits() const override { return 1; } + std::size_t get_num_x_stabilizers() const override { return 1; } + std::size_t get_num_z_stabilizers() const override { return 1; } + + bool declare_feedback; + +public: + explicit feedback_toy_code(bool declare_feedback) + : declare_feedback(declare_feedback) { + operation_encodings.insert(std::make_pair(cudaq::qec::operation::prep0, + cudaq::qec::feedback_toy::prep0)); + operation_encodings.insert( + std::make_pair(cudaq::qec::operation::stabilizer_round, + cudaq::qec::feedback_toy::stabilizer_round)); + m_stabilizers = fromPauliWords({"XX", "ZZ"}); + m_pauli_observables = fromPauliWords({"ZZ"}); + } + + cudaqx::tensor get_inlined_feedback() const override { + if (!declare_feedback) + return code::get_inlined_feedback(); + cudaqx::tensor feedback({2, 2}); + feedback.at({0, 1}) = 1; + return feedback; + } + + cudaqx::tensor get_observable_inlined_feedback_z() const override { + if (!declare_feedback) + return code::get_observable_inlined_feedback_z(); + cudaqx::tensor feedback({1, 2}); + feedback.at({0, 1}) = 1; + return feedback; + } +}; + +// Basis-selection variant of the toy: returns caller-planted observable- +// feedback matrices for the two bases so a test can put a deliberately +// wrong-shape matrix on the basis it expects sample_memory_circuit NOT to +// query. The detector feedback (get_inlined_feedback) stays valid, so the +// underlying Z-memory circuit is otherwise identical to the positive toy. +class feedback_toy_basis_select_code : public feedback_toy_code { + cudaqx::tensor m_obs_z; + cudaqx::tensor m_obs_x; + +public: + feedback_toy_basis_select_code(cudaqx::tensor obs_z, + cudaqx::tensor obs_x) + : feedback_toy_code(/*declare_feedback=*/true), m_obs_z(std::move(obs_z)), + m_obs_x(std::move(obs_x)) {} + + cudaqx::tensor get_observable_inlined_feedback_z() const override { + return m_obs_z; + } + cudaqx::tensor get_observable_inlined_feedback_x() const override { + return m_obs_x; + } +}; + +} // namespace + +TEST(QECCodeTester, checkInlinedFeedbackToyMemoryCircuit) { + feedback_toy_code toy(/*declare_feedback=*/true); + const std::size_t numShots = 20; + const std::size_t numRounds = 4; + auto [syndromes, data] = cudaq::qec::sample_memory_circuit( + toy, cudaq::qec::operation::prep0, numShots, numRounds); + + // 1 first-round boundary + 2 * (numRounds - 1) cross-round + 1 final + // boundary detectors. + ASSERT_EQ(syndromes.rank(), 2); + ASSERT_EQ(syndromes.shape()[0], numShots); + ASSERT_EQ(syndromes.shape()[1], 2 * numRounds); + for (std::size_t shot = 0; shot < numShots; ++shot) + for (std::size_t d = 0; d < syndromes.shape()[1]; ++d) + EXPECT_EQ(syndromes.at({shot, d}), 0) + << "shot " << shot << ", detector " << d; + + // The individual data qubits are randomized by the XX measurement, but the + // ZZ parity must close deterministically every shot. + ASSERT_EQ(data.rank(), 2); + ASSERT_EQ(data.shape()[0], numShots); + ASSERT_EQ(data.shape()[1], 2); + for (std::size_t shot = 0; shot < numShots; ++shot) + EXPECT_EQ(data.at({shot, 0}) ^ data.at({shot, 1}), 0) << "shot " << shot; +} + +TEST(QECCodeTester, checkInlinedFeedbackObservableBasisSelection) { + // The observable-feedback getter must be selected by the memory basis: a + // Z-basis (prep0) memory queries only get_observable_inlined_feedback_z() + // and never get_observable_inlined_feedback_x(). The toy device kernels + // register only prep0, so selection is proved from the Z path in both + // directions. + auto valid_obs = [] { + cudaqx::tensor t({1, 2}); + t.at({0, 1}) = 1; + return t; + }; + // Non-empty and the wrong shape for [num_obs=1 x numCols=2], so any code that + // flattens it trips flatten_feedback_tensor's shape check. + auto wrong_shape_obs = [] { return cudaqx::tensor({3, 5}); }; + + // (a) Valid Z matrix, deliberately wrong-shape X matrix: the Z-memory run + // must succeed exactly as the positive toy test, proving the X getter is + // never queried on the Z path. + { + feedback_toy_basis_select_code toy(valid_obs(), wrong_shape_obs()); + const std::size_t numShots = 20; + const std::size_t numRounds = 4; + auto [syndromes, data] = cudaq::qec::sample_memory_circuit( + toy, cudaq::qec::operation::prep0, numShots, numRounds); + ASSERT_EQ(syndromes.rank(), 2); + ASSERT_EQ(syndromes.shape()[0], numShots); + ASSERT_EQ(syndromes.shape()[1], 2 * numRounds); + for (std::size_t shot = 0; shot < numShots; ++shot) + for (std::size_t d = 0; d < syndromes.shape()[1]; ++d) + EXPECT_EQ(syndromes.at({shot, d}), 0) + << "shot " << shot << ", detector " << d; + ASSERT_EQ(data.shape()[0], numShots); + ASSERT_EQ(data.shape()[1], 2); + for (std::size_t shot = 0; shot < numShots; ++shot) + EXPECT_EQ(data.at({shot, 0}) ^ data.at({shot, 1}), 0) << "shot " << shot; + } + + // (converse) Wrong-shape Z matrix, valid X matrix: the same Z-memory run + // must now throw the shape-validation error, proving the Z getter IS the one + // queried on the Z path (the valid X matrix is irrelevant here). + { + feedback_toy_basis_select_code toy(wrong_shape_obs(), valid_obs()); + try { + cudaq::qec::sample_memory_circuit(toy, cudaq::qec::operation::prep0, + /*numShots=*/20, /*numRounds=*/4); + FAIL() << "expected sample_memory_circuit to reject the wrong-shape " + "Z-basis observable feedback"; + } catch (const std::exception &e) { + EXPECT_NE(std::string(e.what()).find( + "get_observable_inlined_feedback_z() has invalid shape"), + std::string::npos) + << "unexpected failure reason: " << e.what(); + } + } +} + +TEST(QECCodeTester, checkInlinedFeedbackToyNegativeControl) { + // The identical toy without the feedback declaration: the uncorrected + // byproduct makes the record-0 detectors non-deterministic and stim must + // reject the circuit. This proves the declared feedback is doing the work + // in the positive test above. + feedback_toy_code toy(/*declare_feedback=*/false); + try { + cudaq::qec::sample_memory_circuit(toy, cudaq::qec::operation::prep0, + /*numShots=*/20, /*numRounds=*/4); + FAIL() << "expected sample_memory_circuit to reject the feedback-less toy"; + } catch (const std::exception &e) { + // Fail-for-the-right-reason check: stim's determinism analysis must be + // what rejects the circuit (it flags the record-0 cross-round and final + // boundary detectors). + EXPECT_NE(std::string(e.what()).find("non-deterministic detectors"), + std::string::npos) + << "unexpected failure reason: " << e.what(); + } +} + +TEST(QECCodeTester, checkInlinedFeedbackToyDem) { + // The DEM path threads the same feedback declaration through + // dem_from_kernel; with circuit-level noise it must produce a valid model + // (a noiseless model has no error mechanisms and is rejected downstream). + feedback_toy_code toy(/*declare_feedback=*/true); + const std::size_t numRounds = 4; + cudaq::noise_model noise; + noise.add_all_qubit_channel("x", cudaq::qec::two_qubit_depolarization(0.01), + /*num_controls=*/1); + auto dem = cudaq::qec::dem_from_memory_circuit( + toy, cudaq::qec::operation::prep0, numRounds, noise); + // 1 first-round boundary + 2 * (numRounds - 1) cross-round + 1 final + // boundary detectors. + EXPECT_EQ(dem.num_detectors(), 2 * numRounds); + EXPECT_GT(dem.num_error_mechanisms(), 0); + EXPECT_EQ(dem.num_observables(), 1); +} + +TEST(InlinedFeedbackLayout, EmptyTensorsGiveEmptyLayout) { + cudaqx::tensor empty; + auto layout = cudaq::qec::build_inlined_feedback_layout(empty, empty, 2, 1); + EXPECT_TRUE(layout.detector_indices.empty()); + EXPECT_TRUE(layout.detector_offsets.empty()); + EXPECT_TRUE(layout.observable_indices.empty()); + EXPECT_TRUE(layout.observable_offsets.empty()); +} + +TEST(InlinedFeedbackLayout, ToyMatricesCsr) { + // fb = [[0,1],[0,0]]: record 0 heralded by record 1; row 1 empty. + cudaqx::tensor fb({2, 2}); + fb.at({0, 1}) = 1; + // obs_fb = [[0,1]]: observable 0 XORs record 1 every round. + cudaqx::tensor obs_fb({1, 2}); + obs_fb.at({0, 1}) = 1; + auto layout = cudaq::qec::build_inlined_feedback_layout(fb, obs_fb, 2, 1); + EXPECT_EQ(layout.detector_indices, (std::vector{1})); + EXPECT_EQ(layout.detector_offsets, (std::vector{0, 1, 1})); + EXPECT_EQ(layout.observable_indices, (std::vector{1})); + EXPECT_EQ(layout.observable_offsets, (std::vector{0, 1})); +} + +TEST(InlinedFeedbackLayout, MultiRowWeights) { + // 3x3 with row weights 2, 0, 1. + cudaqx::tensor fb({3, 3}); + fb.at({0, 0}) = 1; + fb.at({0, 2}) = 1; + fb.at({2, 1}) = 1; + cudaqx::tensor empty; + auto layout = cudaq::qec::build_inlined_feedback_layout(fb, empty, 3, 0); + EXPECT_EQ(layout.detector_indices, (std::vector{0, 2, 1})); + EXPECT_EQ(layout.detector_offsets, (std::vector{0, 2, 2, 3})); + EXPECT_TRUE(layout.observable_offsets.empty()); +} + +TEST(InlinedFeedbackLayout, AllZeroMatrixGivesZeroWeightRows) { + cudaqx::tensor fb({2, 2}); + cudaqx::tensor empty; + auto layout = cudaq::qec::build_inlined_feedback_layout(fb, empty, 2, 0); + EXPECT_TRUE(layout.detector_indices.empty()); + EXPECT_EQ(layout.detector_offsets, (std::vector{0, 0, 0})); +} + +TEST(InlinedFeedbackLayout, ThrowsOnWrongShape) { + cudaqx::tensor bad({2, 3}); // expected [2 x 2] + cudaqx::tensor empty; + EXPECT_THROW(cudaq::qec::build_inlined_feedback_layout(bad, empty, 2, 0), + std::runtime_error); + cudaqx::tensor bad_obs({2, 2}); // expected [1 x 2] + EXPECT_THROW(cudaq::qec::build_inlined_feedback_layout(empty, bad_obs, 2, 1), + std::runtime_error); +} + +// The entry-point memory-circuit kernel (defined in device/memory_circuit.cpp); +// forward-declared here - like the feedback_toy kernels above - so the +// conformance test below can drive dem_from_kernel directly. +namespace cudaq::qec { +__qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round, + const code::one_qubit_encoding &statePrep, + std::size_t numData, std::size_t numAncx, + std::size_t numAncz, std::size_t numRounds, + const std::vector &x_stabilizers, + const std::vector &z_stabilizers, + const std::vector &obs_matrix_flat, + std::size_t num_observables, + bool measure_in_x_basis, + const std::vector &feedback_flat, + const std::vector &obs_feedback_flat); +} // namespace cudaq::qec + +// Conformance: the device fold inside memory_circuit and the host-side rule in +// build_inlined_feedback_layout must agree on which records compose each +// detector. Build the measurements-to-detectors matrix (M2D) from the kernel +// and compare it, row by row, against the M2D reconstructed from the host +// layout plus the known emission order. If either the device fold or the host +// layout drifts, this test fails. +TEST(InlinedFeedbackLayout, HostLayoutMatchesKernelM2D) { + feedback_toy_code toy(/*declare_feedback=*/true); + const std::size_t numRounds = 4; + const std::size_t numData = 2, numAncx = 1, numAncz = 1; + const std::size_t numCols = numAncx + numAncz; // records per round + const std::size_t num_obs = 1; + + // --- Kernel-derived M2D, mirroring experiments.cpp's sample path. prep0 is a + // Z-basis memory, so measure_in_x_basis = false and the fixed records are the + // Z ancillas (fixed_offset = 0). --- + auto &prep = toy.get_operation( + cudaq::qec::operation::prep0); + auto &stabRound = toy.get_operation( + cudaq::qec::operation::stabilizer_round); + auto parity_x = toy.get_parity_x(); + auto parity_z = toy.get_parity_z(); + std::vector xVec(parity_x.data(), + parity_x.data() + parity_x.size()); + std::vector zVec(parity_z.data(), + parity_z.data() + parity_z.size()); + auto logical_obs = toy.get_observables_z(); + std::vector obs_flat(logical_obs.data(), + logical_obs.data() + logical_obs.size()); + + auto feedback = toy.get_inlined_feedback(); + auto obs_feedback = toy.get_observable_inlined_feedback_z(); + std::vector feedback_flat(feedback.data(), + feedback.data() + feedback.size()); + std::vector obs_feedback_flat( + obs_feedback.data(), obs_feedback.data() + obs_feedback.size()); + + cudaq::noise_model noise; + cudaq::M2DSparseMatrix m2d; + cudaq::M2OSparseMatrix m2o; + cudaq::dem_from_kernel( + cudaq::qec::memory_circuit, &noise, /*options=*/{}, m2d, m2o, stabRound, + prep, numData, numAncx, numAncz, numRounds, xVec, zVec, obs_flat, num_obs, + /*measure_in_x_basis=*/false, feedback_flat, obs_feedback_flat); + + // --- Host-layout-derived expected M2D. Global record index for round r, + // column c is r*numCols + c; the readout for data qubit q is at + // numRounds*numCols + q. --- + auto layout = cudaq::qec::build_inlined_feedback_layout( + feedback, obs_feedback, numCols, num_obs); + auto rec = [numCols](std::size_t round, std::size_t col) { + return round * numCols + col; + }; + auto herald_cols = [](const std::vector &indices, + const std::vector &offsets, + std::size_t row) { + std::vector cols; + if (!offsets.empty()) + cols.assign(indices.begin() + offsets[row], + indices.begin() + offsets[row + 1]); + return cols; + }; + + std::vector> expected; + // Round-0 fixed detector: the fixed record (fixed_offset + 0) of round 0. + expected.push_back({rec(0, 0)}); + // Interior cross-round detectors, emitted each round for j = 0..numCols-1: + // {prev_j, curr_j} plus the layout-row-j heralds taken from the earlier + // round. + for (std::size_t r = 1; r < numRounds; ++r) { + for (std::size_t j = 0; j < numCols; ++j) { + std::set row{rec(r - 1, j), rec(r, j)}; + for (auto k : + herald_cols(layout.detector_indices, layout.detector_offsets, j)) + row.insert(rec(r - 1, k)); + expected.push_back(row); + } + } + // Final boundary detector for the fixed record (row 0): the last-round fixed + // record, the data readouts in the stabilizer support, and the layout-row-0 + // heralds taken from the last round. + { + std::set row{rec(numRounds - 1, 0)}; + for (std::size_t q = 0; q < numData; ++q) + if (zVec[q] != 0) + row.insert(numRounds * numCols + q); + for (auto k : + herald_cols(layout.detector_indices, layout.detector_offsets, 0)) + row.insert(rec(numRounds - 1, k)); + expected.push_back(row); + } + + ASSERT_EQ(m2d.rows.size(), 2 * numRounds); + ASSERT_EQ(m2d.rows.size(), expected.size()); + for (std::size_t d = 0; d < expected.size(); ++d) { + std::set actual(m2d.rows[d].begin(), m2d.rows[d].end()); + EXPECT_EQ(actual, expected[d]) << "detector " << d; + } +}