Add Surface Code

[apassenger]

Hello,

I tried to write a surface code example. I put here, hope it is fine. I will be happy if you can check :) Since the toric code is written as a hypergraph product code (HGP), the non-full rank matrix of the repeating code is used: this is an n×n matrix with rank n-1. If I remove one of the rows - in the code this is the last row - then I have an n-1×n matrix with full rank. I have replaced the Hr matrix with this matrix.

I attached the result in the absence of errors. However, I am not sure for X stabilizer code.
I also did not add z stabilizers too

import numpy as np
import matplotlib.pyplot as plt
from scipy.sparse import hstack, kron, eye, csr_matrix, block_diag
from pymatching import Matching

def repetition_code(n):
    """
    Parity check matrix of a repetition code with length n.
    """
    row_ind, col_ind = zip(*((i, j) for i in range(n) for j in (i, (i+1)%n)))
    #print("row_ind: ",row_ind)
    if row_ind == n-1: print("Last index:",row_ind)
    data = np.ones(2*n, dtype=np.uint8)
    #print("data matrix: ",data)
    
    data_to_array = csr_matrix((data, (row_ind, col_ind))).toarray()
    print("data array: ",data_to_array)
    #data_to_array=data_to_array[:-1,:]
    #data = data_to_array
    #print(data)
    return csr_matrix((data, (row_ind, col_ind)))

def surface_code(L):
        
    Hr = repetition_code(L)
    data_to_array = Hr.toarray()
    print("Hr_to_array:",data_to_array)
    data_to_array=np.array(data_to_array)
    data_to_array=data_to_array[:-1,:]
    print("Hr_to_array_deleted:",data_to_array)
    print("Hr: ",Hr)
    Hr=data_to_array
    print("new Hr: ",Hr)
    return csr_matrix(Hr)
def surface_code_x_stabilisers(L):
    """
    Sparse check matrix for the X stabilisers of a toric code with
    lattice size L, constructed as the hypergraph product of
    two repetition codes.
    """
    Hr=surface_code(L)
    
    H = hstack(
            [kron(Hr, eye(Hr.shape[1])), kron(eye(Hr.shape[0]), Hr.T)],
            dtype=np.uint8
        )
    #print("H :" ,H)    
    H.data = H.data % 2
    H.eliminate_zeros()
    print(H.toarray())
    return csr_matrix(H)
    
    
#a=repetition_code(4)
a=surface_code_x_stabilisers(4)
noise = np.array([0,0,0,0,0])
#z = H@noise % 2
m=Matching(a)
print(m)
#m.draw()
#c = m.decode(z)
m.draw()
#print(c)
plt.show()
print(a)