Descompunerea LU

De la Wikipedia, enciclopedia liberă.
Salt la navigare Salt la căutare

În algebra liniară o descompunere LU , sau descompunerea LUP sau descompunerea Doolittle este o factorizare a unei matrice într-o matrice triunghiulară inferioară , o matrice triunghiulară superioară și o matrice de permutare . Această descompunere este utilizată în analiza numerică pentru a rezolva un sistem de ecuații liniare , pentru a calcula inversul unei matrice sau pentru a calcula determinantul unei matrice.

Definiție

Este o matrice inversabilă . Atunci poate fi descompus ca:

unde este este o matrice de permutare , este o matrice triunghiulară inferioară cu unitatea diagonală ( ) Și este o matrice triunghiulară superioară.

Ideea principală

Descompunerea este similar algoritmului Gauss . În eliminarea gaussiană încercăm să rezolvăm ecuația matricei:

Procesul de eliminare produce o matrice triunghiulară mai înaltă și transformă vectorul în

Atâta timp cât este o matrice triunghiulară mai înaltă, acest sistem de ecuații poate fi ușor rezolvat prin substituție inversă .

În timpul descompunerii , In orice caz, nu este transformat și ecuația poate fi scrisă ca:

deci puteți reutiliza descompunerea dacă doriți să rezolvați același sistem pentru unul diferit .

În cazul mai general, în care factorizarea matricei include și utilizarea schimburilor de rânduri în matrice, este introdusă și o matrice de permutare , iar sistemul devine:

Rezoluția acestui sistem permite determinarea vectorului cautat.

Algoritm

Aplicând serii de transformări ale matricei elementare (adică multiplicări ale stânga) se construiește o matrice triunghiulară superioară care începe de la . Această metodă se numește metoda Gauss . Aceste transformări ale matricei elementare sunt toate transformări liniare de tip combinatorial (al treilea tip listat în „ Matricea elementară ”). Asuma ca fie produsul N transformărilor , atunci matricea triunghiulară superioară este:

Inversul matricei Și:

La fel ca algoritmul Gauss, el folosește doar a treia formă a celor trei tipuri de transformări elementare ale realizării matricei triunghiular superior, se poate deduce că toate sunt triunghiulare inferioare (vezi transformările elementare ale matricei ). Fiind un produs al de asemenea:

este triunghiular inferior.

Matricea este apoi descompusă în produsul de Și :

Aplicații

Matrici inverse

Factorizarea este, de asemenea, utilizat pentru a calcula matricea inversă a . Intr-adevar:

de la care:

Determinant

O altă utilizare a acestei descompuneri este pentru calcularea determinantului matricei . Intr-adevar:

deci prin teorema lui Binet :

Știind că determinantul unei matrici de permutare se menține dacă aceasta corespunde unui număr par de permutări, altfel, și asta , obținem că:

Știind atunci că determinantul unei matrice triunghiulare este dat de produsul termenilor de pe diagonala sa principală , avem că:

dar și știind că termenii de pe diagonala principală a toate sunt , deci putem concluziona în cele din urmă:

unde este indică numărul de schimburi de rânduri efectuate în proces (indicat implicit în matrice ) și termenii Și indicați termenul în linie și coloană respectiv a matricilor Și .

Cod în C

 / * INPUT: A - vector de indicatori către rândurile matricei pătrate de dimensiunea N
* Tol - apelant de toleranță minimă pentru a identifica când matricea este pe cale să degenereze
* IEȘIRE: Matricea A s-a schimbat, conține atât matricile LE, cât și U astfel încât A = (LE) + U și P * A = L * U
* Matricea de permutare nu este stocată în memorie ca o matrice, ci într-un vector
de numere întregi de mărimea N + 1
* conținând indexurile coloanelor în care matricea are ca elemente „1”. Ultimul element P [N] = S + N,
* unde S este numărul de linii schimbate necesar pentru calcularea determinantului, det (P) = (- 1) ^ S
* /
int LUPDecompose ( double ** A , int N , double Tol , int * P ) {

    int i , j , k , imax ; 
    dublu maxA , * ptr , absA ;

    pentru ( i = 0 ; i <= N ; i ++ )
        P [ i ] = i ; // Matrice de permutare unară, P [N] inițializată cu N

    pentru ( i = 0 ; i < N ; i ++ ) {
        maxA = 0,0 ;
        imax = i ;

        pentru ( k = i ; k < N ; k ++ )
            if ((Absa = FABS (A [k] [i]))> Maxa) { 
                maxA = absA ;
                imax = k ;
            }

        if ( maxA < Tol ) returnează 0 ; // Matricea a degenerat

        if ( imax ! = i ) {
            // pivotant P
            j = P [ i ];
            P [ i ] = P [ imax ];
            P [ imax ] = j ;

            // pivotarea liniilor lui A
            ptr = A [ i ];
            A [ i ] = A [ imax ];
            A [ imax ] = ptr ;

            // Numărarea pivoturilor începând de la N pentru calcularea determinantului
            P [ N ] ++ ;
        }

        pentru ( j = i + 1 ; j < N ; j ++ ) {
            A [ j ] [ i ] / = A [ i ] [ i ];

            pentru ( k = i + 1 ; k < N ; k ++ )
                A [ j ] [ k ] - = A [ j ] [ i ] * A [ i ] [ k ];
        }
    }

    retur 1 ; // Descompunerea finalizată
}

