"""
Created on Thu Sep 14 19:11:40 2023
@author: laurentdumas
"""
import numpy as np
def fct_gauss(A, b):
    n = len(b)
    u = np.zeros(n)    
    # Étape d'élimination
    for j in range(n):
         #Choix du pivot
        pivot = A[j , j]
        kp = j      
        # Stratégie de pivot
        for k in range(j, n):
            if abs(A[k, j]) > abs(pivot):
               pivot = A[k, j]
               kp = k        
        # Échange de lignes    
        for l in range(j, n):
            cc=A[j,l]
            A[j,l]=A[kp,l]
            A[kp,l]=cc
        cc=b[j]
        b[j]=b[kp]
        b[kp]=cc      
        for i in range(j+1, n):
            m = A[i, j ] / A[j , j ]
            A[i, :] = A[i, :] - A[j , :] * m
            b[i] = b[i] - m * b[j ]    
    # Remontée
    u[n - 1] = b[n - 1] / A[n - 1, n - 1]
    for j in range(n - 2, -1, -1):
        u[j] = (b[j] - np.dot(A[j, j + 1:], u[j + 1:])) / A[j, j] 
    return u
# Exemple d'utilisation
n = 10
A = np.random.rand(n, n)  # Cas aléatoire
#A=np.array([[1E-19,1],[1,1]])
#b=np.array([[1],[1]])
# Cas de la matrice de Hilbert
#for i in range(n):
#    for j in range(n):
#        A[i, j] = 1 / (i + j + 1)

b = np.dot(A, np.ones(n))
u = fct_gauss(A, b)
error_norm_2 = np.linalg.norm(u - np.ones(n))
print('Erreur norme 2:', error_norm_2)
print('Solution:', u)