Metode Numerice - Curs de PDF
Document Details
Uploaded by Deleted User
Tags
Summary
Acest document conține un curs despre Metode Numerice, acoperind subiecte cum ar fi calculul în virgulă mobilă, rezolvarea sistemelor de ecuații algebrice liniare și neliniare, precum și aproximarea numerică a funcţiilor. Explică concepte şi formule matematice şi conţine exemple.
Full Transcript
METODE NUMERICE – curs 1 Notare: Nota_finală = 0.3 * Nota_laborator + 0.7 * Nota_examen - Realizarea sarcinilor de lucru - Pertinenţa răspunsurilor - Teme (săptămânile 4, 8 și 11) Atenţie!!! Nota_laborator, Nota_examen – mini...
METODE NUMERICE – curs 1 Notare: Nota_finală = 0.3 * Nota_laborator + 0.7 * Nota_examen - Realizarea sarcinilor de lucru - Pertinenţa răspunsurilor - Teme (săptămânile 4, 8 și 11) Atenţie!!! Nota_laborator, Nota_examen – minim 4.50 Condiţie intrare în examen: maxim 1 absenţă la laborator 1 absenţă poate fi recuperată în ultima săptămână METODE NUMERICE – curs 1 Structurarea materiei: Cap. 1 Calculul în virgulă mobilă Cap.2 Rezolvarea sistemelor determinate de ecuaţii algebrice liniare Cap.3 Rezolvarea sistemelor supradeterminate de ecuaţii algebrice liniare Cap.4 Calculul valorilor şi vectorilor proprii Cap.5 Calculul valorilor singulare Cap.6 Rezolvarea ecuaţiilor şi sistemelor de ecuaţii neliniare Cap.7 Aproximarea numerică a funcţiilor Cap.8 Rezolvarea ecuaţiilor diferenţiale METODE NUMERICE – curs 1 Cap. 1 Calculul în virgulă mobilă 1.1 Aritmetica în virgulă mobilă - reprezentarea numerelor în calculatorul numeric depinde de: tipul numerelor (întregi sau reale); structura constructivă a echipamentului de calcul: baza de reprezentare a numerelor; lungimea cuvântului de memorie număr finit de cifre Numerele întregi reprezentabile: - formează o mulţime finită, I; - reprezentarea lor este exacta; - aritmetica cu aceste numere este exactă (cu excepţia operaţiei de împărţire, în general) 𝐼 = 𝑧 ∈ ℤΤ𝑚𝐼 ≤ 𝑧 ≤ 𝑀𝐼 mI, MI – depinde de baza de numeraţie, lungimea cuvântului de memorie, precum şi modul de reprezentare METODE NUMERICE – curs 1 α – baza de numeraţie; tI – număr cifre în baza α ce pot fi reprezentate; 𝑚𝐼 = −𝛼 𝑡𝐼−1 𝑀𝐼 = 𝛼 𝑡𝐼−1 − 1 Exemplu: 𝛼 = 2, 𝑡𝐼 = 16 ⇒ 𝑚𝐼 = −215 ; 𝑀𝐼 = 215 − 1 𝛼 = 2, 𝑡𝐼 = 32 ⇒ 𝑚𝐼 = −231 ; 𝑀𝐼 = 231 − 1 Observaţii: Încercarea de a opera cu numere care nu aparţin domeniului I de reprezentare, determină, la majoritatea calculatoarelor numerice, emiterea unor mesaje de eroare fatală, programele implicate fiind abandonate: “depăşire (binară) inferioară” (dacă z < mI ), respectiv “depăşire (binară) superioară” (dacă z > MI ). uzuală este reprezentarea în baza de numeraţie 2, alocându-se o cifră pentru reprezentarea semnului şi tI – 1 cifre pentru reprezentarea valorii absolute a numărului reprezentare în cod complementar faţă de baza α METODE NUMERICE – curs 1 reprezentare în cod complementar faţă de baza α 𝑧, 𝑧≥0 𝑧𝑐 = ቊ 𝑡𝐼 𝛼 − 𝑧, 𝑧 0 şi semnul “-“ se consideră pentru f < 0 METODE NUMERICE – curs 2 x y c1 c2 c3 F Fig. 2 Principiul rotunjirii simetrice: fl( x ) c 2 , fl( y) c 2 ; x ((c1 c 2 ) / 2, c 2 ), y (c 2 , (c 2 c 3 ) / 2) Observaţie: Neajunsul acestei maniere de reprezentare dacă numărul x este situat la jumătatea distanţei dintre două numere consecutive din F, atunci fl(x) poate lua oricare din cele două valori învecinate. rotunjirea uniformă: f e , | g | 0, 2 f e et , | g | 0, fl( x ) 2 f e et , | g | 0, , ultima cifra f impara 2 f e , | g | 0, , ultima cifra f para 2 - semnul “+” se consideră pentru f > 0 şi semnul “-“ se consideră pentru f < 0; - este adoptată de standardul IEEE METODE NUMERICE – curs 2 Exemplul 1: Se consideră un sistem de numere în virgulă mobilă cu 𝛽 = 10, 𝑡 = 4 şi rotunjire uniformă: 𝑔 < 0.5 54 x 12944 ,9942 0,1294 10 0,4994 10 5 fl(x) 0,1294 10 5 ( 12940 x) 𝑔 > 0.5 x 129551 0,1295 10 6 0,51 10 64 fl(x ) 0,1295 10 6 10 64 ( 129600 x) 𝑔 = 0.5 x 1297 ,5 0,1297 10 4 0,5 10 44 fl(x) 0,1297 10 4 10 44 ( 1298 x) u.c.(f) - impară 𝑔 = 0.5 4 4 x 1296,5 0,1296 10 0,5 10 4 fl(x) 0,1296 10 4 ( 1296 x) u.c.(f) - pară e x | x fl( x ) | x e x / | x | e x / | fl( x ) || x fl( x ) | / | fl( x ) | eroare absolută eroare relativă METODE NUMERICE – curs 2 | x fl( x ) | | f e g e t fl( x ) | 1 et x k k 1 t | fl( x ) | | fl( x ) | (1 / ) e - rotunjirea prin trunchiere 𝑘 = 1 ⇒ 𝜀𝑥 = 𝛽1−𝑡 - cea mai mare spaţiere relativă - rotunjirea simetrică 𝑘 = 1Τ𝛽 ⇒ 𝜀𝑥 = 𝛽 −𝑡 - cea mai mică spaţiere relativă 1.1.3 Operaţii elementare în virgulă mobilă - se definesc operaţiile elementare care au loc cu elementele unei mulţimi F de numere în virgulă mobilă Adunarea/ scăderea - oricare ar fi două numere x şi y din mulţimea G, pentru care există fl(x) şi fl(y) aparţinând mulţimii F, numărului x + y i se asociază fl(x + y) conform algoritmului: METODE NUMERICE – curs 2 Pas 1: se reprezintă intern numerele x şi y prin fl(x) şi, respectiv, fl(y): 𝑓𝑙 𝑥 = 𝑓𝑥 ∙ 𝛽 𝑒𝑥 , 𝑓𝑙 𝑦 = 𝑓𝑦 ∙ 𝛽 𝑒𝑦 Pas 2: dacă 𝑒𝑥 ≠ 𝑒𝑦 numărul cu exponent mai mic se aduce la o formă în care exponentul să fie egal cu cel al celuilalt termen denormalizare deplasarea mantisei spre dreapta, inserând zerouri după virgulă Pas 3: se adună mantisele şi se normalizează rezultatul (dacă este necesar) Pas 4: din rezultat se păstrează t cifre Observaţii: Deplasarea mantisei spre dreapta cu o poziţie determină creşterea exponentului cu o unitate, iar deplasarea spre stânga a mantisei cu o poziţie determină scăderea exponentului cu o unitate. Scăderea se realizează la fel ca adunarea, cu deosebirea că mantisele se scad (scăderea reprezintă o adunare în care scăzătorul are semn schimbat). Adunarea/ scăderea în virgulă mobilă nu este asociativă METODE NUMERICE – curs 2 erori determinate de adunarea/ scăderea în virgulă mobile: omiterea catastrofală apare atunci când se adună doi termeni şi valoarea absolută a unui termen este mai mică decât precizia de reprezentare a celuilalt termen; în acest caz, rezultatul este dat de termenul cu valoare absolută mai mare; neutralizarea termenilor apare atunci când se adună numere cu semne diferite şi cu valori absolute apropiate; în acest caz, în mod eronat, rezultatul este nul. precizia calculelor numerice caracterizată de două mărimi constante (“constante de maşină”): Definiţie: + Epsilonul maşină pentru adunare (notat 𝜀𝑚 ) - cel mai mic număr real pozitiv reprezentabil care + schimbă, prin adunare, unitatea maşinii de calcul: 𝑓𝑙 1 + 𝜀𝑚 >1 − Epsilonul maşină pentru scădere (notat 𝜀𝑚 ) - cel mai mare număr real pozitiv reprezentabil − care schimbă, prin scădere, unitatea maşinii de calcul: 𝑓𝑙 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 j1 k n M k c j (A k ) [* * 0 0 0 k 1,k 0 0 n ,k 0]T 1 j j1 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 M11 M 21 M n11 U A LU L n 1 L M 11 M 21 M n11 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 L1 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 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 Fig. 2.1 Principiul triangularizării cu pivotare parţială a unei matrice: pivotul se găseşte în coloana k, liniile k n; k = 1, …, n - 1 dacă 𝑖𝑘 ≠ 𝑘 se permută liniile k şi 𝑖𝑘 cu ajutorul matricii de permutare de linii 𝑃𝑘 Pk A k A k 1 M k (Pk A k ) METODE NUMERICE – curs 4 Cap. 2 Sisteme determinate de ecuaţii algebrice liniare 2.3 Rezolvarea sistemelor prin triangularizare cu pivotare parţială se consideră sistemul de n ecuaţii algebrice liniare cu n necunoscute: 𝐴 ∙ 𝑥 = 𝑏, 𝐴 ∈ ℝ𝑛×𝑛 , 𝑏 ∈ ℝ𝑛×1 problema de calcul determinarea unei soluţii 𝑥 ∈ ℝ𝑛×1 care să verifice sistemul dat 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 k ],k | max{| a [i,kk] |} k k k i n dacă 𝑖𝑘 ≠ 𝑘 se permută liniile k şi 𝑖𝑘 cu ajutorul matricii de permutare de linii 𝑃𝑘 k Ak = 0 Pk A k A k 1 M k (Pk A k ) ik n | k | i,k | 1, i k 1,..., n Fig. 2.1 Principiul triangularizării cu pivotare parţială a unei matrice: pivotul se găseşte în coloana k, liniile k n; k = 1, …, n - 1 METODE NUMERICE – curs 4 proprietăţi matrice de permutare de linii, Pk: det(Pk ) 1; Pk Pk1 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 4 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 4 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 M11 P2 M 21 Pn 1 M n11 matrice inferior triunghiulară unitate având în fiecare coloană, sub diagonala principală, subvectori Gauss cu liniile permutate METODE NUMERICE – curs 4 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ă PAx Pb L' U x P b P A L' U cPb L' y c yUx 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. METODE NUMERICE – curs 4 2.4 Aplicaţii ale descompunerilor L-U 2.4.1 Calculul determinantului fie o matrice 𝐴 ∈ ℝ𝑛×𝑛 det 𝐴 =? a) calcul bazat pe factorizarea L-U prin triangularizare directă: 𝐴 = 𝐿 ∙ 𝑈 det(A) det(L U) det(L) det(U) L – inferior triunghiulară unitate ⇒ det 𝐿 = 1 U – superior triunghiulară ⇒ det 𝑈 = ς𝑛𝑖=1 𝑢𝑖,𝑖 , 𝑈 = 𝑢𝑖,𝑖 𝑖=1,…,𝑛 n det(A) u i ,i i 1 METODE NUMERICE – curs 4 b) calcul bazat pe factorizarea L-U prin triangularizare cu pivotare parţială: P ∙ 𝐴 = 𝐿′ ∙ 𝑈 A P 1 L' U det(A) det(P 1 L' U) det(P 1 ) det(L' ) det(U) 𝑛 𝑛 det 𝐴 = −1 𝑛𝑝ℓ ∙ ෑ 𝑢𝑖,𝑖 (−1)𝑛𝑝ℓ 1 ෑ 𝑢𝑖,𝑖 𝑖=1 𝑖=1 2.4.2 Calculul inversei unei matrici fie o matrice 𝐴 ∈ ℝ𝑛×𝑛 𝐴−1 =? acest tip de problemă se încadrează în problematica rezolvării ecuaţiilor matriciale de tipul: 𝐴 ∙ 𝑋 = 𝐵, în care 𝐵 = 𝐼𝑛 A 1 X x 1 x k xn METODE NUMERICE – curs 4 a) calcul bazat pe factorizarea L-U prin triangularizare direct a.1) factorizarea 𝐴 = 𝐿 ∙ 𝑈; k a.2) rezolvare sisteme: 𝐴 ∙ 𝑥𝑘 = 𝑒𝑘 , 𝑒𝑘 = 0 ⋯ 0 1 0 ⋯ 0 𝑇 , 𝑘 = 1, 𝑛 - substituţia înainte: 𝐿 ∙ 𝑦𝑘 = 𝑒𝑘 - substituţia inversă: 𝑈 ∙ 𝑥𝑘 = 𝑦𝑘 b) calcul bazat pe triangularizarea cu pivotare parţială: b.1) factorizarea 𝑃 ∙ 𝐴 = 𝐿′ ∙ 𝑈; k b.2) rezolvare sisteme: 𝐴 ∙ 𝑥𝑘 = 𝑒𝑘 , 𝑒𝑘 = 0 ⋯ 0 1 0 ⋯ 0 𝑇 , 𝑘 = 1, 𝑛 - calcul 𝑐 = 𝑃 ∙ 𝑒𝑘 ; - substituţia înainte: 𝐿′ ∙ 𝑦𝑘 = 𝑐; - substituţia inversă: 𝑈 ∙ 𝑥𝑘 = 𝑦𝑘 Concluzie: Inversarea unei matrice necesită un număr mare de operaţii în virgulă mobile în practică, nu se recomandă rezolvarea sistemelor prin metoda bazată pe calculul explicit al inversei matricei sistemului. METODE NUMERICE – curs 4 2.4 Metode iterative 2.4.1 Principiul şi convergenţa metodelor iterative A x b, A nn , b n1 A - matrice nesingulară Metodele iterative construcţia unui şir de vectori convergent la soluţia sistemului: x , x [k] , k 0,1,... [k] x lim x k Relaţia de recurenţă: ANP N - matrice nesingulară ( N P) x b N x P x b Nx [ k 1] Px [k] b, k 0,1,... [ k 1] N 1 P x N 1 b, k 0,1,... [k] x Notaţie: G N 1 P, G nn METODE NUMERICE – curs 4 G are valori proprii în general complexe, care formează mulţimea numită spectrul matricei G (G ) max{| i (G ) |} - rază spectrală a matricei G 1i n Teoremă: Condiţia necesară şi suficientă ca şirul de vectori să fie convergent către soluţia sistemului de ecuaţii este ca matricea G să aibă toate valorile proprii în modul subunitare sau, altfel spus, raza spectrală a matricei G să fie subunitară. Observaţii: Cu cât raza spectrală subunitară a matricei G este mai mică, cu atât viteza de convergenţă a şirului de vectori este mai mare În practică, de multe ori, condiţia necesară şi suficientă prezentată anterior se înlocuieşte printr-o condiţie suficientă, dacă este posibil, şi anume: dacă || G || 1 atunci (G ) 1 norma matricială infinit n || G || max{ | g i , j |} 1 1i n j1 METODE NUMERICE – curs 4 Se consideră: ALDU L - inferior triunghiulară cu diagonala principală nulă D - diagonală având elementele de pe diagonala principală egale cu elementele de pe diagonala principală a matricei A (se presupune nesingulară) U - superior triunghiulară cu diagonala principală nulă 2.4.2 Metoda Jacobi şi metoda Gauss-Seidel metoda Jacobi N D, P (L U) ⇒ D x [ k 1] (L U) x [ k ] b, k 0,1,... ⇕ n a i ,i x [i k 1] a i , j x [jk ] b i , i 1,..., n j1 a i,i 0 j i ⇙ n x [i k 1] (b i / a i ,i ) (a i , j / a i ,i ) x [jk ] , i 1,..., n j1 j i METODE NUMERICE – curs 4 1 1 0, i j G Jacobi N P D (L U) [g i , j ]1i , j n g i, j a i , j / a i ,i , i j Condiţia suficientă care se impune pentru ca metoda Jacobi să fie convergentă este: n n max{ | g i , j |} max{ | a i , j / a i ,i |} 1 1i n j1 1i n j1 j i ⇓ n n | a i , j / a i ,i | 1 | a i ,i | | a i , j | j1 j1 j i j i A - matrice diagonal dominantă pe linii Propoziţie: Dacă matricea A este diagonal dominantă pe linii, atunci metoda Jacobi este convergentă, oricare ar fi estimaţia iniţială a soluţiei sistemului de ecuaţii. CALCULNUMERICE METODE UMERIC ––curs curs44 Observaţie: | a i , j / a i ,i | 1 i 1,..., n , j 1,..., n , i j sunt coeficienţii care multiplică componentele calculate la iteraţii anterioare, deci erorile ce le afectează sunt micşorate pe măsură ce procesul iterativ avansează A - diagonal dominantă pe linii ⇒ procedura este sigur stabilă numeric metoda Gauss-Seidel N L D, P U ⇒ (L D) x [ k 1] U x [ k ] b, k 0,1,... ⇕ i n a i , j x [jk 1] a i , j x [jk ] b i , i 1,..., n j1 ji 1 a i,i 0 ⇙ 𝑖−1 𝑛 [𝑘+1] 𝑏𝑖 𝑎𝑖,𝑗 𝑘+1 [𝑘] 𝑥𝑖 = − ⋅ 𝑥𝑗 − (𝑎𝑖,𝑗 /𝑎𝑖,𝑖 ) ⋅ 𝑥𝑗 , 𝑖 = 1,... , 𝑛 𝑎𝑖,𝑖 𝑎𝑖,𝑖 𝑗=1 𝑗=𝑖+1 METODE NUMERICE – curs 4 Observaţie: Şi în cazul metodei Gauss-Seidel, dacă matricea A este diagonal dominantă pe linii, atunci metoda este convergentă Între razele spectrale subunitare corespunzătoare celor două metode există relaţia: 2 (G Jacobi ) (G Gauss Seidel ) 1 metoda Gauss-Seidel are o viteză de convergenţă mai mare METODE NUMERIC – curs 5 Cap. 3 Rezolvarea numerică a sistemelor supradeterminate de ecuaţii algebrice liniare în sensul celor mai mici pătrate 3.1 Formularea problemei A x b, A mn , b m1 , x n1 mn - sistem supradeterminat Definiţii: m n A are coloanele liniar independente (monică), dacă vectorii coloană ai săi sunt liniar independenţi: n j a j 0m j 0, j 1, , n j1 subspaţiu imagine al matricei A: Im(A) {y m1 / y A x, x n1 } m1 subspaţiu nul (nucleu) al matricei A: N ( A ) {x n1 / A x 0 m } n1 METODE NUMERIC – curs 5 Propoziţie: m n Pentru orice matrice A , m n , următoarele afirmaţii sunt echivalente: A este o matrice monică; rang(A ) n ; N(A) {0 n } ; matricea A T A este pozitiv definită, deci inversabilă (nesingulară). Concluzie: Problema de calcul are o soluţie unică dacă: b 0 m , b Im(A ) în ipoteza că rang(A) = n. Trebuie reformulată problema de calcul pentru a asigura existenţa unei soluţii în caz contrar. Principiul folosit este minimizarea unei funcţii criteriu de tipul: V ( x ) || b A x || * * x - pseudosoluţie a sistemului dacă V ( x ) minim METODE NUMERIC – curs 5 Condiţii de minim: T V ( x ) V ( x ) x {V ( x )} 0 nx1 x 1 x n V - diferenţiabilă not x x {V ( x )} x {V ( x )} x { Tx {V ( x )}} 0 1 1 1 T 1 m 2 V ( x ) || b A x || 2 || r || 2 r r ri 2 2 2 2 2 2 i 1 Problema de calcul devine: * || b A x || 22 min n 1 {|| b A x || 2 2} x* - pseudosoluţie în sensul celor mai mici pătrate x Condiţii de minim: sistem de ecuaţii normale x {V( x )} (1 / 2) 2 (A T A x A T b) 0 n x {V( x )} (1 / 2) 2 (A T A) 0 METODE NUMERIC – curs 5 3.2 Triangularizarea ortogonală a matricelor 3.2.1 Matrici ortogonale Definiţie: Fie o matrice Q mm. Matricea Q se numeşte ortogonală dacă este îndeplinită una din relaţiile: QT Q Im Q T Q 1 Q [q 1 q i q j q m ], q i m1 , i 1,..., m T q i q j 0, i j, i, j 1,..., m; || q i || 22 1 Proprietate o matrice ortogonală păstrază norma vectorială euclidiană: z m1 , || Q z || 22 || z || 22 METODE NUMERIC – curs 5 T matricele Householder: uu U Im 1 T 1 u m1 - vector Houseolder u u || u || 22 2 2 proprietate: U T U U 1 3.2.2 Procedura de triangularizare ortogonală a unei matrici de rang complet Propoziţie: m m Oricare ar fi matricea A mn , m n de rang complet pe coloane, există o matrice Q m n ortogonală, astfel încât matricea A se poate scrie: A Q R unde matricea R are următoarea structură: R1 R 0 ( m n )n unde matricea R 1 nn este superior triunghiulară nesingulară cu rii 0, i 1,..., n METODE NUMERIC – curs 5 Demonstraţia este constructivă şi constituie însuşi algoritmul de triangularizare ortogonală a matricei A (factorizare QR a matricii A): atribuie A 1 A pentru k 1, n execută * determinare matrice Householder U k astfel încât: ( U k A k ) i, k 0, i k 1,..., m atribuie A k 1 U k A k atribuie R A n 1 Fie [1 k 1 k k 1 m ]T m1 coloana k a matricii Ak m m m 2k i2 0 ⇒ U k astfel încât U k [1 k 1 k 0 0]T , 2k 2k i k T U k Im u k u k / k , k || u k || 22 / 2 u k [0 0 u k , k u k 1, k u m ,k ]T m k sign ( k ) i2 , u k ,k k k , u i ,k i , i k 1,..., m, ik k k u k , k , k k. METODE NUMERIC – curs 5 k 0 m 2k i2 0 ⇒ i 0, i k ,..., m ⇒ u k 0m Daca i k ⇒ U k nu există k 0 (| k | impus ) Fie [1 m ]T u u T u u T Uk Im k k k k k k T ⇒ U k u k uk m u j, k j k j k k ⇓ i , i 1,..., k 1 U k i i u i , k , i k ,..., m Daca [ j] [1 j 0 0] T m1 j 1,..., k 1 ⇒ U k [ j] [ j] METODE NUMERIC – curs 5 Tabloul general al transformărilor: U n U n 1 U 1 A R ⇒ UAR U U n U n 1 U 1 ⇒ A Q R Q U 1 U T U1 U 2 U n 3.3 Rezolvarea sistemelor supradeterminate mn m1 n1 U A x b, A , b , x UAx Ub A m n rang( A) n mn R1 R x d; R U A , d U b d1 , d n1 , d ( m n )1 x 1 2 0 ( m n )n d 2 r b A x, V ( x ) || r || 22 || U r || 22 || d 1 R 1 x || 22 || d 2 || 22 ! ⇒ || d 1 R 1 x || 22 0 d 2 0 ( mn ) ⇓ exemplu * R 1 x d1 METODE NUMERIC – curs 5 Cap. 4 Calculul valorilor şi vectorilor proprii 4.1 Formularea problemei nn - Se consideră o matrice reală, pătratică, de ordin n: A Definiţie: nn Oricare ar fi matricea A , un număr în general complex, λ , se numeşte valoare n 1 proprie a matricei A dacă există un vector, x C , x 0 n astfel încât: Ax x x – vector propriu al matricii A asociat valorii proprii λ ( I n A ) x 0 n x 0n I n A - singulară ecuaţie caracteristică p n ( ) det( I n A ) 0 polinom caracteristic p n () n 1 n 1 n 1 n i , i 1,..., n METODE NUMERIC – curs 5 Teoremă de existenţă: Orice matrice pătratică reală, de ordin n, are exact n valori proprii, în general complexe şi nu neapărat distincte, care coincid cu rădăcinile polinomului caracteristic ataşat matricei. Dacă există valori proprii complexe, atunci acestea apar în perechi complex conjugate. Orice matrice A nn are cel puţin un vector propriu. Nu se recomandă calculul numeric al valorilor proprii prin rezolvarea ecuaţiei caracteristice deoarece: - rezolvarea ecuaţiilor polinomiale este o problemă prost condiţionată; - coeficienţii polinomului caracteristic volum mare de calcule erori Metodele practice pentru calculul numeric al valorilor proprii proceduri iterative - matricea A adusă la formă canonică Schur prin transformări ortogonale de asemănare METODE NUMERIC – curs 5 4.2 Forma canonică Schur Definiţie: nn Două matrici, A, A' nn , se numesc ortogonal asemenea, dacă există o matrice ortogonală , Q astfel încât: A' Q T A Q Proprietăţi: Matricile ortogonal asemenea au aceleaşi valori proprii: i (A) 'i (A' ), i 1,..., n Relaţia dintre vectorii proprii ai două matrici ortogonal asemenea: x i (A) Q x 'i (A ' ), i 1,..., n Definiţie: O matrice se spune că este în formă bloc superior triunghiulară, dacă are următoarea structură: T11 T12 T1m 0 T22 T2 m T , Tii pi pi Tii , i 1,..., m - matrici pătratice 0 0 Tmm p i 1, 2, i 1, m formă cvasi-superior triunghiulară METODE NUMERIC – curs 5 Teoremă de existenţă: ~ Oricare ar fi matricea A nn , există o matrice ortogonală Q nn , astfel încât matricea: ~T ~ forma canonică Schur reală a matricei A S Q A Q este în formă cvasi-superior triunghiulară. Blocurile diagonale de ordin întâi ale matricei S reprezintă valorile proprii reale ale matricei A şi ale matricei S, iar blocurile diagonale de ordin doi au valori proprii complex conjugate reprezentând valori proprii complex conjugate ale matricelor A şi S. Observaţii: ~ ~ coloanele matricii Q nn ,q k , k 1,..., n vectori Schur ai matricii A Demonstraţia teoremei este constructivă algoritmul QR 4.3 Algoritmul QR pentru calculul formei canonice Schur Principiu: construcţia unui şir de matrici ortogonal asemenea convergent la forma canonică Schur A 0 A, A 1 ,..., A k ,... A k k S A k a ijk 1 i , j n a ijk k 0 i j 2 a ik1 ,i k 0 i=1,...,n-1 METODE NUMERIC – curs 5 Se defineşte şirul de matrici{A k } k 0 , A 0 A pas QR cu deplasare explicită A k k I n Q k R k , A k 1 R k Q k k I n Q k nn - ortogonală; R k nn - superior triunghiulară k - deplasare (cu rol de accelerare a convergenţei) Propoziţie: Matricele şirului QR sunt ortogonal asemenea: A k 1 Q Tk A k Q k Observaţie: - algoritmul original QR consumator de timp se foloseşte o formă optimizată Forma implementabilă a algoritmului QR cu deplasare explicită parcurge două faze de lucru: faza 1 – pregătitoare (zerorizare elemente de sub sub-diagonala principală) faza 2 – aplicare algoritm QR matricii obţinută la faza 1 se obţine forma canonică Schur METODE NUMERIC – curs 5 Faza 1 - procedură directă de aducere a matricei A la forma superior Hessenberg (H) A= H= 0 c1 c2 ck cn-2 - se folosesc matrici Householder în scopul anulării, coloană cu coloană, a elementelor matricei A situate sub sub-diagonala principală - algoritm: atribuie A 1 A pentru k 1, n 2 execută * determinare matrice Householder astfel încât: ( U k 1 A k ) i ,k 0, i k 2,..., n atribuie A'k 1 U k 1 A k atribuie A k 1 A'k 1 U Tk 1 A'k 1 U k 1 atribuie H A n 1 METODE NUMERIC – curs 5 sinteza matricii Householder, Uk+1 u k 1 [0 0 u k 1, k 1 u n , k 1 ]T n (a [i,kk] ) 2 , T U k 1 I n ( u k 1 u k 1 / k 1 ) k 1 sign (a [kk]1, k ) u k 1, k 1 a [kk]1,k k 1 i k 1 u i,k a [i,kk] , i k 2,..., n, k 1 k 1 u k 1,k 1 tabloul general al transformărilor: U n 1 U 2 A U 2 U n 1 H U UT H U A UT sunt parcurse exact ( n – 2 ) iteraţii METODE NUMERICE – curs 6 4.2 Forma canonică Schur Definiţie: nn Două matrici, A, A' nn , se numesc ortogonal asemenea, dacă există o matrice ortogonală , Q astfel încât: A' Q T A Q Proprietăţi: Matricile ortogonal asemenea au aceleaşi valori proprii: i (A) 'i (A' ), i 1,..., n Relaţia dintre vectorii proprii ai două matrici ortogonal asemenea: x i (A) Q x 'i (A ' ), i 1,..., n Definiţie: O matrice se spune că este în formă bloc superior triunghiulară, dacă are următoarea structură: T11 T12 T1m 0 T22 T2 m T , Tii pi pi Tii , i 1,..., m - matrici pătratice 0 0 Tmm p i 1, 2, i 1, m formă cvasi-superior triunghiulară METODE NUMERICE – curs 6 Exemplu: - matrice cvasi-superior triunghiulară - matrice în formă canonică Schur valori proprii valori proprii reale complex conjugate Teoremă de existenţă: ~ Oricare ar fi matricea A nn , există o matrice ortogonală Q nn , astfel încât matricea: ~T ~ forma canonică Schur reală a matricei A S Q A Q este în formă cvasi-superior triunghiulară. Blocurile diagonale de ordin întâi ale matricei S reprezintă valorile proprii reale ale matricei A şi ale matricei S, iar blocurile diagonale de ordin doi au valori proprii complex conjugate reprezentând valori proprii complex conjugate ale matricelor A şi S. Observaţii: ~ ~ coloanele matricii Q nn ,q k , k 1,..., n vectori Schur ai matricii A Demonstraţia teoremei este constructivă algoritmul QR METODE NUMERICE – curs 6 4.3 Algoritmul QR pentru calculul formei canonice Schur Principiu: construcţia unui şir de matrici ortogonal asemenea convergent la forma canonică Schur A 0 A, A 1 ,..., A k ,... A k k S A k a ijk 1 i , j n a ijk k 0 i j 2 a ik1 ,i k 0 i=1,...,n-1 se defineşte şirul de matrici {A k }k 0 , A 0 A pas QR cu deplasare explicită A k k I n Q k R k , A k 1 R k Q k k I n Q k nn - ortogonală; R k nn - superior triunghiulară k - deplasare (cu rol de accelerare a convergenţei) METODE NUMERICE – curs 6 Propoziţie: Matricele şirului QR sunt ortogonal asemenea: A k 1 Q Tk A k Q k Observaţie: - algoritmul original QR consumator de timp se foloseşte o formă optimizată Forma implementabilă a algoritmului QR cu deplasare explicită parcurge două faze de lucru: faza 1 – pregătitoare (zerorizare elemente de sub sub-diagonala principală) diagonala principală subdiagonala principală faza 2 – aplicare algoritm QR matricii obţinută la faza 1 se obţine forma canonică Schur METODE NUMERICE – curs 6 Faza 1 - procedură directă de aducere a matricei A la forma superior Hessenberg (H) A= H= 0 c1 c2 ck cn-2 - se folosesc matrici Householder în scopul anulării, coloană cu coloană, a elementelor matricei A situate sub sub-diagonala principală - algoritm: atribuie A 1 A pentru k 1, n 2 execută * determinare matrice Householder astfel încât: ( U k 1 A k ) i ,k 0, i k 2,..., n atribuie A'k 1 U k 1 A k atribuie A k 1 A'k 1 U Tk 1 A'k 1 U k 1 atribuie H A n 1 METODE NUMERICE – curs 6 sinteza matricii Householder, Uk+1 u k 1 [0 0 u k 1, k 1 u n , k 1 ]T n (a [i,kk] ) 2 , T U k 1 I n ( u k 1 u k 1 / k 1 ) k 1 sign (a [kk]1, k ) u k 1, k 1 a [kk]1,k k 1 i k 1 u i,k a [i,kk] , i k 2,..., n, k 1 k 1 u k 1,k 1 tabloul general al transformărilor: U n 1 U 2 A U 2 U n 1 H U UT H U A UT sunt parcurse exact ( n – 2 ) iteraţii METODE NUMERICE – curs 6 Faza 2 - procedură iterativă de construcţie a unui şir QR pornind de la matricea H - scopul: anularea unora din elementele de pe sub-diagonala principală a matricei H S11 0 Sii H= S= 0 0 0 Smm -algoritmul QR cu deplasare explicită parcurge următoarele etape: 1. determinare deplasare , ; deplasare explicită (n-1) iteraţii 2. atribuire H H I n ; pas QR cu 3. factorizare H Q R , R - matrice superior triunghiulară, Q - matrice ortogonală; 4. atribuire H R Q ; 5. refacere deplasare H H I n ; 6. corecţie matrice H (efectuare test de decuplare); 7. atribuire S H şi efectuare teste reluare algoritm QR (de la etapa 1). METODE NUMERICE – curs 6 Detaliere etape: Etapa 1 - determinare deplasare - viteza maximă de convergenţă pentru h n ,n Etapa 2 - deplasarea se scade de pe diagonala principală a matricei H - se lucrează în continuare cu această matrice Etapa 3 – factorizare QR a matricei H obţinută la etapa 2 - se zerorizează elementele de pe sub-diagonala principală (n – 1 elemente) - se utilizează matrici de rotaţie plană Givens (ortogonale, nesimetrice) planul (k, k+1) h k ,k rk k+1 R 0 (k)r h k 1,k (k+1)r h k ,k hk+1,k (rk, 0) c k cos (θ) ck dk rk R d k c k h k 1,k k d k sin (θ) rk hk,k METODE NUMERICE – curs 6 I k 1 0( k 1 ) 2 0( k 1 ) ( n k 1 ) - matricea de rotaţie Givens: Pk 0 2( k 1 ) R 0 2( n k 1 ) 0( n k 1 ) ( k 1 ) 0( n k 1 ) 2 I n k 1 nn - se poate demonstra: I k 1 0( k 1 ) 2 0( k 1 ) ( n k 1 ) Pk PkT 0 2( k 1 ) Q 0 2( n k 1 ) 0( n k 1 ) ( k 1 ) 0( n k 1 ) 2 I n k 1 Pk PkT I n n n c 2k d 2k 0 1 0 Q 2 2 0 c k d k 0 1 METODE NUMERICE – curs 6 - tabloul transformărilor de la etapa 3: Pn 1 P1 H R Q T Pn 1 P1 Q P1T PnT1 R – matrice superior triunghiulară Etapa 4 – matricea R se înmulţeşte la dreapta cu matricea Q, rezultatul fiind stocat în H: H R P1T PnT1 Etapa 5 – se adună deplasarea la elementele de pe diagonala prinicipală a matricii H Etapa 6 – test de decuplare Pk h k ,k h k ,k 1 h 'k ,k h 'k ,k 1 PkT h "k ,k h "k ,k 1 h k 1,k h k 1,k 1 0 h 'k 1,k 1 h "k 1,k h "k 1,k 1 k 1,..., n 1 | h k 1, k || h "k 1, k | - inegalitatea strictă elementul din poziţia (k+1, k) devine zero în forma canonică Schur METODE NUMERICE – curs 6 - rol de decuplare a blocurilor diagonale din forma canonică Schur test de decuplare: dacă | h "k 1,k | (| h "k ,k | | h "k 1,k 1 |) atunci atribuie h "k 1, k 0 Etapa 7 – testele care condiţionează reluarea sau nu a algoritmului QR Mi 0 sii 0 sii si,i+1 0 si+1,i s i+1,i+1 0 concluzii privind structura formei canonice Schur: - nu există pe sub-diagonala principală două elemente consecutive nenule - blocuri de ordinul doi pe diagonala matricei S au valori proprii complex conjugate METODE NUMERICE – curs 6 testele care se realizează la această etapă sunt: (T1) dacă i {1, 2,..., n 2} pentru care s i 1,i 0 şi s i 2,i 1 0 atunci matricea S nu este în forma canonică Schur (reluare de la etapa 1) (T2) dacă i care să satisfacă (T1), dar i {1, 2,..., n 1} pentru care s i ,i s i ,i 1 blocul M i are valori proprii reale atunci s i 1,i s i 1,i 1 matricea S nu este în forma canonică Schur (reluare de la etapa 1) Tabloul general al transformărilor din faza a doua a algoritmului QR: S [ Pn[s]1 P1[ s ] ] [ Pn[1]1 P1 ] H [(P1 ) T (Pn[1]1 ) T ] [(P1[ s ] ) T (Pn[ s]1 ) T ] P PT P H PT S METODE NUMERICE – curs 6 În urma aplicării ambelor faze ale formei implementabile a algoritmului QR rezultă: P U A UT PT S ~ ~ QT P U Q UT PT ~ ~ S QT A Q forma implementabilă a algoritmului QR este o procedură stabilă numeric, iar forma canonică Schur calculată pentru o matrice A, notată prin Ŝ, coincide cu forma canonică Schur exactă, S, a matricei A uşor perturbată, Â A E : Ŝ S A E , || E || || A || Exemplu METODE NUMERICE – curs 6 4.4 Calculul valorilor şi vectorilor proprii calculul valorilor proprii inspecţia blocurilor diagonale ale formei canonice Schur, S inspecţia formei canonice Schur elementele aflate pe prima sub-diagonală a matricei S Mi atribuie i 1 0 sii 0 sii si,i+1 cât timp (i n – 1) dacă si+1,i = 0 atunci 0 atribuie i si,i si+1,i si+1,i+1 atribuie i i + 1 0 dacă si+1,i 0 atunci atribuie i,i+1 valorile proprii ale blocului: s i ,i s i,i 1 det( I 2 M i ) 0 Mi s i 1,i s i 1,i 1 atribuie i i + 2 2 (s i ,i s i 1,i 1 ) (s i ,i s i 1,i 1 s i 1,i s i ,i 1 ) 0. METODE NUMERICE – curs 6 calculul vectorilor proprii ~ x i Q r i , i 1,..., n unde S r i i r i , i 1,..., n ~ în practică vectorii Schur ai matricei A vectorii coloană ai matricei Q ~ q i , i 1,..., n ~ - ortogonali Q [~ q1 ~ qi ~ qn] || ~ q i || 22 1, i 1,..., n ~ dacă i s ii (s i 1,i 0) atunci x i q i ~ ~ ~ ~ dacă i , i 1 C (si 1,i 0) atunci x i q i j q i 1 , x i 1 q i j q i 1 METODE NUMERICE – curs 7 Cap. 5 Descompunerea valorilor singulare 5.1 Formularea problemei A m n p min{m, n} r p - matrice deficientă de rang r - rangul matricii A r p r p - matrice de rang complet Teoremă: ~ ~ Oricare ar fi matricea A mn , există două matrici ortogonale U mm şi V nn , astfel încât: ~ ~ ~ ~ UT A V sau A U V T , (1) unde este o matrice pseudo-diagonală, elementele nenule ale acesteia satisfăcând relaţia: diag{1 , , p }, 1 2 p 0. (2) 1 , , p - valori singulare ale matricii A mn - formă canonică diagonală a matricii A relaţia (1) – descompunerea valorilor singulare METODE NUMERICE – curs 7 Demonstraţia teoremei algoritmul de descompunere a valorilor singulare Observaţie: - orice matrice este ortogonal echivalentă bilateral cu o matrice (pseudo)diagonală detaliind structura matricei mn , pot apare următoarele situaţii: mn rang(A) r p n rang(A) r p 1 ×( ) , 1 diag{1 , , n } 0 ( m n )n ( )× ( )×( ) × Exemple: - Matrice de rang complet: - Matrice deficienta de rang: METODE NUMERICE – curs 7 mn rang(A) r p m [1 0 m( n m ) ], 1 diag{1 , , m } rang(A) r p 1 0 r ( n r ) r r , 1 , 1 diag{1 , , r } 0 ( m r )r 0 ( m r )( n r ) Exemple: - Matrice de rang complet: - Matrice deficienta de rang: METODE NUMERICE – curs 7 1 2 r 0, r 1 p 0 rangul matricei A este egal cu rangul matricei , care este egal cu r rangul acestor matrici este egal cu numărul valorilor singulare nenule ~ ~ relaţia (1) poate fi scrisă sub forma: A V U ~ ~ vi - vectorii coloană ai V ~ ui ~ - vectorii coloană ai U A ~v i i ~ ui ~ u i A i ~v i T T i 1,..., p Definiţie: ~ Coloanele matricei V se numesc vectori singulari la dreapta ai matricei A, iar coloanele ~ matricei U se numesc vectori singulari la stânga ai matricei A. METODE NUMERICE – curs 7 , proprietăţi ale decompunerii valorilor singulare: P1. valorile singulare ale matricei A sunt egale cu rădăcina pătratică pozitivă a valorilor proprii ale matricei sau : - reprezintă punctul de plecare în dezvolatarea algoritmului de descompunere P2. dacă A nn este o matrice simetrică şi pozitiv semidefinită, atunci valorile proprii ale matricei A sunt reale, pozitive, iar valorile singulare ale matricei A sunt egale cu valorile ei proprii: i (A ) i ( A ), i 1,..., n METODE NUMERICE – curs 7 ~ ~ P3. scriind matricile U şi V sub forma: ~ ~ ~ U [~ u1 ~ ur ~ u r 1 ~ u m ] [ U1 U 2 ] ~ ~ ~ V [~v 1 ~v r ~ v r 1 ~v n ] [V1 V2 ] ~ 1 0 r( n r ) V1T ~ ~ ~ ~ A U V T [ U1 U2 ] ~ 0 ( m r )r 0 ( m r )( n r ) V2T METODE NUMERICE – curs 7 , ~ P4. primele r coloane ale matricei U , unde r rang(A) , formează o bază ortogonală pentru subspaţiul imagine al matricei A: Im(A) {y / y A x , x n1 }, r r y Ax ~ u i ( i ~ v i x) i ~ T ui i 1 i 1 ~ P5. ultimele n - r coloane ale matricei V , formează o bază ortogonală pentru subspaţiul nul al matricei A: N(A) {x / A x 0 m , x n1 } ~ ~ ~ ~ ~ ~ A V U A [V1 V2 ] [ U 1 U 2 ] ~ ~ ~ A V1 U 1 1 , A V2 0 m( n r ) METODE NUMERICE – curs 7 5.2 Algoritmul SVD Algoritmul SVD (SVD – abr. eng. “Singular Value Decomposition”) construirea unui şir de matrici ortogonal echivalente bilateral, convergent către forma canonică pseudo-diagonală se bazează pe faptul că: i (A) i (A T A ) Principiul de calcul a aplica “mascat” (direct matricei A) algoritmul QR se aplică matricei A, transformările corespunzătoare celor care s-ar aplica matricei B în cadrul algoritmului QR de aducere la forma canonică Schur A B A T A, B B T , B n n se consideră m > n (fără a restrânge generalitatea) Algoritmul comportă două faze de lucru: Faza 1: pregătitoare, de aducere a matricei A la forma superior bidiagonală Faza 2: procedură iterativă, de construcţie a unui şir de matrice ortogonal echivalente bilateral, convergent către forma canonică pseudo-diagonală a matricei A METODE NUMERICE – curs 7 Faza 1 L1 Lk 0 Ln-2 A= J= 0 c1 ck cn U k mm , k 1,..., n - matrici Householder de ordin m care acţionează pe coloane Vk 1 nn , k 1,..., n 2 - matrici Householder de ordin n care acţionează pe linii Tabloul general al transformărilor: U n U k U 1 A V2 Vk 1 Vn 1 J U U n U k U 1 , V V2 Vk 1 Vn 1 UAV J Observaţie: algoritm QR BA AT H - matrice tridiagonală H JT J METODE NUMERICE – curs 7 Faza 2 - zerorizarea elementelor de pe prima supradiagonală a formei superior bidiagonale J 0 0 J= = 0 0 0 - se regăsesc cele şapte etape de la algoritmul QR cu deplasare explicită - se utilizează matrici de rotaţie plană Givens Tabloul general al transformărilor: [R [n]1 R 1[ ] ] [R [n1] 1 R 1 ] J [P1 Pn[1]1 ] [P1[ ] Pn[ 1] ] R P R JP R UAVP ~ ~ ~ UT R U UT A V ~ V VP METODE NUMERICE – curs 7 5.3 Aplicaţii ale descompunerii valorilor singulare - în urma aplicării SVD: ̂ - algoritmul SVD – stabil numeric: Â A E ˆ A A E , || E || 2 || A || 2 ˆ i , i 1,..., p - valori singulare calculate i , i 1,..., p - valori singulare exacte | ˆ i i | ' 5.3.1 Calculul rangului unei matrici ˆ 1 ˆ 2 ˆ p 0 Definiţie: Se numeşte rang efectiv sau numeric, notat cu r̂ , al matricei A, acel index pentru care are loc relaţia de ordine: ˆ r̂ ˆ ' r̂ 1 ˆ r̂ ˆ 1 ' ˆ r̂ 1 ˆ 1 ˆ 1 ' max{m, n} m 𝜎 ≥⋯≥𝜎̂ ≥𝜏>𝜎̂ ≥⋯≥𝜎 ≥0 ' ˆ 1 max{m, n} m ˆ 1 METODE NUMERICE – curs 7 5.3.2 Rezolvarea sistemelor de ecuaţii algebrice liniare în sensul celor mai mici pătrate generalizate A x b, A mn , b m1 , x n1 (3) - dacă A – de rang complet există o pseudosoluţie unică în sensul celor mai mici pătrate (m n) m
