Draw anticommutations in detslice diagrams instead of failing to make…

… the diagram (#736)

- Give better error messages when SparseRevTracker hits an
anticommutation

Fixes #471

src/stim/diagram/detector_slice/detector_slice_set.cc

@@ -83,7 +83,7 @@ bool DetectorSliceSetComputer::process_op_rev(const Circuit &parent, const Circu
 }
 DetectorSliceSetComputer::DetectorSliceSetComputer(
     const Circuit &circuit, uint64_t first_yield_tick, uint64_t num_yield_ticks)
-    : tracker(circuit.count_qubits(), circuit.count_measurements(), circuit.count_detectors()),
+    : tracker(circuit.count_qubits(), circuit.count_measurements(), circuit.count_detectors(), false),
       first_yield_tick(first_yield_tick),
       num_yield_ticks(num_yield_ticks) {
     tick_cur = circuit.count_ticks() + 1;  // +1 because of artificial TICKs at start and end.
@@ -111,6 +111,25 @@ void DetectorSliceSet::write_text_diagram_to(std::ostream &out) const {
     DiagramTimelineAsciiDrawer drawer(num_qubits, false);
     drawer.moment_spacing = 2;
 
+    for (const auto &s : anticommutations) {
+        drawer.reserve_drawing_room_for_targets(s.second);
+        for (const auto &t : s.second) {
+            std::stringstream ss;
+            ss << "ANTICOMMUTED";
+            ss << ":";
+            ss << s.first.second;
+            drawer.diagram.add_entry(AsciiDiagramEntry{
+                AsciiDiagramPos{
+                    drawer.m2x(drawer.cur_moment + 1),
+                    drawer.q2y(t.qubit_value()),
+                    0,
+                    0.5,
+                },
+                ss.str(),
+            });
+        }
+    }
+
     for (const auto &s : slices) {
         drawer.reserve_drawing_room_for_targets(s.second);
         for (const auto &t : s.second) {
@@ -230,6 +249,23 @@ DetectorSliceSet DetectorSliceSet::from_circuit_ticks(
     result.num_ticks = num_ticks;
 
     helper.on_tick_callback = [&]() {
+        // Process anticommutations.
+        for (const auto &[d, g] : helper.tracker.anticommutations) {
+            result.anticommutations[{helper.tick_cur, d}].push_back(g);
+
+            // Stop propagating it backwards if it broke.
+            for (size_t q = 0; q < num_qubits; q++) {
+                if (helper.tracker.xs[q].contains(d)) {
+                    helper.tracker.xs[q].xor_item(d);
+                }
+                if (helper.tracker.zs[q].contains(d)) {
+                    helper.tracker.zs[q].xor_item(d);
+                }
+            }
+        }
+        helper.tracker.anticommutations.clear();
+
+        // Record locations of detectors and observables.
         for (size_t q = 0; q < num_qubits; q++) {
             xs.clear();
             ys.clear();
@@ -861,10 +897,13 @@ void DetectorSliceSet::write_svg_contents_to(
 
     bool haveDrawnCorners = false;
 
-    using tup = std::tuple<uint64_t, stim::DemTarget, SpanRef<const GateTarget>>;
+    using tup = std::tuple<uint64_t, stim::DemTarget, SpanRef<const GateTarget>, bool>;
     std::vector<tup> sorted_terms;
     for (const auto &e : slices) {
-        sorted_terms.push_back({e.first.first, e.first.second, e.second});
+        sorted_terms.push_back({e.first.first, e.first.second, e.second, false});
+    }
+    for (const auto &e : anticommutations) {
+        sorted_terms.push_back({e.first.first, e.first.second, e.second, true});
     }
     std::stable_sort(sorted_terms.begin(), sorted_terms.end(), [](const tup &e1, const tup &e2) -> int {
         int a = (int)std::get<2>(e1).size();
@@ -877,6 +916,22 @@ void DetectorSliceSet::write_svg_contents_to(
         uint64_t tick = std::get<0>(e);
         DemTarget target = std::get<1>(e);
         SpanRef<const GateTarget> terms = std::get<2>(e);
+        bool is_anticommutation = std::get<3>(e);
+        if (is_anticommutation) {
+            for (const auto &anti_target : terms) {
+                auto c = coords(tick + 1, anti_target.qubit_value());
+                out << R"SVG(<circle)SVG";
+                write_key_val(out, "cx", c.xyz[0]);
+                write_key_val(out, "cy", c.xyz[1]);
+                write_key_val(out, "r", scale);
+                write_key_val(out, "fill", "none");
+                write_key_val(out, "stroke", "magenta");
+                write_key_val(out, "stroke-width", scale / 2);
+                out << "/>\n";
+            }
+            continue;
+        }
+
         if (target.is_observable_id()) {
             _draw_observable(
                 out, target.val(), unscaled_coords, coords, tick, terms, pts_workspace, tick < end_tick - 1, scale);

src/stim/diagram/detector_slice/detector_slice_set.h

@@ -45,6 +45,8 @@ struct DetectorSliceSet {
     std::map<uint64_t, std::vector<double>> detector_coordinates;
     /// (tick, DemTarget) -> terms in the slice
     std::map<std::pair<uint64_t, stim::DemTarget>, std::vector<stim::GateTarget>> slices;
+    /// (tick, DemTarget) -> anticommutations in the slice
+    std::map<std::pair<uint64_t, stim::DemTarget>, std::vector<stim::GateTarget>> anticommutations;
 
     /// Args:
     ///     circuit: The circuit to make a detector slice diagram from.

src/stim/diagram/detector_slice/detector_slice_set.test.cc

@@ -398,3 +398,54 @@ TEST(inv_space_fill_transform, inv_space_fill_transform) {
     ASSERT_EQ(inv_space_fill_transform({0, 0}), 0);
     ASSERT_EQ(inv_space_fill_transform({4, 55.5}), 339946);
 }
+
+TEST(detector_slice_set, from_circuit_with_errors) {
+    CoordFilter empty_filter;
+    auto slice_set = DetectorSliceSet::from_circuit_ticks(
+        stim::Circuit(R"CIRCUIT(
+            TICK
+            R 0
+            TICK
+            R 0
+            TICK
+            MXX 0 1
+            DETECTOR rec[-1]
+        )CIRCUIT"),
+        0,
+        5,
+        {&empty_filter});
+    ASSERT_EQ(
+        slice_set.anticommutations,
+        (std::map<std::pair<uint64_t, stim::DemTarget>, std::vector<stim::GateTarget>>{
+            {{2, DemTarget::relative_detector_id(0)}, {GateTarget::x(0)}},
+        }));
+    ASSERT_EQ(
+        slice_set.slices,
+        (std::map<std::pair<uint64_t, stim::DemTarget>, std::vector<stim::GateTarget>>{
+            {{3, DemTarget::relative_detector_id(0)}, {GateTarget::x(0), GateTarget::x(1)}},
+        }));
+}
+
+TEST(circuit_diagram_timeline_text, anticommuting_detector_circuit) {
+    auto circuit = Circuit(R"CIRCUIT(
+        TICK
+        R 0
+        TICK
+        R 0
+        TICK
+        MXX 0 1
+        M 2
+        DETECTOR rec[-1]
+        DETECTOR rec[-2]
+    )CIRCUIT");
+    CoordFilter empty_filter;
+    std::stringstream ss;
+    ss << DetectorSliceSet::from_circuit_ticks(circuit, 0, 10, &empty_filter);
+    ASSERT_EQ("\n" + ss.str() + "\n", R"DIAGRAM(
+q0: -------ANTICOMMUTED:D1--X:D1-
+                            |
+q1: ------------------------X:D1-
+
+q2: -Z:D0--Z:D0-------------Z:D0-
+)DIAGRAM");
+}

src/stim/diagram/timeline/timeline_svg_drawer.test.cc

@@ -362,3 +362,22 @@ TEST(diagram_timeline_svg_drawer, test_circuit_all_ops_detslice) {
         circuit, ss, 0, UINT64_MAX, DiagramTimelineSvgDrawerMode::SVG_MODE_TIME_DETECTOR_SLICE, {&empty_filter});
     expect_string_is_identical_to_saved_file(ss.str(), "circuit_all_ops_detslice.svg");
 }
+
+TEST(diagram_timeline_svg_drawer, anticommuting_detector_circuit) {
+    CoordFilter empty_filter;
+    std::stringstream ss;
+    auto circuit = Circuit(R"CIRCUIT(
+        TICK
+        R 0
+        TICK
+        R 0
+        TICK
+        MXX 0 1
+        M 2
+        DETECTOR rec[-1]
+        DETECTOR rec[-2]
+    )CIRCUIT");
+    DiagramTimelineSvgDrawer::make_diagram_write_to(
+        circuit, ss, 0, UINT64_MAX, DiagramTimelineSvgDrawerMode::SVG_MODE_TIME_DETECTOR_SLICE, {&empty_filter});
+    expect_string_is_identical_to_saved_file(ss.str(), "anticommuting_detslice.svg");
+}

src/stim/simulators/sparse_rev_frame_tracker.cc

@@ -197,31 +197,48 @@ SparseUnsignedRevFrameTracker::SparseUnsignedRevFrameTracker(
       anticommutations() {
 }
 
+void SparseUnsignedRevFrameTracker::fail_due_to_anticommutation(const CircuitInstruction &inst) {
+    std::stringstream ss;
+    ss << "While running backwards through the circuit, during reverse-execution of the instruction\n";
+    ss << "    " << inst << "\n";
+    ss << "the following detecting region vs dissipation anticommutations occurred\n";
+    for (auto &[d, g] : anticommutations) {
+        ss << "    " << d << " vs " << g << "\n";
+    }
+    ss << "Therefore invalid detectors/observables are present in the circuit.\n";
+    throw std::invalid_argument(ss.str());
+}
+
 void SparseUnsignedRevFrameTracker::handle_xor_gauge(
-    SpanRef<const DemTarget> sorted1, SpanRef<const DemTarget> sorted2) {
+    SpanRef<const DemTarget> sorted1,
+    SpanRef<const DemTarget> sorted2,
+    const CircuitInstruction &inst,
+    GateTarget location) {
     if (sorted1 == sorted2) {
         return;
     }
-    if (fail_on_anticommute) {
-        throw std::invalid_argument("A detector or observable anticommuted with a dissipative operation.");
-    }
     SparseXorVec<DemTarget> dif;
     dif.xor_sorted_items(sorted1);
     dif.xor_sorted_items(sorted2);
     for (const auto &d : dif) {
-        anticommutations.insert(d);
+        anticommutations.insert({d, location});
+    }
+    if (fail_on_anticommute) {
+        fail_due_to_anticommutation(inst);
     }
 }
 
-void SparseUnsignedRevFrameTracker::handle_gauge(SpanRef<const DemTarget> sorted) {
+void SparseUnsignedRevFrameTracker::handle_gauge(SpanRef<const DemTarget> sorted,
+    const CircuitInstruction &inst,
+    GateTarget location) {
     if (sorted.empty()) {
         return;
     }
-    if (fail_on_anticommute) {
-        throw std::invalid_argument("A detector or observable anticommuted with a dissipative operation.");
-    }
     for (const auto &d : sorted) {
-        anticommutations.insert(d);
+        anticommutations.insert({d, location});
+    }
+    if (fail_on_anticommute) {
+        fail_due_to_anticommutation(inst);
     }
 }
 
@@ -302,19 +319,19 @@ void SparseUnsignedRevFrameTracker::undo_ZCZ_single(GateTarget c, GateTarget t)
 void SparseUnsignedRevFrameTracker::handle_x_gauges(const CircuitInstruction &dat) {
     for (size_t k = dat.targets.size(); k-- > 0;) {
         auto q = dat.targets[k].qubit_value();
-        handle_gauge(xs[q].range());
+        handle_gauge(xs[q].range(), dat, GateTarget::x(q));
     }
 }
 void SparseUnsignedRevFrameTracker::handle_y_gauges(const CircuitInstruction &dat) {
     for (size_t k = dat.targets.size(); k-- > 0;) {
         auto q = dat.targets[k].qubit_value();
-        handle_xor_gauge(xs[q].range(), zs[q].range());
+        handle_xor_gauge(xs[q].range(), zs[q].range(), dat, GateTarget::y(q));
     }
 }
 void SparseUnsignedRevFrameTracker::handle_z_gauges(const CircuitInstruction &dat) {
     for (size_t k = dat.targets.size(); k-- > 0;) {
         auto q = dat.targets[k].qubit_value();
-        handle_gauge(zs[q].range());
+        handle_gauge(zs[q].range(), dat, GateTarget::z(q));
     }
 }
 void SparseUnsignedRevFrameTracker::undo_MPP(const CircuitInstruction &target_data) {

src/stim/simulators/sparse_rev_frame_tracker.h

@@ -41,7 +41,7 @@ struct SparseUnsignedRevFrameTracker {
     /// If true, an exception is raised if anticommutation is detected.
     bool fail_on_anticommute;
     /// Where anticommuting dets and obs are stored.
-    std::set<DemTarget> anticommutations;
+    std::set<std::pair<DemTarget, GateTarget>> anticommutations;
 
     SparseUnsignedRevFrameTracker(
         uint64_t num_qubits,
@@ -62,8 +62,9 @@ struct SparseUnsignedRevFrameTracker {
     void undo_gate(const CircuitInstruction &inst);
     void undo_gate(const CircuitInstruction &op, const Circuit &parent);
 
-    void handle_xor_gauge(SpanRef<const DemTarget> sorted1, SpanRef<const DemTarget> sorted2);
-    void handle_gauge(SpanRef<const DemTarget> sorted);
+    void handle_xor_gauge(SpanRef<const DemTarget> sorted1, SpanRef<const DemTarget> sorted2, const CircuitInstruction &inst, GateTarget location);
+    void handle_gauge(SpanRef<const DemTarget> sorted, const CircuitInstruction &inst, GateTarget location);
+    void fail_due_to_anticommutation(const CircuitInstruction &inst);
     void undo_classical_pauli(GateTarget classical_control, GateTarget target);
     void undo_ZCX_single(GateTarget c, GateTarget t);
     void undo_ZCY_single(GateTarget c, GateTarget t);

src/stim/simulators/sparse_rev_frame_tracker.test.cc

@@ -706,8 +706,11 @@ TEST(SparseUnsignedRevFrameTracker, tracks_anticommutation) {
     rev.undo_circuit(circuit);
     ASSERT_EQ(
         rev.anticommutations,
-        (std::set<DemTarget>{
-            DemTarget::relative_detector_id(0), DemTarget::relative_detector_id(2), DemTarget::observable_id(2)}));
+        (std::set<std::pair<DemTarget, GateTarget>>{
+            {DemTarget::relative_detector_id(0), GateTarget::x(2)},
+            {DemTarget::relative_detector_id(2), GateTarget::x(2)},
+            {DemTarget::observable_id(2), GateTarget::x(1)},
+            {DemTarget::observable_id(2), GateTarget::x(2)}}));
 
     SparseUnsignedRevFrameTracker rev2(circuit.count_qubits(), circuit.count_measurements(), circuit.count_detectors());
     ASSERT_THROW({ rev.undo_circuit(circuit); }, std::invalid_argument);
@@ -730,3 +733,15 @@ TEST(SparseUnsignedRevFrameTracker, MZZ) {
     ASSERT_TRUE(rev.zs[1].empty());
     ASSERT_EQ(rev.zs[2].sorted_items, (std::vector<DemTarget>{DemTarget::relative_detector_id(0)}));
 }
+
+TEST(SparseUnsignedRevFrameTracker, fail_anticommute) {
+    Circuit circuit(R"CIRCUIT(
+        RX 0 1 2
+        M 2
+        DETECTOR rec[-1]
+    )CIRCUIT");
+
+    SparseUnsignedRevFrameTracker rev(
+        circuit.count_qubits(), circuit.count_measurements(), circuit.count_detectors(), true);
+    ASSERT_THROW({ rev.undo_circuit(circuit); }, std::invalid_argument);
+}

src/stim/stabilizers/flow.inl

@@ -341,8 +341,8 @@ std::vector<bool> check_if_circuit_has_unsigned_stabilizer_flows(
             result[t.val()] = false;
         }
     }
-    for (const auto &anti : rev.anticommutations) {
-        result[anti.val()] = false;
+    for (const auto &[dem_target, gate_target] : rev.anticommutations) {
+        result[dem_target.val()] = false;
     }
 
     return result;