/ * INPUT: matrice A, P completate cu funcția LUPDecompose; b - vector; N - dimensiune
* IEȘIRE: x - vectorul soluției lui A * x = b
* /
void LUPSolve ( double ** A , int * P , double * b , int N , double * x ) {

    for ( int i = 0 ; i < N ; i ++ ) {
        x [ i ] = b [ P [ i ]];

        pentru ( int k = 0 ; k < i ; k ++ )
            x [ i ] - = A [ i ] [ k ] * x [ k ];
    }

    pentru ( int i = N - 1 ; i > = 0 ; i - ) {
        pentru ( int k = i + 1 ; k < N ; k ++ )
            x [ i ] - = A [ i ] [ k ] * x [ k ];

        x [ i ] = x [ i ] / A [ i ] [ i ];
    }
}

/ * INPUT: matrice A, P completate cu funcția LUPDecompose; N - dimensiune
* IEȘIRE: IA este inversul matricei inițiale
* /
void LUPInvert ( dublu ** A , int * P , int N , dublu ** IA ) {
  
    for ( int j = 0 ; j < N ; j ++ ) {
        for ( int i = 0 ; i < N ; i ++ ) {
            if ( P [ i ] == j ) 
                IA [ i ] [ j ] = 1,0 ;
            altceva
                IA [ i ] [ j ] = 0,0 ;

            pentru ( int k = 0 ; k < i ; k ++ )
                IA [ i ] [ j ] - = A [ i ] [ k ] * IA [ k ] [ j ];
        }

        pentru ( int i = N - 1 ; i > = 0 ; i - ) {
            pentru ( int k = i + 1 ; k < N ; k ++ )
                IA [ i ] [ j ] - = A [ i ] [ k ] * IA [ k ] [ j ];

            IA [ i ] [ j ] = IA [ i ] [ j ] / A [ i ] [ i ];
        }
    }
}

/ * INPUT: matrice A, P completate cu funcția LUPDecompose; N - dimensiune.
* IEȘIRE: Funcția returnează valoarea determinantului matricei inițiale
* /
dublu LUPDeterminant ( dublu ** A , int * P , int N ) {

    dublu det = A [ 0 ] [ 0 ];

    pentru ( int i = 1 ; i < N ; i ++ )
        det * = A [ i ] [ i ];

    if (( P [ N ] - N ) % 2 == 0 )
        return det ; 
    altceva
        retur - det ;
}

Cod C #

 public class SystemOfLinearEquations
    {
        public double [] SolveUsingLU ( double [,] matrix , double [] rightPart , int n )
        {
            // descompunerea matricei
            double [,] lu = new double [ n , n ];
            suma dublă = 0 ;
            for ( int i = 0 ; i < n ; i ++)
            {
                pentru ( int j = i ; j < n ; j ++)
                {
                    suma = 0 ;
                    pentru ( int k = 0 ; k < i ; k ++)
                        sum + = lu [ i , k ] * lu [ k , j ];
                    lu [ i , j ] = matrice [ i , j ] - sumă ;
                }
                pentru ( int j = i + 1 ; j < n ; j ++)
                {
                    suma = 0 ;
                    pentru ( int k = 0 ; k < i ; k ++)
                        sum + = lu [ j , k ] * lu [ k , i ];
                    lu [ j , i ] = ( 1 / lu [ i , i ]) * ( matrice [ j , i ] - sumă );
                }
            }
            
            // lu = L + UI
            // Calculul soluțiilor lui Ly = b
            double [] y = new double [ n ];
            for ( int i = 0 ; i < n ; i ++)
            {
                suma = 0 ;
                pentru ( int k = 0 ; k < i ; k ++)
                    sum + = lu [ i , k ] * y [ k ];
                y [ i ] = rightPart [ i ] - sumă ;
            }
            // Calculul soluțiilor lui Ux = y
            double [] x = new double [ n ];
            pentru ( int i = n - 1 ; i > = 0 ; i -)
            {
                suma = 0 ;
                pentru ( int k = i + 1 ; k < n ; k ++)
                    sum + = lu [ i , k ] * x [ k ];
                x [ i ] = ( 1 / lu [ i , i ]) * ( y [ i ] - sumă );
            }
            returnează x ;
        }
}

