Matrix

HNF(A)

Computes row-style Hermite Normal Form of a integer matrix A.

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def HNF(A: MatrixInt) -> Tuple[MatrixInt, SquareMatrixInt, int]:
    r'''
    Computes row-style Hermite Normal Form of a integer matrix A.

    Args:

    Returns:

    '''
    H = A.copy()
    m, n = H.shape
    p = min(m,n)
    k, j = 0, 0

    U = np.identity(m, dtype=int)
    detU = 1


    while j != p:
        # Choose pivot
        col = H[k:, j]
        non_zero = col[col != 0]
        if len(non_zero) == 0:
            j += 1
            k += 1
            continue
        min_val = np.min(np.abs(non_zero))
        i0 = np.where(np.abs(col) == min_val)[0][0] + k
        if i0 > k:
            H[[k, i0]] = H[[i0, k]]
            detU *= -1
            U = elementary_row_swap(m, k, i0) @ U

        if H[k,j] < 0:
            H[k] = -H[k]
            detU *= -1
            U = elementary_row_mul(m, k, -1) @ U

        # Reduce Rows
        b = H[k,j]
        for i in range(k+1, m):
            q = round(H[i,j] / b)
            H[i] -= q * H[k]
            U = elementary_row_add(m, i, k, -q) @ U

        # Check if column is done
        if np.all(H[k+1:, j] == 0):
            j += 1
            k += 1

    # Final reductions
    k = 0
    for j in range(p):
        if H[k,j] < 0:
            H[k] = -H[k]
            U = elementary_row_mul(m, k, -1) @ U
            detU *= -1

        b = H[k,j]
        if b == 0: continue
        for i in range(k):
            q = H[i,j] // b
            H[i] -= q * H[k]
            U = elementary_row_add(m, i, k, -q) @ U

        k += 1

    return H, U, detU

cofactor(A, i, j)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def cofactor(A: SquareMatrixInt, i: int, j: int) -> int:
    r'''

    Args:

    Returns:

    '''
    return minor(A, i, j) * ((-1) ** (i + 1 + j + 1))

cofactor_matrix(A)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def cofactor_matrix(A: SquareMatrixInt) -> SquareMatrixInt:
    r'''

    Args:

    Returns:

    '''
    n = A.shape[0]
    C = np.zeros((n,n), dtype=int)
    for i in range(n):
        for j in range(n):
            C[i,j] = cofactor(A, i, j)
    return C

det(A)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def det(A: SquareMatrixInt) -> int:
    r'''

    Args:

    Returns:

    '''
    H, U, detU = HNF(A)
    return np.prod(np.diagonal(H)) * detU

dual_q_ary_lattice_basis(A, modulus)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def dual_q_ary_lattice_basis(A: MatrixInt, modulus: int) -> SquareMatrixModInt:
    r'''

    Args:

    Returns:

    '''
    m, n = A.shape
    Y = mod_left_kernel(A, modulus)
    assert Y.shape == (m - n, m)

    q = modulus
    return np.block([[Y], [np.zeros((m-n, m-n), dtype=int), q * np.identity(n, dtype=int)]])

elementary_row_add(n, i, j, s)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def elementary_row_add(n: int, i: int, j: int, s: int):
    r'''

    Args:

    Returns:

    '''
    E = np.identity(n, int)
    E[i, j] = s
    return E

elementary_row_mul(n, i, s)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def elementary_row_mul(n: int, i: int, s: int):
    r'''

    Args:

    Returns:

    '''
    E = np.identity(n, int)
    E[i, i] = s
    return E

elementary_row_swap(n, i, j)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def elementary_row_swap(n: int, i: int, j: int):
    r'''

    Args:

    Returns:

    '''
    E = np.identity(n, int)
    E[[i,j]] = E[[j, i]]
    return E

left_kernel(A)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def left_kernel(A: MatrixInt) -> MatrixInt|None:
    r'''

    Args:

    Returns:

    '''
    H, U, _ = HNF(A)
    r = 0
    for row in H[::-1]:
        if np.all(row == 0):
            r += 1
        else:
            break
    if r == 0:
        return None
    return U[-r::]

matrix_modinv(A, modulus)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def matrix_modinv(A: SquareMatrixInt, modulus: int) -> SquareMatrixModInt:
    r'''

    Args:

    Returns:

    '''
    C = cofactor_matrix(A) % modulus
    det_inv = integer_ring.modinv(det(A), modulus)

    return (det_inv * C.T) % modulus

