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

Author: Strilanc

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."
+            ),
         }