Skip to content

Commit

Permalink
Clean up stim.Tableau.from_stabilizers a bit (#716)
Browse files Browse the repository at this point in the history
- Switch to streaming computations of the gates
- Improve the anticommutation check from O(n^2) additional work to O(n)
additional work
- Fix a sizing error and test that it's fixed
- Cut costs on slowest test
  • Loading branch information
Strilanc authored Mar 15, 2024
1 parent 1cdf45f commit 2161fc9
Show file tree
Hide file tree
Showing 6 changed files with 95 additions and 92 deletions.
4 changes: 2 additions & 2 deletions src/stim/mem/simd_word.test.cc
Original file line number Diff line number Diff line change
Expand Up @@ -33,14 +33,14 @@ union WordOr64 {

TEST_EACH_WORD_SIZE_W(simd_word_pick, popcount, {
WordOr64 v;
auto n = sizeof(simd_word<W>) * 8;
auto n = sizeof(simd_word<W>) * 2;

for (size_t expected = 0; expected <= n; expected++) {
std::vector<uint64_t> bits{};
for (size_t i = 0; i < n; i++) {
bits.push_back(i < expected);
}
for (size_t reps = 0; reps < 100; reps++) {
for (size_t reps = 0; reps < 10; reps++) {
std::shuffle(bits.begin(), bits.end(), INDEPENDENT_TEST_RNG());
for (size_t i = 0; i < n; i++) {
v.p[i >> 6] = 0;
Expand Down
26 changes: 13 additions & 13 deletions src/stim/simulators/tableau_simulator.test.cc
Original file line number Diff line number Diff line change
Expand Up @@ -488,24 +488,24 @@ bool vec_sim_corroborates_measurement_process(

TEST_EACH_WORD_SIZE_W(TableauSimulator, measurement_vs_vector_sim, {
auto rng = INDEPENDENT_TEST_RNG();
for (size_t k = 0; k < 10; k++) {
for (size_t k = 0; k < 5; k++) {
auto state = Tableau<W>::random(2, rng);
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {1}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 1}));
}
for (size_t k = 0; k < 10; k++) {
for (size_t k = 0; k < 5; k++) {
auto state = Tableau<W>::random(4, rng);
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 1}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {2, 1}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 1, 2, 3}));
}
{
auto state = Tableau<W>::random(12, rng);
auto state = Tableau<W>::random(8, rng);
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 1, 2, 3}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 10, 11}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {11, 5, 7}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 6, 7}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {7, 3, 4}));
ASSERT_TRUE(vec_sim_corroborates_measurement_process(state, {0, 1, 2, 3, 4, 5, 6, 7}));
}
})

