Metode Numerice - Curs 3 PDF

METODE NUMERICE – curs 3 Cap. 2 Sisteme determinate de ecuaţii algebrice liniare 2.1 Formularea problemei  Metodele numerice de rezolvare a unui sistem determinat de ecuaţii algebrice liniare metode directe: metode iterative...

METODE NUMERICE – curs 3 Cap. 2 Sisteme determinate de ecuaţii algebrice liniare 2.1 Formularea problemei  Metodele numerice de rezolvare a unui sistem determinat de ecuaţii algebrice liniare metode directe: metode iterative (indirecte):  sistem echivalent, direct rezolvabil prin mijloace  construirea unui şir de aproximaţii pentru elementare  eliminarea progresivă a soluţia sistemului, convergent la aceasta: necunoscutelor (eliminare gaussiană) 𝑘→∞ 𝑥 [𝑘] 𝑥  prin transformări elementare de echivalenţă, se 𝑘≥0 aduce matricea A la anumite forme tipice: u 11   u 1n   l11 0  0  lim || x  x ||  0 [k]  0         k  U L          0      calculele se opresc la un index de iterare  0  0 u nn   n1 l   l nn  [s], când este îndeplinită condiţia: [s 1] x ||   impus [s ] matrice superior matrice inferior || x triunghiulară triunghiulară  procedura de transformare se numeşte triangularizare  METODE NUMERICE – curs 3 2.2 Rezolvarea sistemelor prin triangularizare directă 2.2.1 Principiul metodei  se consideră sistemul de n ecuaţii algebrice liniare cu n necunoscute: 𝐴 ∙ 𝑥 = 𝑏, 𝐴 ∈ ℝ𝑛×𝑛 , 𝑏 ∈ ℝ𝑛×1  notaţie: A  [a ij ]1i n 1 j n  matricile: [a 11 ] , …, [a ij ]1i n 1 - submatrici principale (minori principali directori) ale lui A 1 j n 1 Teoremă: Dacă matricea 𝐴 ∈ ℝ𝑛×𝑛 are toate submatricele principale inversabile (nesingulare), atunci există matricele 𝐿, 𝐷, 𝑈 ∈ ℝ𝑛×𝑛 astfel încât: A LDU unde L este o matrice inferior triunghiulară, D este o matrice diagonală şi U este o matrice superior triunghiulară.  METODE NUMERICE – curs 3 Observaţii: 1. uzuale sunt factorizările: 𝐴 = 𝐿 ∙ 𝑈 (factorizare L-U); 2. demonstraţia teoremei enunţate este constructivă, constituind însuşi algoritmul de descompunere L-U a matricei A; 3. algoritmul de descompunere  procedeul de eliminare gaussiană a necunoscutelor - matricea A este adusă la forma superior triunghiulară în urma unui şir de transformări de asemănare  se “acumulează” într-o matrice inferior triunghiulară, cu elementele de pe diagonala principală egale cu 1 (descompunere Doolittle); 4. considerând descompunerea L-U a matricei sistemului A, rezolvarea sistemului implică două subetape: a) substituţie înainte: rezolvarea sistemului 𝐿 ∙ 𝑦 = 𝑏 - din aproape în aproape: 𝑦1 din prima ecuaţie  𝑦2 din a doua ecuaţie  … 𝑦𝑛 din ultima ecuaţie (𝑦 = 𝑦𝑖 𝑖=1,𝑛 ) b) substituţie inversă: rezolvarea sistemului U∙ 𝑥 = 𝑦 - din aproape în aproape: 𝑥𝑛 din ultima ecuaţie  𝑥𝑛−1 din penultima ecuaţie  … 𝑥1 din prima ecuaţie (𝑥 = 𝑥𝑖 𝑖=1,𝑛 )  METODE NUMERICE – curs 3 Propoziţie: Dacă matricea A admite o descompunere L-U, atunci această descompunere este unică.  procedura de triangularizare directă necesită un număr de operaţii în virgulă mobilă de 3 ordinul lui 𝑛 ൗ3;  dacă matricea A este simetrică (𝐴 = 𝐴𝑇 ) şi pozitiv definită (∀𝑥 ∈ ℝ𝑛×1 , 𝑥 ≠ 0𝑛×1 , 𝑥 𝑇 ∙ 𝐴 ∙ 𝑥 > 0 şi  𝑥 𝑇 ∙ 𝐴 ∙ 𝑥 = 0 ⇔ 𝑥 = 0𝑛×1 ), A se descompune sub forma 𝐴 = 𝐿 ∙ 𝐿𝑇 (descompunerea Cholesky). 2.2.2 Procedura de triangularizare directă a unei matrice  algoritmul triangularizării directe: atribuie 𝐴1 ⟵ 𝐴 pentru 𝑘 = 1, 𝑛 − 1 execută  determinare matrice 𝑀𝑘 astfel încât 𝐴𝑘+1 = 𝑀𝑘 ∙ 𝐴𝑘 să aibă elementele: [𝑘+1] [𝑘+1] [𝑘] 𝑎𝑖,𝑘 = 0, 𝑖 = 𝑘 + 1, … , 𝑛 şi 𝑎𝑖,𝑗 = 𝑎𝑖,𝑗 , 𝑖 = 1, … , 𝑛; 𝑗 = 1, … , 𝑘 − 1 atribuie 𝐴𝑘+1 ⟵ 𝑀𝑘 ∙ 𝐴𝑘  atribuie 𝑈 ⟵ 𝐴𝑛  METODE NUMERICE – curs 3  algoritmul parcurge (n – 1) etape în care se zerorizează elementele unei coloane, situate sub diagonala principală, lăsând nemodificate coloanele transformate la etapele anterioare [𝑘]  notaţie: 𝜉 = 𝑐𝑘 𝐴𝑘 ; 𝑐𝑘 𝐴𝑘 - coloana k a matricii 𝐴𝑘 = 𝑎𝑖,𝑗 𝑖,𝑗=1,…,𝑛   [1   k  k 1   n ]T  [a 1[ k,k]  a [kk,k] a [kk]1,k  a [nk,k] ]T 𝑇  𝑀𝑘 se construieşte încât 𝑀𝑘 ∙ 𝜉 = 𝜉1 ⋯ 𝜉𝑘 0 ⋯ 0  se consideră vectorul de multiplicatori: 𝑚𝑘 = 0 ⋯ 0 𝜇𝑘+1,𝑘 ⋯ 𝜇𝑛,𝑘 𝑇 Definiţie: Matricea 𝑀𝑘 se numeşte matrice de transformare elementară de ordin n şi indice k sau matrice Gauss şi este definită prin: Mk  In  mk  ek T unde 𝑒𝑘 este coloana k a matricii unitate de ordin n, 𝐼𝑛. CALCUL NUMERIC – curs 3  proprietăţi ale matricilor Gauss: sunt inferior triunghiulare; elementele diagonalei principale sunt egale cu 1; sunt nesingulare, având inversa: M k1  I n  m k  e k T  efectul aplicării matricei 𝑀𝑘 asupra vectorului 𝜉: M k    (I n  m k  e k )      m k  e k      m k   k T T  [ 1   k  k 1   k 1,k   k   n   n ,k   k ]T trebuie zerorizate  presupunând 𝜉𝑘 ≠ 0, se aleg multiplicatorii Gauss:  i,k   i /  k , i  k  1,..., n pivot M k    [1   k 0  0]T  METODE NUMERICE – curs 3 Observaţii: 1. În practică, pe calculator, etapa k descrisă mai sus se poate realiza testând condiţia: | a [kk,k] |  [𝑘] în loc de a verifica 𝑎𝑘,𝑘 ≠ 0, unde 𝜀 este o constantă impusă, de valoare mică sau foarte mica (de exemplu, epsilonul-maşină); aceasta se realizează datorită faptului că, dacă în aritmetica exactă pivotul este nul, în aritmetica virgulei mobile, datorită erorilor de calcul, această situaţie este echivalentă cu: | a k ,k |  ; [k] 2. Când pivotul este în modul mai mic sau egal cu 𝜀, eliminarea gaussiană eşuează; aceasta corespunde situaţiei când matricea iniţială A are submatricea principală de ordin k singulară, deci conform teormei enunţate anterior, descompunerea L-U a matricei A nu există.  efectul aplicării transformării 𝑀𝑘 asupra celorlaltor coloane ale matricei 𝐴𝑘 : 𝑇 - se consideră un vector 𝜂 ≠ 𝜉, cu 𝜂 = 𝜂1 ⋯ 𝜂𝑛 M k    [1  k k 1   k 1,k  k  n   n ,k  k ]T METODE NUMERICE – curs 3 Concluzii: 1. Matricea 𝑀𝑘 lasă nemodificate primele k – 1 coloane ale matricei 𝐴𝑘 :   c j (A k )  [*  * 0  0  0]T , j  k, j  1,..., k  1 1 j j1 k n M k  c j (A k )  [*  * 0  0 0   k 1,k  0  0   n ,k  0]T 1 j j1 k k 1 n  c j (A k ) 2. Matricea 𝑀𝑘 transformă coloana k a matricei 𝐴𝑘 zerorizând liniile 𝑘 + 1, ⋯ , 𝑛; 3. Matricea 𝑀𝑘 transformă coloanele 𝑘 + 1, ⋯ , 𝑛 ale matricei 𝐴𝑘 în liniile 𝑘 + 1, ⋯ , 𝑛:   c j (A k ), j  k  1,..., n M k  c j (A k )  [*  * a [kk]1, j   k 1,k  a [kk, ]j  a [nk, ]j   n,k  a [kk, ]j ]T (notaţia * semnificând faptul că elementele implicate rămân nemodificate) METODE NUMERICE – curs 3  tabloul general al transformărilor: M n 1    M 2  M 1  A  U A  M11  M 21    M n11  U A LU L n 1 L  M 11  M 21    M n11  I n   m k  e k T k 1 Observaţie: 1. Matricea L este inferior triunghiulară unitate şi conţine în fiecare coloană, sub elementul unitar de pe diagonala principală, subvectorii Gauss; 2. Prima sub-etapă de rezolvare a sistemului 𝐴 ∙ 𝑥 = 𝑏 este substituţia înainte aplicată sistemului de ecuaţii 𝐿 ∙ 𝑦 = 𝑏 : y  L1  b  M n 1    M 2  M1  b  METODE NUMERICE – curs 3 Concluzie: La triangularizarea simplă, unde elementele matricei A se modifică corespunzător relaţiei: a [ijk 1]  a [ijk ]   ik  a [kjk ] , i  k  1,..., n; j  k,..., n multiplicatorii 𝜇𝑖𝑘 pot avea, în principiu, orice valoare - dacă aceste valori sunt mari sau foarte mari, atunci: pot apare fenomenele de omitere catastrofală şi/sau de neutralizare a termenilor; [𝑘] amplifică erorile prezente în termenii 𝑎𝑘𝑗. triangularizarea simplă este instabilă numeric  stabilizarea algoritmului triangularizării unei matrici  multiplicatori Gauss cu modul subunitar  METODE NUMERICE – curs 3 𝟐 −3 5 Exemplu: 𝐴 = 𝐴1 = 4 7 −3 1 9 4 iteraţia 1: 0 0 1 0 0 𝟐 −3 5 𝑚1 = 4Τ𝟐 = 2 𝑀1 = −2 1 0 𝐴2 = 𝑀1 ∙ 𝐴1 = 0 𝟏𝟑 −13 1Τ𝟐 0.5 −0.5 0 1 0 10.5 1.5 iteraţia 2: 0 0 1 0 0 𝑚2 = 0 = 0 𝑀2 = 0 1 0 10.5Τ𝟏𝟑 0.8077 0 −0.8077 1 𝟐 −3 5 𝐴3 = 𝑀2 ∙ 𝐴2 = 0 𝟏𝟑 −13 = 𝑈 0 0 12 1 0 0 𝐿= 𝑀1−1 ∙ 𝑀2−1 = 2 1 0 0.5 0.8077 1 METODE NUMERICE – curs 3 2.3 Rezolvarea sistemelor prin triangularizare cu pivotare parţială  principiu: la pasul k se caută pivotul printre elementele din coloana k, pornind de la elementul de pe diagonala principală în jos, alegându-se elementul care are cea mai mare valoare în modul: | a [i kk ],k | max{| a [i,kk] |}   k k i  n k Ak = 0  ik n | k  dacă 𝑖𝑘 ≠ 𝑘  se permută liniile k şi 𝑖𝑘 cu ajutorul matricii de permutare de linii 𝑃𝑘 Pk  A k A k 1  M k  (Pk  A k ) |  i,k | 1, i  k  1,..., n METODE NUMERICE – curs 3  proprietăţi matrice de permutare de linii, Pk: det(Pk )  1; Pk  Pk1 Teoremă: Dacă matricea A este nesingulară, atunci există o matrice numită matrice generală de permutare de linii, astfel încât: P  A  L'  U în care U este o matrice superior triunghiulară şi Lˈ este o matrice inferior triunghiulară unitate cu elementele | l i, j | 1, i  j  Demonstraţia  constructivă, constituind însuşi algoritmul de triangularizare cu pivotare parţială a unei matrici METODE NUMERICE – curs 3  algoritmul de triangularizare cu pivotare parţială: atribuie 𝐴1 ⟵ 𝐴 pentru 𝑘 = 1, 𝑛 − 1 execută [𝑘] [𝑘] [𝑘] determinare 𝜉𝑘 ← 𝑎𝑖𝑘 ,𝑘 unde 𝑎𝑖𝑘 ,𝑘 = max 𝑎𝑖,𝑘 𝑖=𝑘,𝑛 dacă (𝑖𝑘 ≠ 𝑘) atunci determină 𝑃𝑘 altfel atribuie 𝑃𝑘 ← 𝐼𝑛  calcul 𝑃𝑘 ∙ 𝐴𝑘 determinare matrice 𝑀𝑘 astfel încât 𝐴𝑘+1 = 𝑀𝑘 ∙ 𝑃𝑘 ∙ 𝐴𝑘 să aibă elementele: [𝑘+1] [𝑘+1] [𝑘] 𝑎𝑖,𝑘 = 0, 𝑖 = 𝑘 + 1, … , 𝑛 şi 𝑎𝑖,𝑗 = 𝑎𝑖,𝑗 , 𝑖 = 1, … , 𝑛; 𝑗 = 1, … , 𝑘 − 1 atribuie 𝐴𝑘+1 ⟵ 𝑀𝑘 ∙ 𝑃𝑘 ∙ 𝐴𝑘  atribuie 𝑈 ⟵ 𝐴𝑛  METODE NUMERICE – curs 3  tabloul general al transformărilor: M n 1  Pn 1    M 2  P2  M1  P1  A  U - folosind 𝑃𝑘 = 𝑃𝑘−1 (M n 1  Pn 1    M1  P1  P1  P2    Pn 1 )  (Pn 1    P2  P1 )  A  U −1 P 𝐿′ matrice generală de permutare de linii P  A  L'  U L'  Pn 1    P2  M11  P2  M 21    Pn 1  M n11 matrice inferior triunghiulară unitate având în fiecare coloană, sub diagonala principală, subvectori Gauss cu liniile permutate METODE NUMERICE – curs 3  etapele rezolvării sistemului 𝐴 ∙ 𝑥 = 𝑏, A ∈ ℝ𝑛×𝑛 , 𝑏 ∈ ℝ𝑛×1: a) factorizarea L_U a matricii A prin triangularizare cu pivotare parţială: P  A  L'  U b) calculul vectorului: 𝑐 = 𝑃 ∙ 𝑏 c) rezolvare sistem 𝐿′ ∙ 𝑦 = 𝑐 prin substituţie înainte d) Rezolvare sistem 𝑈 ∙ 𝑥 = 𝑦 prin substituţie inversă PAx Pb L'  U  x  P  b P  A  L'  U cPb L'  y  c yUx  Observaţie: Dacă algoritmul de triangularizare cu pivotare parţială eşuează, în sensul că pivotul găsit la o anumită etapă [k] este nul sau foarte mic în modul, aceasta corespunde situaţiei când în aritmetica reală primele k coloane ale matricei A sunt liniar dependente.