Quantum Error Correction with Circuit-level Noise Modeling

This example builds upon the previous code-capacity noise model example. In the circuit-level noise modeling experiment, we have many of the same components from the CUDA-Q QEC library: QEC codes, decoders, and noisy data. The primary difference here, is that we can begin to run CUDA-Q kernels to generate noisy data, rather than just generating random bitstring to represent our errors.

Along with the stabilizers, parity check matrices, and logical observables, the QEC code type also has an encoding map. This map allows codes to define logical gates in terms of gates on the underlying physical qubits. These encodings operate on the qec.patch type, which represents three registers of physical qubits making up a logical qubit. A data qubit register, an X-stabilizer ancilla register, and a Z-stabilizer ancilla register.

The most notable encoding stored in the QEC map, is how the qec.operation.stabilizer_round, which encodes a cudaq.kernel which stores the gate-level information for how to do a stabilizer measurement. These stabilizer rounds are the gate-level way to encode the parity check matrix of a QEC code into quantum circuits.

This example walks through how to use the CUDA-Q QEC library to perform a quantum memory experiment simulation. These experiments model how well QEC cycles, or rounds of stabilizer measuments, can protect the information encoded in a logical qubit. If noise is turned off, then the information is protected indefinitely. Here, we will model depolarization noise after each CX gate, and track how many logical errors occur.

CUDA-Q QEC Implementation

Here’s how to use CUDA-Q QEC to perform a circuit-level noise model experiment in both Python and C++:

import numpy as np
import cudaq
import cudaq_qec as qec

# Get a QEC code
cudaq.set_target("stim")
steane = qec.get_code("steane")

# Get the parity check matrix of a code
# Can get the full code, or for CSS codes
# just the X or Z component
H = steane.get_parity()
print(f"H:\n{H}")
observables = steane.get_pauli_observables_matrix()
Lz = steane.get_observables_z()
print(f"observables:\n{observables}")
print(f"Lz:\n{Lz}")

nShots = 3
nRounds = 4

# error probabily
p = 0.01
noise = cudaq.NoiseModel()
noise.add_all_qubit_channel("x", qec.TwoQubitDepolarization(p), 1)

# prepare logical |0> state, tells the sampler to do z-basis experiment
statePrep = qec.operation.prep0
# our expected measurement in this state is 0
expected_value = 0

# sample the steane memory circuit with noise on each cx gate
# reading out the syndromes after each stabilizer round (xor'd against the previous)
# and readout out the data qubits at the end of the experiment
syndromes, data = qec.sample_memory_circuit(steane, statePrep, nShots, nRounds,
                                            noise)
print("From sample function:\n")
print("syndromes:\n", syndromes)
print("data:\n", data)

# Get a decoder
decoder = qec.get_decoder("single_error_lut", H)
nLogicalErrors = 0

# Logical Mz each shot (use Lx if preparing in X-basis)
logical_measurements = (Lz @ data.transpose()) % 2
# only one logical qubit, so do not need the second axis
logical_measurements = logical_measurements.flatten()
print("LMz:\n", logical_measurements)

# initialize a Pauli frame to track logical flips
# through the stabilizer rounds
pauli_frame = np.array([0, 0], dtype=np.uint8)
for shot in range(0, nShots):
    print("shot:", shot)
    for syndrome in syndromes:
        print("syndrome:", syndrome)
        # decode the syndrome
        convergence, result = decoder.decode(syndrome)
        data_prediction = np.array(result, dtype=np.uint8)

        # see if the decoded result anti-commutes with the observables
        print("decode result:", data_prediction)
        decoded_observables = (observables @ data_prediction) % 2
        print("decoded_observables:", decoded_observables)

        # update pauli frame
        pauli_frame = (pauli_frame + decoded_observables) % 2
        print("pauli frame:", pauli_frame)

    # after pauli frame has tracked corrections through the rounds
    # apply the pauli frame correction to the measurement, and see
    # if this matches the state we intended to prepare
    # We prepared |0>, so we check if logical measurement Mz + Pf_X = 0
    corrected_mz = (logical_measurements[shot] + pauli_frame[0]) % 2
    print("Expected value:", expected_value)
    print("Corrected value:", corrected_mz)
    if (corrected_mz != expected_value):
        nLogicalErrors += 1

# Count how many shots the decoder failed to correct the errors
print("Number of logical errors:", nLogicalErrors)

/*******************************************************************************
 * Copyright (c) 2024 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.    *
 ******************************************************************************/

// Compile and run with
// nvq++ --enable-mlir -lcudaq-qec circuit_level_noise.cpp -o circuit_level
// ./circuit_level

#include "cudaq.h"
#include "cudaq/qec/decoder.h"
#include "cudaq/qec/experiments.h"
#include "cudaq/qec/noise_model.h"

