Expose memory_circuit kernels

libs/qec/include/cudaq/qec/realtime/decoding.h

@@ -16,10 +16,21 @@
 // offload decoding work to our QEC decoders in real time
 // (within qubit coherence times)
 //
-// The design here is as follows: We declare but do not
-// implement the API. Then we allow users to specify concrete
-// implementations of the API via the target specification passed to
-// nvq++.
+// The design here is as follows: this header only declares the API; what a
+// kernel's calls bind to is resolved per binary.
+//
+//  - A target decoding shim (cudaq-qec-realtime-decoding-{simulation,
+//    quantinuum,...}), when linked, provides the real implementations. They
+//    win kernel-registry resolution because the shim registers after its
+//    libcudaq-qec dependency.
+//  - Without a shim, the calls bind to the default no-op implementations
+//    that libcudaq-qec registers (lib/device/memory_circuit.cpp): enqueues
+//    and resets do nothing, and corrections come back all-zero.
+//
+// In short: link a shim and these calls decode; link nothing extra and they
+// are no-ops. The host-side entry points additionally no-op while a circuit
+// is being analyzed for its detector error model, so kernels containing these
+// calls are safe to pass to dem_from_memory_circuit.
 
 namespace cudaq::qec::decoding {
 // CUDA-Q QEC Realtime Decoding API (declarations)

libs/qec/lib/CMakeLists.txt

@@ -122,6 +122,14 @@ if (CUDAQX_QEC_USE_FLOAT)
    target_compile_definitions(${LIBRARY_NAME} PUBLIC -DCUDAQX_QEC_FLOAT_TYPE=float)
 endif()
 
+# Host (g++) translation units -- this library's own and every consumer's, via
+# the target's usage requirements -- must see kernel declarations only when
+# including cudaq/qec/device/memory_circuit.h; the kernel *implementations* are
+# compiled solely by nvq++ device-code TUs (cudaqx_add_device_code does not
+# inherit compile definitions, so those still receive them).
+target_compile_definitions(${LIBRARY_NAME} PUBLIC
+  CUDAQ_QEC_MEMORY_CIRCUIT_DECLARATIONS_ONLY)
+
 target_include_directories(${DECODERS_LIBRARY_NAME}
   PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}
   PUBLIC

libs/qec/lib/code.cpp

@@ -7,8 +7,6 @@
  ******************************************************************************/
 #include "cudaq/qec/code.h"
 
-#include "device/memory_circuit.h"
-
 INSTANTIATE_REGISTRY(cudaq::qec::code, const cudaqx::heterogeneous_map &)
 
 namespace cudaq::qec {

libs/qec/lib/device/CMakeLists.txt

@@ -9,4 +9,6 @@
 cudaqx_add_device_code(cudaq-qec
   SOURCES 
     memory_circuit.cpp
+  COMPILER_FLAGS
+    -DCUDAQ_QEC_DISABLE_REALTIME_DECODING
 )

libs/qec/lib/device/memory_circuit.cpp

@@ -5,94 +5,44 @@
  * This source code and the accompanying materials are made available under    *
  * the terms of the Apache License 2.0 which accompanies this distribution.    *
  ******************************************************************************/
-#include "memory_circuit.h"
-#include <numeric>
 
-namespace cudaq::qec {
-
-__qpu__ void memory_circuit(const code::stabilizer_round &stabilizer_round,
-                            const code::one_qubit_encoding &statePrep,
-                            std::size_t num_data, std::size_t numAncx,
-                            std::size_t numAncz, std::size_t num_rounds,
-                            const std::vector<std::size_t> &x_stabilizers,
-                            const std::vector<std::size_t> &z_stabilizers,
-                            const std::vector<std::size_t> &obs_matrix_flat,
-                            std::size_t num_observables,
-                            bool measure_in_x_basis) {
-  // Allocate the data and ancilla qubits
-  cudaq::qvector data(num_data), xstab_anc(numAncx), zstab_anc(numAncz);
-
-  // Create the logical patch
-  patch logical(data, xstab_anc, zstab_anc);
-
-  // Prepare the initial state
-  statePrep({data, xstab_anc, zstab_anc});
-
-  // The "off-basis" detectors will be non-deterministic after the first
-  // stabilizer round.
-  auto final_syndrome = stabilizer_round(logical, x_stabilizers, z_stabilizers);
-  std::size_t num_fixed_measurements =
-      measure_in_x_basis ? xstab_anc.size() : zstab_anc.size();
-  std::size_t fixed_offset =
-      measure_in_x_basis ? final_syndrome.size() - num_fixed_measurements : 0;
-  for (std::size_t i = 0; i < num_fixed_measurements; ++i) {
-    cudaq::detector(final_syndrome[fixed_offset + i]);
-  }
-
-  // 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);
-    final_syndrome = syndrome;
-  }
-
-  if (measure_in_x_basis) {
-    h(data);
-  }
-  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<cudaq::measure_result> 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);
-  }
-
-  // For each stabilizer, form detectors from data qubit readout connected with
-  // final stabilizer round.
-  const std::vector<size_t> &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<cudaq::measure_result> 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);
-  }
+// The translation unit that compiles the `memory_circuit` device kernels into
+// libcudaq-qec. This TU must exist: it is what turns the header-defined kernels
+// into linkable symbols + JIT kernel registrations for the sample/DEM entry
+// points in experiments.cpp (and for applications, which only include the
+// declarations in cudaq/qec/device/memory_circuit.h).
+#include "cudaq/qec/device/memory_circuit.h"
+
+namespace cudaq::qec::decoding {
+
+// Default no-op implementations of the realtime decoding API (see
+// cudaq/qec/realtime/decoding.h), so that binaries which do not link a target
+// decoding shim still link and run: enqueues and resets do nothing, and
+// corrections come back all-zero.
+//
+// Deliberately defined in the SAME translation unit as the kernels above:
+// nvq++ inlines these empty bodies into the library's registered copy of the
+// kernels, making that copy deterministically decoding-free (it serves the
+// sample/DEM entry points). Realtime applications compile their own kernel
+// copy (by including memory_circuit.h in one nvq++ TU) and link a
+// decoding shim, whose implementations their copy binds to instead.
+
+__qpu__ void
+enqueue_syndromes(std::uint64_t decoder_id,
+                  const std::vector<cudaq::measure_result> &syndromes,
+                  std::uint64_t tag) {}
+
+__qpu__ void enqueue_syndromes_test(std::uint64_t decoder_id,
+                                    const std::vector<bool> &syndromes,
+                                    std::uint64_t tag) {}
+
+__qpu__ std::vector<bool> get_corrections(std::uint64_t decoder_id,
+                                          std::uint64_t return_size,
+                                          bool reset) {
+  std::vector<bool> result(return_size);
+  return result;
 }
 
-} // namespace cudaq::qec
+__qpu__ void reset_decoder(std::uint64_t decoder_id) {}
+
+} // namespace cudaq::qec::decoding