Expand Down Expand Up @@ -613,7 +613,7 @@ TEST_EACH_WORD_SIZE_W(TableauSimulator, correlated_error, {
expected);

int hits[3]{};
size_t n = 10000;
size_t n = 5000;
for (size_t k = 0; k < n; k++) {
auto sample = TableauSimulator<W>::sample_circuit(
Circuit(R"circuit(
Expand All @@ -628,8 +628,8 @@ TEST_EACH_WORD_SIZE_W(TableauSimulator, correlated_error, {
hits[2] += sample[2];
}
ASSERT_TRUE(0.45 * n < hits[0] && hits[0] < 0.55 * n);
ASSERT_TRUE((0.125 - 0.05) * n < hits[1] && hits[1] < (0.125 + 0.05) * n);
ASSERT_TRUE((0.28125 - 0.05) * n < hits[2] && hits[2] < (0.28125 + 0.05) * n);
ASSERT_TRUE((0.125 - 0.08) * n < hits[1] && hits[1] < (0.125 + 0.08) * n);
ASSERT_TRUE((0.28125 - 0.08) * n < hits[2] && hits[2] < (0.28125 + 0.08) * n);
})

TEST_EACH_WORD_SIZE_W(TableauSimulator, quantum_cannot_control_classical, {
Expand Down Expand Up @@ -865,12 +865,12 @@ TEST_EACH_WORD_SIZE_W(TableauSimulator, peek_bloch, {

TEST_EACH_WORD_SIZE_W(TableauSimulator, paulis, {
auto rng = INDEPENDENT_TEST_RNG();
TableauSimulator<W> sim1(INDEPENDENT_TEST_RNG(), 500);
TableauSimulator<W> sim2(INDEPENDENT_TEST_RNG(), 500);
sim1.inv_state = Tableau<W>::random(500, rng);
TableauSimulator<W> sim1(INDEPENDENT_TEST_RNG(), 300);
TableauSimulator<W> sim2(INDEPENDENT_TEST_RNG(), 300);
sim1.inv_state = Tableau<W>::random(300, rng);
sim2.inv_state = sim1.inv_state;

sim1.paulis(PauliString<W>(500));
sim1.paulis(PauliString<W>(300));
ASSERT_EQ(sim1.inv_state, sim2.inv_state);
sim1.paulis(PauliString<W>(5));
ASSERT_EQ(sim1.inv_state, sim2.inv_state);
Expand Down
120 changes: 55 additions & 65 deletions src/stim/stabilizers/conversions.inl
Original file line number Diff line number Diff line change
Expand Up @@ -558,10 +558,9 @@ Tableau<W> stabilizers_to_tableau(
num_qubits = std::max(num_qubits, e.num_qubits);
}

simd_bit_table<W> buf_xs(num_qubits, stabilizers.size());
simd_bit_table<W> buf_zs(num_qubits, stabilizers.size());
simd_bit_table<W> buf_xs(stabilizers.size(), num_qubits);
simd_bit_table<W> buf_zs(stabilizers.size(), num_qubits);
simd_bits<W> buf_signs(stabilizers.size());
simd_bits<W> buf_workspace(stabilizers.size());
for (size_t k = 0; k < stabilizers.size(); k++) {
memcpy(buf_xs[k].u8, stabilizers[k].xs.u8, stabilizers[k].xs.num_u8_padded());
memcpy(buf_zs[k].u8, stabilizers[k].zs.u8, stabilizers[k].zs.num_u8_padded());
Expand All @@ -570,19 +569,6 @@ Tableau<W> stabilizers_to_tableau(
buf_xs = buf_xs.transposed();
buf_zs = buf_zs.transposed();

for (size_t k1 = 0; k1 < stabilizers.size(); k1++) {
for (size_t k2 = k1 + 1; k2 < stabilizers.size(); k2++) {
if (!stabilizers[k1].ref().commutes(stabilizers[k2])) {
std::stringstream ss;
ss << "Some of the given stabilizers anticommute.\n";
ss << "For example:\n ";
ss << stabilizers[k1];
ss << "\nanticommutes with\n";
ss << stabilizers[k2] << "\n";
throw std::invalid_argument(ss.str());
}
}
}
Circuit elimination_instructions;

size_t used = 0;
Expand All @@ -597,11 +583,6 @@ Tableau<W> stabilizers_to_tableau(

// Check for incompatible / redundant stabilizers.
if (pivot == num_qubits) {
for (size_t q = 0; q < num_qubits; q++) {
if (buf_xs[q][k]) {
throw std::invalid_argument("Some of the given stabilizers anticommute.");
}
}
if (buf_signs[k]) {
throw std::invalid_argument("Some of the given stabilizers contradict each other.");
}
Expand All @@ -620,19 +601,21 @@ Tableau<W> stabilizers_to_tableau(
CircuitInstruction instruction{g, {}, &t};
elimination_instructions.safe_append(instruction);
size_t q = pivot;
simd_bits_range_ref<W> xs1 = buf_xs[q];
simd_bits_range_ref<W> zs1 = buf_zs[q];
simd_bits_range_ref<W> ss = buf_signs;
switch (g) {
case GateType::H_YZ:
buf_xs[q] ^= buf_zs[q];
buf_workspace = buf_zs[q];
buf_workspace.invert_bits();
buf_workspace &= buf_xs[q];
buf_signs ^= buf_workspace;
ss.for_each_word(xs1, zs1, [](auto &s, auto &x, auto &z) {
x ^= z;
s ^= z.andnot(x);
});
break;
case GateType::H:
buf_xs[q].swap_with(buf_zs[q]);
buf_workspace = buf_zs[q];
buf_workspace &= buf_xs[q];
buf_signs ^= buf_workspace;
ss.for_each_word(xs1, zs1, [](auto &s, auto &x, auto &z) {
std::swap(x, z);
s ^= x & z;
});
break;
default:
throw std::invalid_argument("Unrecognized gate type.");
Expand All @@ -649,47 +632,34 @@ Tableau<W> stabilizers_to_tableau(
elimination_instructions.safe_append(instruction);
size_t q1 = targets[0].qubit_value();
size_t q2 = targets[1].qubit_value();
simd_bits_range_ref<W> x1 = buf_xs[q1];
simd_bits_range_ref<W> z1 = buf_zs[q1];
simd_bits_range_ref<W> x2 = buf_xs[q2];
simd_bits_range_ref<W> z2 = buf_zs[q2];
simd_bits_range_ref<W> ss = buf_signs;
simd_bits_range_ref<W> xs1 = buf_xs[q1];
simd_bits_range_ref<W> zs1 = buf_zs[q1];
simd_bits_range_ref<W> xs2 = buf_xs[q2];
simd_bits_range_ref<W> zs2 = buf_zs[q2];
switch (g) {
case GateType::XCX:
buf_workspace = x1;
buf_workspace ^= x2;
buf_workspace &= z1;
buf_workspace &= z2;
buf_signs ^= buf_workspace;
x1 ^= z2;
x2 ^= z1;
ss.for_each_word(xs1, zs1, xs2, zs2, [](auto &s, auto &x1, auto &z1, auto &x2, auto &z2) {
s ^= (x1 ^ x2) & z1 & z2;
x1 ^= z2;
x2 ^= z1;
});
break;
case GateType::XCY:
x1 ^= x2;
x1 ^= z2;
x2 ^= z1;
z2 ^= z1;
buf_workspace = x1;
buf_workspace |= x2;
buf_workspace.invert_bits();
buf_workspace &= z1;
buf_workspace &= z2;
buf_signs ^= buf_workspace;
buf_workspace = z2;
buf_workspace.invert_bits();
buf_workspace &= z1;
buf_workspace &= x1;
buf_workspace &= x2;
buf_signs ^= buf_workspace;
ss.for_each_word(xs1, zs1, xs2, zs2, [](auto &s, auto &x1, auto &z1, auto &x2, auto &z2) {
x1 ^= x2 ^ z2;
x2 ^= z1;
z2 ^= z1;
s ^= x1.andnot(z1) & x2.andnot(z2);
s ^= x1 & z1 & z2.andnot(x2);
});
break;
case GateType::XCZ:
z2 ^= z1;
x1 ^= x2;
buf_workspace = z2;
buf_workspace ^= x1;
buf_workspace.invert_bits();
buf_workspace &= x2;
buf_workspace &= z1;
buf_signs ^= buf_workspace;
ss.for_each_word(xs1, zs1, xs2, zs2, [](auto &s, auto &x1, auto &z1, auto &x2, auto &z2) {
z2 ^= z1;
x1 ^= x2;
s ^= (z2 ^ x1).andnot(z1 & x2);
});
break;
default:
throw std::invalid_argument("Unrecognized gate type.");
Expand Down Expand Up @@ -717,6 +687,26 @@ Tableau<W> stabilizers_to_tableau(
used++;
}

// All stabilizers will have been mapped into Z products, if they commuted.
for (size_t q = 0; q < num_qubits; q++) {
if (buf_xs[q].not_zero()) {
for (size_t k1 = 0; k1 < stabilizers.size(); k1++) {
for (size_t k2 = k1 + 1; k2 < stabilizers.size(); k2++) {
if (!stabilizers[k1].ref().commutes(stabilizers[k2])) {
std::stringstream ss;
ss << "Some of the given stabilizers anticommute.\n";
ss << "For example:\n ";
ss << stabilizers[k1];
ss << "\nanticommutes with\n";
ss << stabilizers[k2] << "\n";
throw std::invalid_argument(ss.str());
}
}
}
throw std::invalid_argument("The given stabilizers commute but the solver failed in a way that suggests they anticommute. Please report this as a bug.");
}
}

if (used < num_qubits) {
if (!allow_underconstrained) {
throw std::invalid_argument(
Expand Down
15 changes: 14 additions & 1 deletion src/stim/stabilizers/conversions.test.cc
Original file line number Diff line number Diff line change
Expand Up @@ -217,7 +217,7 @@ TEST_EACH_WORD_SIZE_W(conversions, stabilizer_state_vector_to_circuit_basic, {
TEST_EACH_WORD_SIZE_W(conversions, stabilizer_state_vector_to_circuit_fuzz_round_trip, {
auto rng = INDEPENDENT_TEST_RNG();
for (const auto &little_endian : std::vector<bool>{false, true}) {
for (size_t n = 0; n < 10; n++) {
for (size_t n = 0; n < 5; n++) {
// Pick a random stabilizer state.
TableauSimulator<W> sim(INDEPENDENT_TEST_RNG(), n);
sim.inv_state = Tableau<W>::random(n, rng);
Expand Down Expand Up @@ -663,6 +663,19 @@ TEST_EACH_WORD_SIZE_W(conversions, stabilizer_to_tableau_detect_anticommutation,
ASSERT_THROW({ stabilizers_to_tableau<W>(input_stabilizers, false, false, false); }, std::invalid_argument);
})

TEST_EACH_WORD_SIZE_W(conversions, stabilizer_to_tableau_size_affecting_redundancy, {
std::vector<stim::PauliString<W>> input_stabilizers;
input_stabilizers.push_back(PauliString<W>::from_str("X_"));
input_stabilizers.push_back(PauliString<W>::from_str("_X"));
for (size_t k = 0; k < 150; k++) {
input_stabilizers.push_back(PauliString<W>::from_str("__"));
}
auto t = stabilizers_to_tableau<W>(input_stabilizers, true, true, false);
ASSERT_EQ(t.num_qubits, 2);
ASSERT_EQ(t.zs[0], PauliString<W>::from_str("X_"));
ASSERT_EQ(t.zs[1], PauliString<W>::from_str("_X"));
})

TEST(conversions, independent_to_disjoint_xyz_errors) {
double out_x;
double out_y;
Expand Down
20 changes: 10 additions & 10 deletions src/stim/stabilizers/tableau.test.cc
Original file line number Diff line number Diff line change
Expand Up @@ -722,15 +722,15 @@ TEST_EACH_WORD_SIZE_W(tableau, expand_pad_equals, {

TEST_EACH_WORD_SIZE_W(tableau, transposed_access, {
auto rng = INDEPENDENT_TEST_RNG();
size_t n = 1000;
size_t n = W > 256 ? 1000 : 400;
Tableau<W> t(n);
auto m = t.xs.xt.data.num_bits_padded();
t.xs.xt.data.randomize(m, rng);
t.xs.zt.data.randomize(m, rng);
t.zs.xt.data.randomize(m, rng);
t.zs.zt.data.randomize(m, rng);
for (size_t inp_qubit = 0; inp_qubit < 1000; inp_qubit += 99) {
for (size_t out_qubit = 0; out_qubit < 1000; out_qubit += 99) {
for (size_t inp_qubit = 0; inp_qubit < n; inp_qubit += 99) {
for (size_t out_qubit = 0; out_qubit < n; out_qubit += 99) {
bool bxx = t.xs.xt[inp_qubit][out_qubit];
bool bxz = t.xs.zt[inp_qubit][out_qubit];
bool bzx = t.zs.xt[inp_qubit][out_qubit];
Expand Down Expand Up @@ -885,26 +885,26 @@ TEST_EACH_WORD_SIZE_W(tableau, transposed_xz_input, {

TEST_EACH_WORD_SIZE_W(tableau, direct_sum, {
auto rng = INDEPENDENT_TEST_RNG();
auto t1 = Tableau<W>::random(260, rng);
auto t2 = Tableau<W>::random(270, rng);
auto t1 = Tableau<W>::random(160, rng);
auto t2 = Tableau<W>::random(170, rng);
auto t3 = t1;
t3 += t2;
ASSERT_EQ(t3, t1 + t2);

PauliString<W> p1 = t1.xs[5];
p1.ensure_num_qubits(260 + 270, 1.0);
p1.ensure_num_qubits(160 + 170, 1.0);
ASSERT_EQ(t3.xs[5], p1);

std::string p2 = t2.xs[6].str();
std::string p3 = t3.xs[266].str();
std::string p3 = t3.xs[166].str();
ASSERT_EQ(p2[0], p3[0]);
p2 = p2.substr(1);
p3 = p3.substr(1);
for (size_t k = 0; k < 260; k++) {
for (size_t k = 0; k < 160; k++) {
ASSERT_EQ(p3[k], '_');
}
for (size_t k = 0; k < 270; k++) {
ASSERT_EQ(p3[260 + k], p2[k]);
for (size_t k = 0; k < 170; k++) {
ASSERT_EQ(p3[160 + k], p2[k]);
}
})

Expand Down
2 changes: 1 addition & 1 deletion src/stim/stabilizers/tableau_iter.test.cc
Original file line number Diff line number Diff line change
Expand Up @@ -135,7 +135,7 @@ TEST_EACH_WORD_SIZE_W(tableau_iter, iter_tableau, {
ASSERT_EQ(n1, 6);
ASSERT_EQ(s1, 24);
ASSERT_EQ(n2, 720);
ASSERT_EQ(n3, 1451520);
// ASSERT_EQ(n3, 1451520); // Note: disabled because it takes 2-3 seconds.
})

TEST_EACH_WORD_SIZE_W(tableau_iter, iter_tableau_distinct, {
Expand Down

0 comments on commit 2161fc9

Please sign in to comment.