Cod în Matlab

 funcţie X = SolveLinearSystem ( A, b )
    n = lungimea ( b );
    X = zerouri ( n , 1 );
    y = zerouri ( n , 1 );
    % Descompunerea matricei prin metoda Doolittle
    pentru the = 1 : 1 : n
        pentru j = 1 : 1 :( i - 1 )
            alfa = A ( i , j );
            pentru k = 1 : 1 :( j - 1 )
                alfa = alfa - A ( i , k ) * A ( k , j );
            Sfârșit
            A ( i , j ) = alfa / A ( j , j );
        Sfârșit
        pentru j = i : 1 : n
            alfa = A ( i , j );
            pentru k = 1 : 1 :( i - 1 )
                alfa = alfa - A ( i , k ) * A ( k , j );
            Sfârșit
            A ( i , j ) = alfa ;
        Sfârșit
    Sfârșit
    % A = L + UI
    % calculul soluțiilor Ly = b
    pentru the = 1 : 1 : n
        alfa = 0 ;
        pentru k = 1 : 1 : i
            alfa = alfa + A ( i , k ) * y ( k );
        Sfârșit
        y ( i ) = b ( i ) - alfa ;
    Sfârșit
    % calculul soluțiilor lui Ux = y
    pentru the = n :( - 1 ): 1
        alfa = 0 ;
        pentru k = ( i + 1 ): 1 : n
            alfa = alfa + A ( i , k ) * x ( k );
        Sfârșit
        x ( i ) = ( y ( i ) - alfa ) / A ( i , i );
    Sfârșit    
Sfârșit

Bibliografie

  • ( EN ) Bunch, James R.; Hopcroft, John (1974), „Factorizarea și inversarea triunghiulară prin multiplicarea rapidă a matricei”, Matematica calculului 28: 231–236, ISSN 0025-5718.
  • ( EN ) Cormen, Thomas H.; Leiserson, Charles E.; Rivest, Ronald L.; Stein, Clifford (2001), Introducere în algoritmi , MIT Press și McGraw-Hill, ISBN 978-0-262-03293-3 .
  • ( EN ) Golub, Gene H.; Van Loan, Charles F. (1996), Matrix Computations (ediția a treia) , Baltimore: Johns Hopkins, ISBN 978-0-8018-5414-9 .
  • (EN) Horn, Roger A.; Johnson, Charles R. (1985), Matrix Analysis , Cambridge University Press, ISBN 0-521-38632-2 . A se vedea secțiunea 3.5.
  • (EN) Householder, Alston S. (1975), Theory of Matrices in Numerical Analysis, New York: Dover Publications, MR 0378371.
  • ( EN ) Okunev, Pavel; Johnson, Charles R. (1997), Condiții necesare și suficiente pentru existența factorizării LU a unei matrice arbitrare , arXiv: math.NA/0506382.
  • ( EN ) Poole, David (2006), Algebra liniară: o introducere modernă (ediția a doua) , Canada: Thomson Brooks / Cole, ISBN 0-534-99845-3 .
  • ( EN ) Apăsați, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (2007), „Secțiunea 2.3”, Rețete numerice: arta calculelor științifice (ediția a 3-a) , New York: Cambridge University Press, ISBN 978-0-521-88068-8 .
  • (EN) Trefethen, Lloyd N.; Bau, David (1997), Algebra liniară numerică , Philadelphia: Society for Industrial and Applied Mathematics, ISBN 978-0-89871-361-9 .

Elemente conexe

linkuri externe

Matematica Portalul de matematică : accesați intrările Wikipedia care se ocupă de matematică