libs/qec/lib/experiments.cpp

@@ -7,7 +7,7 @@
  ******************************************************************************/
 
 #include "cudaq/qec/experiments.h"
-#include "device/memory_circuit.h"
+#include "cudaq/qec/device/memory_circuit.h"
 #include "cudaq/algorithms/dem.h"
 #include "cudaq/qec/dem_sampling.h"
 #include "cudaq/qec/pcm_utils.h"
@@ -316,14 +316,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,
+                         /*enqueue_decoder_id=*/std::int64_t{-1});
 
   // 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, /*enqueue_decoder_id=*/std::int64_t{-1});
 
   // mzTable[shot, meas_idx] = raw 0/1 outcome; shape (numShots,
   // numMeasPerShot). Measurement layout per shot: numRounds*numCols ancilla,
@@ -435,11 +437,10 @@ sample_memory_circuit(const code &code, std::size_t numShots,
 namespace details {
 /// @brief Given a memory circuit setup, generate a DEM. This is the main driver
 /// function that all of the function overloads invoke.
-cudaq::qec::detector_error_model
-dem_from_memory_circuit(const code &code, operation statePrep,
-                        std::size_t numRounds, cudaq::noise_model &noise,
-                        bool keep_x_stabilizers, bool keep_z_stabilizers,
-                        bool decompose_errors) {
+cudaq::qec::detector_error_model dem_from_memory_circuit(
+    const code &code, operation statePrep, std::size_t numRounds,
+    cudaq::noise_model &noise, bool keep_x_stabilizers, bool keep_z_stabilizers,
+    bool decompose_errors) {
   if (!keep_x_stabilizers && !keep_z_stabilizers)
     throw std::runtime_error("dem_from_memory_circuit error - no stabilizers "
                              "to keep.");
@@ -483,7 +484,8 @@ dem_from_memory_circuit(const code &code, operation statePrep,
   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,
+      /*enqueue_decoder_id=*/std::int64_t{-1});
   auto dem = cudaq::qec::dem_from_stim_text(dem_text, decompose_errors);
 
   const auto numXStabs = code.get_num_x_stabilizers();

libs/qec/lib/realtime/realtime_decoding.cpp

@@ -7,6 +7,7 @@
  ******************************************************************************/
 
 #include "realtime_decoding.h"
+#include "common/ExecutionContext.h"
 #include "cudaq/qec/decoder.h"
 #include "cudaq/qec/logger.h"
 #include "cudaq/qec/pcm_utils.h"
@@ -334,8 +335,15 @@ set_syndrome_capture_callback(void (*callback)(const uint8_t *, size_t)) {
   g_syndrome_capture_callback = callback;
 }
 
+static bool in_dem_analysis_context() {
+  auto *ctx = cudaq::getExecutionContext();
+  return ctx && ctx->name == "dem";
+}
+
 void enqueue_syndromes(std::size_t decoder_id, uint8_t *syndromes,
                        std::uint64_t syndrome_length, std::uint64_t tag) {
+  if (in_dem_analysis_context())
+    return;
   if (decoder_id >= g_decoders.size()) {
     throw std::invalid_argument(
         fmt::format("Decoder {} not found", decoder_id));
@@ -407,6 +415,13 @@ void enqueue_syndromes(std::size_t decoder_id, uint8_t *syndromes,
 
 void get_corrections(std::size_t decoder_id, uint8_t *corrections,
                      std::uint64_t correction_length, bool reset) {
+  if (in_dem_analysis_context()) {
+    // Circuit analysis (e.g. DEM extraction): return zeroed corrections rather
+    // than driving a decoder that isn't configured.
+    if (corrections && correction_length > 0)
+      std::memset(corrections, 0, correction_length);
+    return;
+  }
   CUDA_QEC_INFO("Entered get_corrections function decoder_id={}, "
                 "correction_length={}, reset={}",
                 decoder_id, correction_length, reset);
@@ -458,6 +473,8 @@ void get_corrections(std::size_t decoder_id, uint8_t *corrections,
 }
 
 void reset_decoder(std::size_t decoder_id) {
+  if (in_dem_analysis_context())
+    return;
   CUDA_QEC_INFO("Entered reset_decoder for decoder_id={}", decoder_id);
   if (decoder_id >= g_decoders.size()) {
     throw std::invalid_argument(

libs/qec/unittests/test_qec.cpp

@@ -15,6 +15,7 @@
 
 #include "cuda-qx/core/library_utils.h"
 #include "cudaq/qec/codes/surface_code.h"
+#include "cudaq/qec/device/memory_circuit.h"
 #include "cudaq/qec/experiments.h"
 #include "cudaq/qec/pcm_utils.h"
 #include "cudaq/qec/plugin_loader.h"
@@ -2192,3 +2193,18 @@ TEST(PluginLoaderTester, checkCleanupPluginsEdgeCases) {
   // not do anything.
   cudaq::qec::cleanup_plugins(cudaq::qec::PluginType::CODE);
 }
+
+TEST(QECCodeTester, checkMemoryCircuitMultiReadout) {
+  // memory_circuit_multi lays out logical qubit i's data readout at
+  // [i * numData, (i+1) * numData); the helper reshapes that flat buffer.
+  std::vector<int> flat{0, 1, 2, 3, 4, 5};
+  auto per_logical = cudaq::qec::memory_circuit_multi_readout(flat, 2, 3);
+  ASSERT_EQ(per_logical.size(), 2u);
+  EXPECT_EQ(per_logical[0], (std::vector<int>{0, 1, 2}));
+  EXPECT_EQ(per_logical[1], (std::vector<int>{3, 4, 5}));
+
+  // Degenerate single-logical case is the identity reshape.
+  auto single = cudaq::qec::memory_circuit_multi_readout(flat, 1, 6);
+  ASSERT_EQ(single.size(), 1u);
+  EXPECT_EQ(single[0], flat);
+}