Méthode de gauss pour la résolution d'un système linéaire

Le tableau ci-dessous énumère trois méthodes directes populaires, chacune d'entre elles utilisant des opérations élémentaires pour produire sa propre forme finale d'équations faciles à résoudre.

•  \(U\) : Matrice triangulaire supérieure
•  \(L\) : Matrice triangulaire inférieure
•  \(I\) : Matrice identité

Élimination de Gauss

L'élimination de Gauss est la méthode la plus familière pour résoudre un système équations linéaires. Elle se compose de deux parties : la phase d'élimination et la phase de substitutions.

•   La fonction de la phase d'élimination est de transformer le Système sous la forme \(Ux = c\).
•  Le système est ensuite résolu par substitution.

\begin{align*} 4x_1-2x_2 +3x_3& = 11 \tag{a}\\ -2x_1+4x_2 -2x_3& = -16 \tag{b}\\ x_1-2x_2 +4x_3& = 17 \tag{c} \end{align*}

Phase d'élimination

La phase d'élimination n'utilise qu'une seule des opérations élémentaires—Multiplier une équation (disons l'équation j) par une constante \(\lambda\) et la soustraire d'une autre équation (équation i). \begin{equation} Eq.(i) \leftarrow Eq.(i) - \lambda \times Eq.(j) \tag{1} \end{equation}

L'équation à soustraire, à savoir l'équation (j), est appelée l'équation du pivot.
Nous commençons l'élimination en prenant l'équation (a) comme équation pivot et en choisissant les multiplicateurs \(\lambda\) de manière à éliminer \(x_1\) dans les équations (b) et (c) : \begin{align*} Eq.(b) \leftarrow Eq.(b) - (-0.5) \times Eq.(a) \\ Eq.(c) \leftarrow Eq.(c) - (0.25) \times Eq.(a) \end{align*} Après cette transformation, les équations deviennent : \begin{align*} 4x_1-2x_2 +3x_3& = 11 \tag{a}\\ 3x_2 -1.5x_3& = -10.5 \tag{b}\\ -1.5x_2 +3.75x_3& = 14.25 \tag{c} \end{align*}

Maintenant, nous choisissons (b) comme équation de pivot et éliminons $x_2$ de (c) : \begin{align*} Eq.(c) \leftarrow Eq.(c) - (-0.5) \times Eq.(b) \end{align*} ce qui donne les équations suivantes : \begin{align*} 4x_1-2x_2 +3x_3& = 11 \tag{a}\\ 3x_2 -1.5x_3& = -10.5 \tag{b}\\ 3x_3& = 9 \tag{c} \end{align*}

Comme indiqué précédemment, la matrice de coefficients augmentés est un instrument plus pratique pour effectuer les calculs. Ainsi, les équations originales seraient écrites comme : \begin{equation} \left[ \begin{matrix} 4& -2& 1\\ -2& 4& -2\\ 1&-2&4 \end{matrix} \left| \, \begin{matrix} 11 \\ -16 \\ 17 \\ \end{matrix} \right. \right] \tag{2} \end{equation} et les équations équivalentes produites par le premier et le second passage de l'élimination de Gauss seraient les suivantes : \begin{equation} \left[ \begin{matrix} 4& -2& 1\\ 0& 3& -1.5\\ 0&-1.5&3.75 \end{matrix} \left| \, \begin{matrix} 11 \\ -10.5 \\ 14.25 \\ \end{matrix} \right. \right] \tag{3} \end{equation} \begin{equation} \left[ \begin{matrix} 4& -2& 1\\ 0& 3& -1.5\\ 0&0&3 \end{matrix} \left| \, \begin{matrix} 11 \\ -10.5 \\ 9 \\ \end{matrix} \right. \right] \tag{4} \end{equation}

Algorithme

Supposons que les k premières lignes de A ont déjà été transformées en forme triangulaire supérieure. Par conséquent, l'équation de pivot actuelle est la kème équation, et toutes les équations en dessous doivent encore être transformées. \begin{equation} \left[ \begin{matrix} A_{11}&A_{12}&A_{13}&\cdots&A_{1k}&\cdots&A_{1j}&\cdots&A_{1n}\\ 0&A_{22}&A_{23}&\cdots&A_{2k}&\cdots&A_{2j}&\cdots&A_{2n}\\ 0&0&A_{33}&\cdots&A_{3k}&\cdots&A_{3j}&\cdots&A_{3n}\\ \vdots&\vdots&\vdots&&\vdots&&\vdots&&\vdots\\ 0&0&0&\cdots&A_{kk}&\cdots&A_{kj}&\cdots&A_{kn}\\ \vdots&\vdots&\vdots&&\vdots&&\vdots&&\vdots\\ 0&0&0&\cdots&A_{ik}&\cdots&A_{ij}&\cdots&A_{in}\\ \vdots&\vdots&\vdots&&\vdots&&\vdots&&\vdots\\ 0&0&0&\cdots&A_{nk}&\cdots&A_{nj}&\cdots&A_{nn}\\ \end{matrix} \left| \, \begin{matrix} b_1 \\ b_2 \\ b_3 \\ \vdots \\ b_k \\ \vdots \\ b_i \\ \vdots \\ b_n \\ \end{matrix} \right. \right] \tag{5} \end{equation}

