algo/statiliser.py

class Statiliser

Statiliser is a class that generates a basis of states stabilized by a given set of stabilizers. Stabilizers are represented as strings of Paulis (e.g., "ZZZII", "IIZZZ", "XIXXI", "IXXIX").

generate uses gradient-free optimization to find states that are stabilized by the provided stabilizers, and then applies the Gram-Schmidt process to ensure the resulting states form an orthonormal basis.

name	type	default
basis	Optional[pt.Tensor]	None
num_states	int	None
stabilisers	List[pt.Tensor]	None
sz	int	None

init

Create a generator for the subspace stabilised by stabilisers.

def __init__(self: Any, stabilisers: Sequence[Union[str, pt.Tensor]]) -> Any

Implementation

def __init__(self, stabilisers: Sequence[Union[str, pt.Tensor]]):
    stabilisers = [S(s) if isinstance(s, str) else s for s in stabilisers]
    self.stabilisers = list(stabilisers)
    self.sz = int(pt.log2(pt.tensor(stabilisers[0].shape[0])).item())
    self.num_states = 2 ** (self.sz - len(stabilisers))
    self.basis = None

generate

Generate an orthonormal basis for the stabilised subspace.

def generate(self: Any, mode: str, tol: float, minimal: bool) -> pt.Tensor

Implementation

def generate(self, mode: str='real', tol: float=1e-06, minimal: bool=True) -> pt.Tensor:
    basis: List[pt.Tensor] = []
    factor = 2 if mode == 'complex' else 1
    for _ in range(self.num_states):
        res = minimize(self._fun, x0=pt.rand(2 ** self.sz * factor).numpy(), args=(mode, int(minimal)), method='Powell', tol=tol).x
        state = pt.tensor(res, dtype=pt.float32 if mode == 'real' else pt.cfloat)
        if mode != 'real':
            state = self._to_complex(state)
        state = state / pt.norm(state)
        basis.append(state)
    basis_t = GramSchmidt(basis)
    basis_t = basis_t.to(dtype=pt.float16 if mode == 'real' else pt.cfloat)
    self.basis = basis_t
    return basis_t

S

Convert a Pauli string (e.g. "XIIZ") into its operator tensor.

def S(string: str) -> pt.Tensor

Implementation

def S(string: str) -> pt.Tensor:
    paulis = [_paulis[i] for i in string]
    return _S(*paulis)

GramSchmidt

Orthonormalise a list/sequence of vectors using Gram-Schmidt.

Vectors with (near-)zero residual norm are dropped.

def GramSchmidt(vectors: Sequence[pt.Tensor]) -> pt.Tensor

Implementation

def GramSchmidt(vectors: Sequence[pt.Tensor]) -> pt.Tensor:
    ortho: List[pt.Tensor] = []
    for v in vectors:
        w = v - sum((v @ u.conj() * u for u in ortho))
        if pt.norm(w) > 1e-08:
            ortho.append(w / pt.norm(w))
    return pt.stack(ortho)