Improve handling of disjoint errors in circuit-to-dem conversion (#617)

- Analytically solve the bulk of the PAULI_CHANNEL_1 disjoint-to-independent problem, with a fallback to Newton-Rhapson approximation, with a fallback to just leaving it disjoint
- Merge undistinguished cases of a disjoint error channel using addition instead of Bernoulli addition
- Add `HERALDED_PAULI_CHANNEL_1` gate
- Pull `dev/regen_crumble_html.sh` script out of cpp string generation script

Fixes #602

dev/canvas_with_texture_for_3d_diagrams.html

@@ -97,13 +97,25 @@
     let {x, y} = pickRect(i);
     ctx.fillStyle = 'black';
     ctx.fillRect(x * 32, y * 32, 32, 32);
-    ctx.fillStyle = 'white';
+    ctx.fillStyle = 'yellow';
     ctx.font = '6pt serif';
     ctx.fillText('HERALD', x * 32 + 16 - ctx.measureText('HERALD').width / 2, y * 32 + 11);
     ctx.font = '7pt serif';
     ctx.fillText('ERASE', x * 32 + 16 - ctx.measureText('ERASE').width / 2, y * 32 + 18);
 }
 
+function drawHeraldPauliError1(ctx, i) {
+    let {x, y} = pickRect(i);
+    ctx.fillStyle = 'black';
+    ctx.fillRect(x * 32, y * 32, 32, 32);
+    ctx.fillStyle = 'yellow';
+    ctx.font = '6pt serif';
+    ctx.fillText('HERALD', x * 32 + 16 - ctx.measureText('HERALD').width / 2, y * 32 + 3);
+    ctx.font = '9pt serif';
+    ctx.fillText('Pauli', x * 32 + 16 - ctx.measureText('Pauli').width / 2, y * 32 + 11);
+    ctx.fillText('Err1', x * 32 + 16 - ctx.measureText('Err1').width / 2, y * 32 + 21);
+}
+
 /**
  * @param {!CanvasRenderingContext2D} ctx
  */
@@ -190,6 +202,7 @@
     drawRect(ctx, "M_ZZ", "#44F", "#000", n++);
     drawRect(ctx, "M_PAD", "#888", "#000", n++);
     drawHeraldErase(ctx, n++);
+    drawHeraldPauliError1(ctx, n++);
 }
 
 draw(document.getElementById('cv').getContext('2d'))

dev/regen_crumble_html.sh

@@ -0,0 +1,13 @@
+#!/bin/bash
+set -e
+
+cd "$( dirname "${BASH_SOURCE[0]}" )"
+cd "$(git rev-parse --show-toplevel)"
+
+cat glue/crumble/crumble.html | grep -v "^<script";
+echo "<script>";
+# HACK: this temp file is to work around https://github.com/rollup/rollup/issues/5097
+rollup glue/crumble/main.js > tmp_crumble.tmp
+uglifyjs -c -m --mangle-props --toplevel < tmp_crumble.tmp;
+rm tmp_crumble.tmp
+echo "</script>";

dev/regen_crumble_to_cpp_string_write_to_stdout.sh

@@ -5,25 +5,15 @@ set -e
 cd "$( dirname "${BASH_SOURCE[0]}" )"
 cd "$(git rev-parse --show-toplevel)"
 
-{
-    echo '#include "stim/diagram/crumble_data.h"';
-    echo '';
-    echo 'std::string stim_draw_internal::make_crumble_html() {';
-    echo '    std::string result;'
-    {
-        cat glue/crumble/crumble.html | grep -v "^<script";
-        echo "<script>";
-        # HACK: this temp file is to work around https://github.com/rollup/rollup/issues/5097
-        rollup glue/crumble/main.js > tmp_crumble.tmp
-        uglifyjs -c -m --mangle-props --toplevel < tmp_crumble.tmp;
-        rm tmp_crumble.tmp
-        echo "</script>";
-    } | python -c '
+echo '#include "stim/diagram/crumble_data.h"'
+echo '';
+echo 'std::string stim_draw_internal::make_crumble_html() {'
+echo '    std::string result;'
+dev/regen_crumble_html.sh | python -c '
 import sys
 for line in sys.stdin:
     for k in range(0, len(line), 1024):
         part = line[k:k+1024]
         print(f"""    result.append(R"CRUMBLE_PART({part})CRUMBLE_PART");""")';
-    echo '    return result;'
-    echo '}';
-}
+echo '    return result;'
+echo '}'