int main() {
  // Choose a QEC code
  auto steane = cudaq::qec::get_code("steane");

  // Access the parity check matrix
  auto H = steane->get_parity();
  std::cout << "H:\n";
  H.dump();

  // Access the logical observables
  auto observables = steane->get_pauli_observables_matrix();
  auto Lz = steane->get_observables_z();

  // Data qubits the logical Z observable is supported on
  std::cout << "Lz:\n";
  Lz.dump();

  // Observables are stacked as Z over X for mat-vec multiplication
  std::cout << "Obs:\n";
  observables.dump();

  // How many shots to run the experiment
  int nShots = 3;
  // For each shot, how many rounds of stabilizer measurements
  int nRounds = 4;

  // can set seed for reproducibility
  // cudaq::set_random_seed(1337);
  cudaq::noise_model noise;

  // Add a depolarization noise channel after each cx gate
  noise.add_all_qubit_channel(
      "x", cudaq::qec::two_qubit_depolarization(/*probability*/ 0.01),
      /*numControls*/ 1);

  // Perform a noisy z-basis memory circuit experiment
  auto [syndromes, data] = cudaq::qec::sample_memory_circuit(
      *steane, cudaq::qec::operation::prep0, nShots, nRounds, noise);

  // With noise, many syndromes will flip each QEC cycle, these are the
  // syndrome differences from the previous cycle.
  std::cout << "syndromes:\n";
  syndromes.dump();

  // With noise, Lz will sometimes be flipped
  std::cout << "data:\n";
  data.dump();

  // Use z-measurements on data qubits to determine the logical mz
  // In an x-basis experiment, use Lx.
  auto logical_mz = Lz.dot(data.transpose()) % 2;
  std::cout << "logical_mz each shot:\n";
  logical_mz.dump();

  // Select a decoder
  auto decoder = cudaq::qec::get_decoder("single_error_lut", H);

  // Initialize a pauli_frame to track the logical errors
  cudaqx::tensor<uint8_t> pauli_frame({observables.shape()[0]});

  // Start a loop to count the number of logical errors
  size_t numLerrors = 0;
  for (size_t shot = 0; shot < nShots; ++shot) {
    std::cout << "shot: " << shot << "\n";

    for (size_t round = 0; round < nRounds - 1; ++round) {
      std::cout << "round: " << round << "\n";

      // Access one row of the syndrome tensor
      size_t count = shot * (nRounds - 1) + round;
      size_t stride = syndromes.shape()[1];
      cudaqx::tensor<uint8_t> syndrome({stride});
      syndrome.borrow(syndromes.data() + stride * count);
      std::cout << "syndrome:\n";
      syndrome.dump();

      // Decode the syndrome
      auto [converged, v_result] = decoder->decode(syndrome);
      cudaqx::tensor<uint8_t> result_tensor;
      cudaq::qec::convert_vec_soft_to_tensor_hard(v_result, result_tensor);
      std::cout << "decode result:\n";
      result_tensor.dump();

      // See if the decoded result anti-commutes with observables
      auto decoded_observables = observables.dot(result_tensor);
      std::cout << "decoded observable:\n";
      decoded_observables.dump();

      // update from previous stabilizer round
      pauli_frame = (pauli_frame + decoded_observables) % 2;
      std::cout << "pauli frame:\n";
      pauli_frame.dump();
    }

    // prep0 means we expected to measure out 0.
    uint8_t expected_mz = 0;
    // Apply the pauli frame correction to our logical measurement
    uint8_t corrected_mz = (logical_mz.at({0, shot}) + pauli_frame.at({0})) % 2;

    // Check if Logical_mz + pauli_frame_X = 0?
    std::cout << "Corrected readout: " << +corrected_mz << "\n";
    std::cout << "Expected readout: " << +expected_mz << "\n";
    if (corrected_mz != expected_mz)
      numLerrors++;
    std::cout << "\n";
  }

  std::cout << "numLogicalErrors: " << numLerrors << "\n";
}

Compile and run with

nvq++ --enable-mlir -lcudaq-qec circuit_level_noise.cpp -o circuit_level_noise
./circuit_level_noise
  1. QEC Code and Decoder types:

    • As in the code capacity example, our central objects are the qec.code and qec.decoder types.

  2. Clifford simulation backend:

    • As the size of QEC circuits can grow quite large, Clifford simulation is often the best tool for these simulations.

    • cudaq.set_target("stim") selects the highly performant Stim simulator as the simulation backend.

  3. Noise model:

    • To add noisy gates we use the cudaq.NoiseModel type.

    • CUDA-Q supports the generation of arbitrary noise channels, but here we use a qec.TwoQubitDepolarization channel to add a depolarization channel.

    • This is added to the CX gate by adding it to the X gate with 1 control.

    • This noisy gate is added to every qubit via that noise.add_all_qubit_channel function.

  4. Getting circuit-level noisy data:

    • The qec.code is the first input parameter here, as the code’s stabilizer_round determines the circuits executed.

    • Each memory circuit runs for an input number of nRounds, which specifies how many stabilizer_round kernels are ran.

    • After nRounds the data qubits are measured and the run is over.

    • This is performed nShots number of times.

    • During a shot, each syndrome is xor’d against the preceding syndrome, so that we can track a sparser flow of data showing which round each parity check was violated.

    • This means we will get a total of nShots * (nRounds - 1) syndromes to decode and analyze.

  5. Data qubit measurements:

    • The data qubits are only read out after the end of each shot, so there are nShots worth of data readouts.

    • The basis of the data qubit measurements depends on the state preparation used.

    • Z-basis readout when preparing the logical |0> or logical |1> state with the qec.operation.prep0 or qec.operation.prep1 kernels.

    • X-basis readout when preparing the logical |+> or logical |-> state with the qec.operation.prepp or qec.operation.prepm kernels.

  6. Logical Errors:

    • From here, the decoding procedure is again similar to the code capacity case, expect for we use a pauli frame to track errors that happen each QEC cycle.

    • The final values of the pauli frame tell us how our logical state flipped during the experiment, and what needs to be done to correct it.

    • We compare our known initial state (corrected by the Pauli frame), against our measured data qubits to determine if a logical error occurred.

The CUDA-Q QEC library thus provides a platform for numerical QEC experiments. The qec.code can be used to analyze a variety of QEC codes (both library or user provided), with a variety of decoders (both library or user provided). The CUDA-Q QEC library also provides tools to speed up the automation of generating noisy data and syndromes.