import numpy as np
rng = np.random.default_rng(0)  #generator
#----------- Input ---------------------------------
n = 6  # dimensiunea matricei
A = rng.integers(1,11, size = (n,n))
A = np.triu(A,-1)
A[n-1,:] = rng.integers(1,11, size = (1,n))
A =A.astype('float')
b =  rng.integers(-10,11, size = (n,1)).astype('float')
#-------------------------------------------------------------
print('-----------A------------')
print(A)
print('-----------b------------')
print(b)

def EGPP(A):
    n = len(A)

    p = np.zeros(len(A), dtype = int)

    for k in range(n - 1):
        if abs(A[k + 1][k]) > abs(A[k][k]):
            ik = k + 1
            p[k] = ik
            for i in range(k, n):
                A[ik][i], A[k][i] = A[k][i], A[ik][i]

        A[k + 1][k] /= A[k][k]

        for j in range(k + 1, n):
            A[k + 1][j] -= A[k + 1][j] * A[k][j]

    return A, p

def S_SUP_TR(U, b):
    n = len(U)

    rezultat = np.zeros(n, dtype = float)

    for i in range(n - 1, -1, -1):
        s = b[i]
        for j in range(i + 1, n):
            s -= U[i][j] * rezultat[j]
        rezultat[i] = s / U[i][i]

    return rezultat

def SL_GPP(A, b):
    U, p = EGPP(A)

    n = len(A)

    for k in range(n - 1):
        b[k], b[p[k]] = b[p[k]], b[k]
        for i in range(k + 1, n):
            b[i] -= A[i][k] * b[k]
    x = S_SUP_TR(U, b)

    return x

solutie = SL_GPP(A, b)
print(solutie)

inv = np.linalg.inv(A)

if solutie is inv * b:
    print("Solutia gasita e corecta")
else:
    print("Solutia gasita e gresita.........")