Soit la ième ligne une ligne typique sous l'équation de pivot qui doit être transformée, ce qui signifie que l'élément \(A_{ik}\) doit être éliminé. Nous pouvons y parvenir en multipliant la ligne pivot par \(\lambda = \frac{A_{ik}} {A_{kk}}\) et en la soustrayant de la ième ligne. \begin{equation} A_{ij} \leftarrow A_{ij} - \lambda A_{kj},\, j=k,k+1,\cdots,n \tag{6} \end{equation} \begin{equation} b_i \leftarrow b_i - \lambda b_k \tag{7} \end{equation} Pour transformer la matrice de coefficients entière en forme triangulaire supérieure, k et i dans les équations. (2 et 3) doit avoir les valeurs \(k = 1, 2,\cdots ,n-1\) (choisit la ligne pivot), \(i = k +1, k + 2,\cdots , n\) (choisit la ligne à transformer).

    # pour chaque pivot
    for k in range(0,n-1):
        # si le pivot égal zéro
        # on cherche un pivot différent de zero dans les équations suivantes
        if A[k,k]==0:
            lpivot=-1 # stocker l'indice du ligne du pivot
            for L in range(k+1,n):
                if A[L,k] !=0:
                    lpivot=L
                    break
                    
            if lpivot!=-1:
                # échange l'équation k avec lpivot
                A[[k, lpivot]] = A[[lpivot, k]]

            # le système n'admit pas de solution
            else:
                return None

        for i in range(k+1,n):
            if A[i,k] != 0.0:
                lam = A[i,k]/A[k,k]
                A[i,k:n+1] = A[i,k:n+1] - lam*A[k,k:n+1]
                                    

Après élimination de Gauss, la matrice de coefficients augmentés a la forme : $$ \left[ A \left| \, b \right. \right] = \left[ \begin{matrix} A_{11}&A_{12}&A_{13}&\cdots&A_{1n}&\\ 0&A_{22}&A_{23}&\cdots&A_{2n}&\\ 0&0&A_{23}&\cdots&A_{3n}&\\ \vdots&\vdots&\vdots&\ddots&\vdots&\\ 0&0&0&\cdots&A_{nn}& \end{matrix} \left| \, \begin{matrix} b_1 \\ b_2 \\ b_3 \\ \vdots \\ b_n \\ \end{matrix} \right. \right] $$ La dernière équation, \(A_{nn}x_n = b_n\), est résolue en premier, ce qui donne : \begin{equation} x_n=b_n / A_{nn} \tag{8} \end{equation}

Phase de substitution

Les inconnues peuvent maintenant être calculées par substitution. Résoudre les équations. (c), (b) et (a) dans cet ordre, nous obtenons : \begin{align*} x_3&=9/3=3\\ x_2&=(-10.5+1.5x_3)/3=(-10.5+1.5\times 3)/3=-2\\ x_1&=(11+2x_2-x_3)/4=(11+2\times (-2)-3)/4=1 \end{align*}

Algorithme

Considérons maintenant l'étape de la substitution, où $x_n, x_{n-1},\cdots , x_{k+1}$ ont déjà été calculés (dans cet ordre), et nous sommes sur le point de déterminer $x_k$ à partir de la kième équation : \begin{align*} A_{k,k}x_{k} + A_{k,k+1}x_{k+1}+\cdots+A_{kn}x_{n}=b_k \\ \end{align*} La solution est : \begin{equation} x_k=\left(b_k-\sum_{j=k+1}^{n} A_{kj}x_j\right)\frac{1}{A_{kk}}, \, \, k=n-1,n-2,\cdots,1 \tag{9} \end{equation}

# phase de substitution
for k in range(n-1,-1,-1):
    # On utilisation le produit scalaire des vecteurs
    # Ici on multiple le vecteur A[k,k+1:n] par x[k+1:n]
    # Fonction np.dot permet de calculer le produit de 2 vecteur
    x[k]=(A[k,n] - np.dot(A[k,k+1:n],x[k+1:n]))/A[k,k]
                                    

Exemple

Tout d'abord, nous créons une fonction qui prend la matrice A et le vecteur b, puis elle renvoie la matrice augmentée

A=[[4,-2,1],[-2,4,-2],[1,-2,4]]
b=[11,-16,17]
def augumentee(A,b):
    n=len(A)
    C=np.zeros((len(b),len(b)+1))
    C[:,:n]=A
    C[:,n]=b
    return C
    def gaussElimin(A):
    n=len(A)
    # pour chaque pivot
    for k in range(0,n-1):
        # si le pivot égal zéro
        # on cherche un pivot différent de zero dans les équations suivantes
        if A[k,k]==0:
            lpivot=-1 # stocker l'indice du ligne du pivot
            for L in range(k+1,n):
                if A[L,k] !=0:
                    lpivot=L
                    break
                    
            if lpivot!=-1:
                # échange l'équation k avec lpivot
                A[[k, lpivot]] = A[[lpivot, k]]

            # le système n'admit pas de solution
            else:
                return None

        for i in range(k+1,n):
            if A[i,k] != 0.0:
                lam = A[i,k]/A[k,k]
                A[i,k:n+1] = A[i,k:n+1] - lam*A[k,k:n+1]
    

    #phase de substitution
    x=[0]*len(A)
    #ann=A[-1,-2]  et bn=A[-1,-1] puisque A est la matrice augumentée
    x[-1]=A[-1,-2]/A[-1,-1]
    for k in range(n-1,-1,-1):
        # On utilisation le produit scalaire des vecteurs
        # Ici on multiple le vecteur A[k,k+1:n] par x[k+1:n]
        # Fonction np.dot permet de calculer le produit de 2 vecteur
        x[k]=(A[k,n] - np.dot(A[k,k+1:n],x[k+1:n]))/A[k,k]


# C c'est la matrice des coefficients augumentées
C=augumentee(A,b)
sol=gaussElimin(C)
if sol is None:
    print('Le système n\'admit pas de solution')
else:
    print(sol)