doc/gates.md

@@ -52,6 +52,7 @@
     - [E](#E)
     - [ELSE_CORRELATED_ERROR](#ELSE_CORRELATED_ERROR)
     - [HERALDED_ERASE](#HERALDED_ERASE)
+    - [HERALDED_PAULI_CHANNEL_1](#HERALDED_PAULI_CHANNEL_1)
     - [PAULI_CHANNEL_1](#PAULI_CHANNEL_1)
     - [PAULI_CHANNEL_2](#PAULI_CHANNEL_2)
     - [X_ERROR](#X_ERROR)
@@ -1911,6 +1912,21 @@ record. When it doesn't fire, nothing happens to the target qubit and a
 is erased to the maximally mixed state by applying X_ERROR(0.5) and
 Z_ERROR(0.5).
 
+CAUTION: when converting a circuit with this error into a detector
+error model, this channel is split into multiple potential effects.
+In the context of a DEM, these effects are considered independent.
+This is an approximation, because independent effects can be combined.
+The effect of this approximation, assuming a detector is declared
+on the herald, is that it appears this detector can be cancelled out
+by two of the (originally disjoint) heralded effects firing together.
+Sampling from the DEM instead of the circuit can thus produce unheralded
+errors, even if the circuit noise model only contains heralded errors.
+These issues occur with probability p^2, where p is the probability of a
+heralded error, since two effects that came from the same heralded error
+must occur together to cancel out the herald detector. This also means
+a decoder configured using the DEM will think there's a chance of unheralded
+errors even if the circuit the DEM came from only uses heralded errors.
+
 Parens Arguments:
 
     A single float (p) specifying the chance of the noise firing.
@@ -1954,6 +1970,75 @@ Examples:
     DETECTOR rec[-4]
     DETECTOR rec[-5]
 
+<a name="HERALDED_PAULI_CHANNEL_1"></a>
+### The 'HERALDED_PAULI_CHANNEL_1' Instruction
+
+A heralded error channel that applies biased noise.
+
+This error records a bit into the measurement record, indicating whether
+or not the herald fired. How likely it is that the herald fires, and the
+corresponding chance of each possible error effect (I, X, Y, or Z) are
+configured by the parens arguments of the instruction.
+
+CAUTION: when converting a circuit with this error into a detector
+error model, this channel is split into multiple potential effects.
+In the context of a DEM, these effects are considered independent.
+This is an approximation, because independent effects can be combined.
+The effect of this approximation, assuming a detector is declared
+on the herald, is that it appears this detector can be cancelled out
+by two of the (originally disjoint) heralded effects firing together.
+Sampling from the DEM instead of the circuit can thus produce unheralded
+errors, even if the circuit noise model only contains heralded errors.
+These issues occur with probability p^2, where p is the probability of a
+heralded error, since two effects that came from the same heralded error
+must occur together to cancel out the herald detector. This also means
+a decoder configured using the DEM will think there's a chance of unheralded
+errors even if the circuit the DEM came from only uses heralded errors.
+
+Parens Arguments:
+
+    This instruction takes four arguments (pi, px, py, pz). The
+    arguments are disjoint probabilities, specifying the chances
+    of heralding with various effects.
+
+    pi is the chance of heralding with no effect (a false positive).
+    px is the chance of heralding with an X error.
+    py is the chance of heralding with a Y error.
+    pz is the chance of heralding with a Z error.
+
+Targets:
+
+    Qubits to apply heralded biased noise to.
+
+Pauli Mixture:
+
+    1-pi-px-py-pz: record 0, apply I
+               pi: record 1, apply I
+               px: record 1, apply X
+               py: record 1, apply Y
+               pz: record 1, apply Z
+
+Examples:
+
+    # With 10% probability perform a phase flip of qubit 0.
+    HERALDED_PAULI_CHANNEL_1(0, 0, 0, 0.1) 0
+    DETECTOR rec[-1]  # Include the herald in detectors available to the decoder
+
+    # With 20% probability perform a heralded dephasing of qubit 0.
+    HERALDED_PAULI_CHANNEL_1(0.1, 0, 0, 0.1) 0
+    DETECTOR rec[-1]
+
+    # Subject a Bell Pair to heralded noise.
+    MXX 0 1
+    MZZ 0 1
+    HERALDED_PAULI_CHANNEL_1(0.01, 0.02, 0.03, 0.04) 0 1
+    MXX 0 1
+    MZZ 0 1
+    DETECTOR rec[-1] rec[-5]  # Did ZZ stabilizer change?
+    DETECTOR rec[-2] rec[-6]  # Did XX stabilizer change?
+    DETECTOR rec[-3]    # Did the herald on qubit 1 fire?
+    DETECTOR rec[-4]    # Did the herald on qubit 0 fire?
+
 <a name="PAULI_CHANNEL_1"></a>
 ### The 'PAULI_CHANNEL_1' Instruction
 
@@ -1962,9 +2047,9 @@ A single qubit Pauli error channel with explicitly specified probabilities for e
 Parens Arguments:
 
     Three floats specifying disjoint Pauli case probabilities.
-    px: Probability of applying an X operation.
-    py: Probability of applying a Y operation.
-    pz: Probability of applying a Z operation.
+    px: Disjoint probability of applying an X error.
+    py: Disjoint probability of applying a Y error.
+    pz: Disjoint probability of applying a Z error.
 
 Targets:
 

doc/python_api_reference_vDev.md

@@ -12510,8 +12510,8 @@ def gate_data(
         >>> stim.gate_data('cnot').is_two_qubit_gate
         True
         >>> gate_dict = stim.gate_data()
-        >>> len(gate_dict)
-        65
+        >>> len(gate_dict) > 50
+        True
         >>> gate_dict['MX'].produces_measurements
         True
     """

doc/stim.pyi

@@ -9741,8 +9741,8 @@ def gate_data(
         >>> stim.gate_data('cnot').is_two_qubit_gate
         True
         >>> gate_dict = stim.gate_data()
-        >>> len(gate_dict)
-        65
+        >>> len(gate_dict) > 50
+        True
         >>> gate_dict['MX'].produces_measurements
         True
     """

file_lists/benchmark_files

@@ -14,6 +14,7 @@ src/stim/simulators/dem_sampler.perf.cc
 src/stim/simulators/error_analyzer.perf.cc
 src/stim/simulators/frame_simulator.perf.cc
 src/stim/simulators/tableau_simulator.perf.cc
+src/stim/stabilizers/conversions.perf.cc
 src/stim/stabilizers/pauli_string.perf.cc
 src/stim/stabilizers/tableau.perf.cc
 src/stim/stabilizers/tableau_iter.perf.cc

file_lists/source_files_no_main

@@ -82,3 +82,4 @@ src/stim/simulators/matched_error.cc
 src/stim/simulators/sparse_rev_frame_tracker.cc
 src/stim/simulators/transform_without_feedback.cc
 src/stim/simulators/vector_simulator.cc
+src/stim/stabilizers/conversions.cc

file_lists/test_files

@@ -17,6 +17,7 @@ src/stim/cmd/command_sample_dem.test.cc
 src/stim/dem/detector_error_model.test.cc
 src/stim/diagram/ascii_diagram.test.cc
 src/stim/diagram/base64.test.cc
+src/stim/diagram/coord.test.cc
 src/stim/diagram/detector_slice/detector_slice_set.test.cc
 src/stim/diagram/graph/match_graph_3d_drawer.test.cc
 src/stim/diagram/graph/match_graph_svg_drawer.test.cc

glue/cirq/stimcirq/_stim_to_cirq.py

@@ -485,6 +485,9 @@ def handler(
             "HERALDED_ERASE": not_impl(
                 "Converting HERALDED_ERASE to cirq is not supported."
             ),
+            "HERALDED_PAULI_CHANNEL_1": not_impl(
+                "Converting HERALDED_PAULI_CHANNEL_1 to cirq is not supported."
+            ),
         }
 
 

glue/python/src/stim/__init__.pyi

@@ -9741,8 +9741,8 @@ def gate_data(
         >>> stim.gate_data('cnot').is_two_qubit_gate
         True
         >>> gate_dict = stim.gate_data()
-        >>> len(gate_dict)
-        65
+        >>> len(gate_dict) > 50
+        True
         >>> gate_dict['MX'].produces_measurements
         True
     """

src/stim/circuit/circuit.cc

@@ -764,7 +764,7 @@ const Circuit Circuit::aliased_noiseless_circuit() const {
     for (const auto &op : operations) {
         auto flags = GATE_DATA.items[op.gate_type].flags;
         if (flags & GATE_PRODUCES_RESULTS) {
-            if (op.gate_type == GateType::HERALDED_ERASE) {
+            if (op.gate_type == GateType::HERALDED_ERASE || op.gate_type == GateType::HERALDED_PAULI_CHANNEL_1) {
                 // Replace heralded errors with fixed MPAD.
                 result.target_buf.ensure_available(op.targets.size());
                 auto &tail = result.target_buf.tail;
@@ -796,7 +796,7 @@ Circuit Circuit::without_noise() const {
     for (const auto &op : operations) {
         auto flags = GATE_DATA.items[op.gate_type].flags;
         if (flags & GATE_PRODUCES_RESULTS) {
-            if (op.gate_type == GateType::HERALDED_ERASE) {
+            if (op.gate_type == GateType::HERALDED_ERASE || op.gate_type == GateType::HERALDED_PAULI_CHANNEL_1) {
                 // Replace heralded errors with fixed MPAD.
                 result.target_buf.ensure_available(op.targets.size());
                 auto &tail = result.target_buf.tail;

src/stim/circuit/circuit.test.cc

@@ -1647,6 +1647,7 @@ Circuit stim::generate_test_circuit_with_all_operations() {
         Y_ERROR(0.02) 1
         Z_ERROR(0.03) 2
         HERALDED_ERASE(0.04) 3
+        HERALDED_PAULI_CHANNEL_1(0.01, 0.02, 0.03, 0.04) 6
         TICK
 
         # Collapsing Gates

src/stim/circuit/gate_data.h

@@ -73,7 +73,7 @@ constexpr inline uint16_t gate_name_to_hash(const char *v, size_t n) {
         auto c3 = v[3] | 0x20;
         auto c5 = v[5] | 0x20;
         result ^= c3 * 61;
-        result += c5 * 27;
+        result += c5 * 77;
     }
     return result & 0x1FF;
 }
@@ -82,7 +82,7 @@ constexpr inline uint16_t gate_name_to_hash(const char *c) {
     return gate_name_to_hash(c, std::char_traits<char>::length(c));
 }
 
-constexpr const size_t NUM_DEFINED_GATES = 66;
+constexpr const size_t NUM_DEFINED_GATES = 67;
 
 enum GateType : uint8_t {
     NOT_A_GATE = 0,
@@ -132,6 +132,7 @@ enum GateType : uint8_t {
     ELSE_CORRELATED_ERROR,
     // Heralded noise channels
     HERALDED_ERASE,
+    HERALDED_PAULI_CHANNEL_1,
     // Pauli gates
     I,
     X,

src/stim/circuit/gate_data.pybind.cc

@@ -80,8 +80,8 @@ void stim_pybind::pybind_gate_data_methods(pybind11::module &m, pybind11::class_
                 >>> stim.gate_data('cnot').is_two_qubit_gate
                 True
                 >>> gate_dict = stim.gate_data()
-                >>> len(gate_dict)
-                65
+                >>> len(gate_dict) > 50
+                True
                 >>> gate_dict['MX'].produces_measurements
                 True
         )DOC")