minor(A, i, j)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def minor(A: SquareMatrixInt, i: int, j: int) -> int:
    r'''

    Args:

    Returns:

    '''
    return det(np.delete(np.delete(A, i, axis=0), j, axis=1))

mod_REF(A, modulus)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def mod_REF(A: MatrixInt, modulus: int) -> Tuple[MatrixModInt, SquareMatrixModInt]:
    r'''

    Args:

    Returns:

    '''
    m, n = A.shape
    inv = lambda a: integer_ring.modinv(a, modulus)

    M = A.copy()
    M %= modulus

    U = np.identity(m, dtype=int)

    for j in range(min(m,n)):
        #print(f"j= {j}")
        #print(M)
        col = M[ : ,j]
        nonzero = np.nonzero(col[j:])[0]
        if len(nonzero) == 0:
            # no pivot, in this column
            continue
        pivot_i = nonzero[0] + j

        #print("current column:")
        #print(col)
        #print(f"pivot index: {pivot_i}")
        #print(f"pivot value: {M[pivot_i, j]}")


        if pivot_i != j:
            row_swap(M, pivot_i, j)
            U = (elementary_row_swap(m, pivot_i, j) @ U) % modulus

        #print("After swap")
        #print(M)


        pivot_inv = inv(M[j,j])
        row_scale(M, j, pivot_inv)
        U = (elementary_row_mul(m, j, pivot_inv) @ U) % modulus

        M %= modulus

        #print("After scale")
        #print(M)


        for i in range(j + 1, m):
            if M[i, j] != 0:
                U = (elementary_row_add(m, i, j, -M[i,j]) @ U) % modulus
                row_add(M, i, j, -M[i,j])
                M %= modulus

        #print("After reduce")
        #print(M)
        #print("===============")

    return M % modulus, U % modulus

mod_RREF(A, modulus)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def mod_RREF(A: MatrixInt, modulus: int) -> Tuple[MatrixModInt, SquareMatrixModInt]:
    r'''

    Args:

    Returns:

    '''
    M, U = mod_REF(A, modulus)
    m, n = M.shape

    h = m - 1
    for j in range(n - 1, -1, -1):
        # check if current column is a pivot column
        if M[h, j] != 1:
            # skip zero rows at the bottom
            if M[h, j] == 0:
                h -= 1
            continue

        # reduce column's entries below pivot
        for i in range(h - 1, -1, -1):
            coeff = M[i,j]
            row_add(M, i, h, -coeff)
            U = elementary_row_add(m, i, h, -coeff) @ U

        # move to the next row
        h -= 1


    return M % modulus, U % modulus

mod_left_kernel(A, modulus)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def mod_left_kernel(A: MatrixInt, modulus: int) -> MatrixInt|None:
    r'''

    Args:

    Returns:

    '''
    G, U = mod_RREF(A, modulus)
    r = 0
    for row in G[::-1]:
        if np.all(row == 0):
            r += 1
        else:
            break
    if r == 0:
        return None
    return U[-r::]

nullity(A)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def nullity(A: MatrixInt) -> int:
    r'''

    Args:

    Returns:

    '''
    H, U, _ = HNF(A)
    r = 0
    for row in H[::-1]:
        if np.all(row == 0):
            r += 1
        else:
            break
    return r

q_ary_lattice_basis(A, modulus)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def q_ary_lattice_basis(A: MatrixInt, modulus: int) -> SquareMatrixInt:
    r'''

    Args:

    Returns:

    '''
    m, n = A.shape
    assert n >= m
    # A = (A1 | A2)
    A1 = A[ : ,:m] # m x m
    A2 = A[ : ,m:] # m x (n - m)

    B11 = np.identity(m, dtype=int) # m x m
    B12 = (matrix_modinv(A1, modulus) @ A2) % modulus # m x (n - m)
    B21 = np.zeros((n - m, m), dtype=int)
    B22 = modulus * np.identity(n - m, dtype=int)

    print(B11)
    print()
    print(B12)
    print()
    print(B21)
    print()
    print(B22)

    return np.block([[B11, B12], [B21, B22]])

row_add(M, i, k, s)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def row_add(M: MatrixInt, i: int, k: int, s: int):
    r'''

    Args:

    Returns:

    '''
    M[i] += s * M[k]

row_scale(M, i, s)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def row_scale(M: MatrixInt, i: int, s: int):
    r'''

    Args:

    Returns:

    '''
    M[i] *= s

row_swap(M, i, j)

Args:

Returns:

Source code in src/lbpqc/primitives/matrix.py

def row_swap(M: MatrixInt, i: int, j: int):
    r'''

    Args:

    Returns:

    '''
    M[[i,j]] = M[[j,i]]