Fix stim.Circuit.to_tableau ignoring SPP and SPP_DAG (#847)

- Fix `SPP` not being classified as a unitary gate
- Fix `SPP_DAG` not being classified as a unitary gate
- Fix some code/tests assuming a gate being unitary implies it has a
non-target-dependent unitary matrix

Fixes #846

doc/gates.md

@@ -3043,7 +3043,7 @@ Decomposition (into H, S, CX, M, R):
 
 
 <a name="SPP"></a>
-### The 'SPP' Instruction
+### The 'SPP' Gate
 
 The generalized S gate. Phases the -1 eigenspace of Pauli product observables by i.
 
@@ -3109,7 +3109,7 @@ Decomposition (into H, S, CX, M, R):
 
 
 <a name="SPP_DAG"></a>
-### The 'SPP_DAG' Instruction
+### The 'SPP_DAG' Gate
 
 The generalized S_DAG gate. Phases the -1 eigenspace of Pauli product observables by -i.
 

src/stim/circuit/circuit_pybind_test.py

@@ -864,7 +864,6 @@ def test_likeliest_error_sat_problem():
         DETECTOR rec[-1] rec[-2]
     """)
     sat_str = c.likeliest_error_sat_problem(quantization=100)
-    print(sat_str)
     assert sat_str == 'p wcnf 2 4 401\n18 -1 0\n100 -2 0\n401 -1 0\n401 2 0\n'
 
 
@@ -1902,3 +1901,14 @@ def test_pop():
 def test_circuit_create_with_odd_cx():
     with pytest.raises(ValueError, match="0, 1, 2"):
         stim.Circuit("CX 0 1 2")
+
+
+def test_to_tableau():
+    assert stim.Circuit().to_tableau() == stim.Tableau(0)
+    assert stim.Circuit("QUBIT_COORDS 0").to_tableau() == stim.Tableau(1)
+    assert stim.Circuit("I 0").to_tableau() == stim.Tableau(1)
+    assert stim.Circuit("H 0").to_tableau() == stim.Tableau.from_named_gate("H")
+    assert stim.Circuit("CX 0 1").to_tableau() == stim.Tableau.from_named_gate("CX")
+    assert stim.Circuit("SPP Z0").to_tableau() == stim.Tableau.from_named_gate("S")
+    assert stim.Circuit("SPP X0").to_tableau() == stim.Tableau.from_named_gate("SQRT_X")
+    assert stim.Circuit("SPP_DAG Y0*Y1").to_tableau() == stim.Tableau.from_named_gate("SQRT_YY_DAG")

src/stim/cmd/command_help.cc

@@ -290,7 +290,7 @@ void print_stabilizer_generators(Acc &out, const Gate &gate) {
 }
 
 void print_bloch_vector(Acc &out, const Gate &gate) {
-    if (!(gate.flags & GATE_IS_UNITARY) || (gate.flags & GATE_TARGETS_PAIRS)) {
+    if (!(gate.flags & GATE_IS_UNITARY) || !(gate.flags & GATE_IS_SINGLE_QUBIT_GATE)) {
         return;
     }
 
@@ -343,7 +343,7 @@ void print_bloch_vector(Acc &out, const Gate &gate) {
 }
 
 void print_unitary_matrix(Acc &out, const Gate &gate) {
-    if (!(gate.flags & GATE_IS_UNITARY)) {
+    if (!gate.has_known_unitary_matrix()) {
         return;
     }
     auto matrix = gate.unitary();

src/stim/gates/gate_data_pauli_product.cc

@@ -99,7 +99,7 @@ S 1
             .id = GateType::SPP,
             .best_candidate_inverse_id = GateType::SPP_DAG,
             .arg_count = 0,
-            .flags = (GateFlags)(GATE_TARGETS_PAULI_STRING | GATE_TARGETS_COMBINERS),
+            .flags = (GateFlags)(GATE_TARGETS_PAULI_STRING | GATE_TARGETS_COMBINERS | GATE_IS_UNITARY),
             .category = "P_Generalized Pauli Product Gates",
             .help = R"MARKDOWN(
 The generalized S gate. Phases the -1 eigenspace of Pauli product observables by i.
@@ -174,7 +174,7 @@ CX 2 1
             .id = GateType::SPP_DAG,
             .best_candidate_inverse_id = GateType::SPP,
             .arg_count = 0,
-            .flags = (GateFlags)(GATE_TARGETS_PAULI_STRING | GATE_TARGETS_COMBINERS),
+            .flags = (GateFlags)(GATE_TARGETS_PAULI_STRING | GATE_TARGETS_COMBINERS | GATE_IS_UNITARY),
             .category = "P_Generalized Pauli Product Gates",
             .help = R"MARKDOWN(
 The generalized S_DAG gate. Phases the -1 eigenspace of Pauli product observables by -i.

src/stim/gates/gates.cc

@@ -267,6 +267,10 @@ std::array<float, 4> Gate::to_axis_angle() const {
     return {rx, ry, rz, acosf(rs) * 2};
 }
 
+bool Gate::has_known_unitary_matrix() const {
+    return (flags & GateFlags::GATE_IS_UNITARY) && (flags & (GateFlags::GATE_IS_SINGLE_QUBIT_GATE | GateFlags::GATE_TARGETS_PAIRS));
+}
+
 std::vector<std::vector<std::complex<float>>> Gate::unitary() const {
     if (unitary_data.size() != 2 && unitary_data.size() != 4) {
         throw std::out_of_range(std::string(name) + " doesn't have 1q or 2q unitary data.");

src/stim/gates/gates.h

@@ -258,7 +258,7 @@ struct Gate {
 
     template <size_t W>
     std::vector<Flow<W>> flows() const {
-        if (flags & GateFlags::GATE_IS_UNITARY) {
+        if (has_known_unitary_matrix()) {
             auto t = tableau<W>();
             if (flags & GateFlags::GATE_TARGETS_PAIRS) {
                 return {
@@ -285,6 +285,12 @@ struct Gate {
     bool is_symmetric() const;
     GateType hadamard_conjugated(bool ignoring_sign) const;
 
+    /// Determines if the gate has a specified unitary matrix.
+    ///
+    /// Some unitary gates, such as SPP, don't have a specified matrix because the
+    /// matrix depends crucially on the targets.
+    bool has_known_unitary_matrix() const;
+
     /// Converts a single qubit unitary gate into an euler-angles rotation.
     ///
     /// Returns:

src/stim/gates/gates.pybind.cc

@@ -50,7 +50,7 @@ pybind11::object gate_tableau(const Gate &self) {
     return pybind11::none();
 }
 pybind11::object gate_unitary_matrix(const Gate &self) {
-    if (self.flags & GATE_IS_UNITARY) {
+    if (self.has_known_unitary_matrix()) {
         auto r = self.unitary();
         auto n = r.size();
         std::complex<float> *buffer = new std::complex<float>[n * n];

src/stim/gates/gates.test.cc

@@ -144,7 +144,7 @@ TEST_EACH_WORD_SIZE_W(gate_data, decompositions_are_correct, {
 
 TEST_EACH_WORD_SIZE_W(gate_data, unitary_inverses_are_correct, {
     for (const auto &g : GATE_DATA.items) {
-        if (g.flags & GATE_IS_UNITARY) {
+        if (g.has_known_unitary_matrix()) {
             auto g_t_inv = g.tableau<W>().inverse(false);
             auto g_inv_t = GATE_DATA[g.best_candidate_inverse_id].tableau<W>();
             EXPECT_EQ(g_t_inv, g_inv_t) << g.name;

src/stim/gates/gates_test.py

@@ -59,14 +59,17 @@ def test_gate_data_repr():
 def test_gate_data_inverse():
     for v in stim.gate_data().values():
         assert v.is_unitary == (v.inverse is not None)
-        if v.is_unitary:
-            assert np.allclose(v.unitary_matrix.conj().T, v.inverse.unitary_matrix, atol=1e-6)
+        matrix = v.unitary_matrix
+        if matrix is not None:
+            assert v.is_unitary
+            assert np.allclose(matrix.conj().T, v.inverse.unitary_matrix, atol=1e-6)
             assert v.inverse == v.generalized_inverse
 
     assert stim.gate_data('H').inverse == stim.gate_data('H')
     assert stim.gate_data('S').inverse == stim.gate_data('S_DAG')
     assert stim.gate_data('M').inverse is None
     assert stim.gate_data('CXSWAP').inverse == stim.gate_data('SWAPCX')
+    assert stim.gate_data('SPP').inverse == stim.gate_data('SPP_DAG')
 
     assert stim.gate_data('S').generalized_inverse == stim.gate_data('S_DAG')
     assert stim.gate_data('M').generalized_inverse == stim.gate_data('M')

src/stim/py/stim_pybind_test.py

@@ -136,6 +136,13 @@ def test_main_write_to_file():
             assert "Generated repetition_code" in f.read()
 
 
+def test_main_help(capsys):
+    assert stim.main(command_line_args=["help"]) == 0
+    captured = capsys.readouterr()
+    assert captured.err == ""
+    assert 'Available stim commands' in captured.out
+
+
 def test_main_redirects_stdout(capsys):
     assert stim.main(command_line_args=[
         "gen",

src/stim/simulators/error_analyzer.test.cc

@@ -307,7 +307,7 @@ TEST_EACH_WORD_SIZE_W(ErrorAnalyzer, unitary_gates_match_frame_simulator, {
         data.push_back(GateTarget::qubit(k));
     }
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             e.undo_gate({gate.id, {}, data});
             f.do_gate({gate.inverse().id, {}, data});
             for (size_t q = 0; q < 16; q++) {

src/stim/simulators/frame_simulator.test.cc

@@ -97,7 +97,7 @@ bool is_bulk_frame_operation_consistent_with_tableau(const Gate &gate) {
 
 TEST_EACH_WORD_SIZE_W(FrameSimulator, bulk_operations_consistent_with_tableau_data, {
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             EXPECT_TRUE(is_bulk_frame_operation_consistent_with_tableau<W>(gate)) << gate.name;
         }
     }

src/stim/simulators/graph_simulator.test.cc

@@ -257,7 +257,7 @@ TEST(graph_simulator, all_unitary_gates_work) {
     SpanRef<GateTarget> t1 = t2;
     t1.ptr_end--;
     for (const auto &gate : GATE_DATA.items) {
-        if (!(gate.flags & GATE_IS_UNITARY)) {
+        if (!gate.has_known_unitary_matrix()) {
             continue;
         }
         Circuit circuit;

src/stim/simulators/sparse_rev_frame_tracker.test.cc

@@ -144,7 +144,7 @@ static std::vector<GateTarget> qubit_targets(const std::vector<uint32_t> &target
 TEST_EACH_WORD_SIZE_W(SparseUnsignedRevFrameTracker, fuzz_all_unitary_gates_vs_tableau, {
     auto rng = INDEPENDENT_TEST_RNG();
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             size_t n = (gate.flags & GATE_TARGETS_PAIRS) ? 2 : 1;
             SparseUnsignedRevFrameTracker tracker_gate(n + 3, 0, 0);
             for (size_t q = 0; q < n; q++) {

src/stim/simulators/tableau_simulator.test.cc

@@ -311,7 +311,7 @@ TEST_EACH_WORD_SIZE_W(TableauSimulator, unitary_gates_consistent_with_tableau_da
     auto t = Tableau<W>::random(10, rng);
     TableauSimulator<W> sim(INDEPENDENT_TEST_RNG(), 10);
     for (const auto &gate : GATE_DATA.items) {
-        if (!(gate.flags & GATE_IS_UNITARY)) {
+        if (!gate.has_known_unitary_matrix()) {
             continue;
         }
         sim.inv_state = t;

src/stim/stabilizers/pauli_string_ref.test.cc

@@ -71,7 +71,7 @@ TEST_EACH_WORD_SIZE_W(pauli_string, do_instruction_agrees_with_tableau_sim, {
     sim.inv_state = Tableau<W>::random(sim.inv_state.num_qubits, sim.rng);
 
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             check_pauli_string_do_instruction_agrees_with_tableau_sim<W>(gate, sim);
         }
     }

src/stim/stabilizers/tableau.test.cc

@@ -184,15 +184,15 @@ TEST_EACH_WORD_SIZE_W(tableau, str, {
 
 TEST_EACH_WORD_SIZE_W(tableau, gate_tableau_data_vs_unitary_data, {
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             EXPECT_TRUE(tableau_agrees_with_unitary<W>(gate.tableau<W>(), gate.unitary())) << gate.name;
         }
     }
 })
 
 TEST_EACH_WORD_SIZE_W(tableau, inverse_data, {
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             auto &inv_gate = gate.inverse();
             auto tab = gate.tableau<W>();
             auto inv_tab = inv_gate.tableau<W>();
@@ -1163,7 +1163,7 @@ TEST_EACH_WORD_SIZE_W(tableau, unitary_big_endian, {
 
 TEST_EACH_WORD_SIZE_W(tableau, unitary_vs_gate_data, {
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             std::vector<std::complex<float>> flat_expected;
             for (const auto &row : gate.unitary()) {
                 flat_expected.insert(flat_expected.end(), row.begin(), row.end());

src/stim/util_top/stabilizers_vs_amplitudes.test.cc

@@ -10,7 +10,7 @@ using namespace stim;
 
 TEST_EACH_WORD_SIZE_W(conversions, unitary_to_tableau_vs_gate_data, {
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             EXPECT_EQ(unitary_to_tableau<W>(gate.unitary(), true), gate.tableau<W>()) << gate.name;
         }
     }
@@ -20,7 +20,7 @@ TEST_EACH_WORD_SIZE_W(conversions, tableau_to_unitary_vs_gate_data, {
     VectorSimulator v1(2);
     VectorSimulator v2(2);
     for (const auto &gate : GATE_DATA.items) {
-        if (gate.flags & GATE_IS_UNITARY) {
+        if (gate.has_known_unitary_matrix()) {
             auto actual = tableau_to_unitary<W>(gate.tableau<W>(), true);
             auto expected = gate.unitary();
             v1.state.clear();