Are there any packages that can calculate stabilizer tableau of a QECC

Question

I'm experimenting with some small quantum error correcting codes (QECC). For example $[[5,1,3]]$, $[[8,3,3]]$ or toric codes $[[2d^2,2,d]]$ ($d=2,3,\cdots$). The last one being defined by redundant stabilizers. What packages can take in a set of $m'$ stabilizers and produce the tableau? ($m'$ could be larger than $m=n-k$ in case stabilizers are not independent).

Here are the details for the $[[8,3,3]]$ code : the code is in standard form; $1=X$, $2=Z$, $3=XZ$;

[[3,0,2,2,1,1,3,0],
[2,3,0,2,1,3,0,1],
[2,0,3,0,1,2,1,3],
[2,2,0,1,0,1,3,3],
[2,2,2,2,2,2,2,2],
[2,2,0,2,0,2,0,0],
[2,0,2,2,0,0,2,0],
[0,2,2,2,0,0,0,2],
[2,0,0,0,0,0,0,0],
[0,2,0,0,0,0,0,0],
[0,0,2,0,0,0,0,0],
[0,0,0,2,0,0,0,0],
[0,0,0,0,1,0,0,0],
[0,2,2,0,1,1,0,0],
[2,0,0,2,1,0,1,0],
[0,0,2,2,1,0,0,1]]

first 5 rows of the matrix above are stabilizers; next 3 are logical Z; next 5 are destabilizers; last 3 logical X.

Answer 1

You can use stim for this, although you do have to write the stabilizer projection procedure for yourself.

Write some methods to project a system into the +1 eigenstate of several stabilizers:

from typing import List

import stim

def find_compatible_tableau(stabilizers: List[stim.PauliString]) -> stim.Tableau:
    num_qubits = max(len(e) for e in stabilizers)
    sim = stim.TableauSimulator()
    # Start the target qubits in a state that overlaps all stabilizers.
    for q in range(num_qubits):
        sim.h(q)
        sim.cnot(q, q + num_qubits + 1)
    # Project into each stabilizer's +1 eigenbasis.
    for s in stabilizers:
        project_stabilizer(sim, s, ancilla=num_qubits)
    # Discard ancillary qubits.
    sim.set_num_qubits(num_qubits)

    # Simulator happens to track the inverse tableau.
    # Invert it to get the normal one.
    return sim.current_inverse_tableau()**-1


def project_stabilizer(sim: stim.TableauSimulator, stabilizer: stim.PauliString, ancilla: int):
    assert ancilla >= len(stabilizer)
    sim.reset(ancilla)
    sim.h(ancilla)
    for q, p in enumerate(stabilizer):
        if p == 1:
            sim.cnot(ancilla, q)
        elif p == 2:
            sim.cy(ancilla, q)
        elif p == 3:
            sim.cz(ancilla, q)
    if stabilizer.sign == -1:
        sim.z(ancilla)
    sim.h(ancilla)
    returned_true, kickback = sim.measure_kickback(ancilla)
    if returned_true:
        if kickback is None:
            raise ValueError("Contradictory stabilizers.")
        sim.do(kickback)

Use it on your problem:

solved_tableau = find_compatible_tableau(stabilizers=[
    stim.PauliString("+ZZ_"),
    stim.PauliString("+_ZZ"),
    stim.PauliString("-XXX"),
])

print(repr(solved_tableau))

And voilà:

stim.Tableau.from_conjugated_generators(
    xs=[
        stim.PauliString("-X__"),
        stim.PauliString("+_X_"),
        stim.PauliString("+__Z"),
    ],
    zs=[
        stim.PauliString("+Z_Z"),
        stim.PauliString("+_ZZ"),
        stim.PauliString("-XXX"),
    ],
)

The stabilizers are a bit re-arranged, but the table is consistent with the ones that were asked for.