Metoda deflației Householder#

Rugăm cititorul să treacă mai întâi atât prin noțiunile teoretice, cât și prin factorizarea QR folosind metoda Householder.

Prin deflație înțelegem transformarea unei matrice \(A\in\mathbb{R}^{n\times n}\) de dimensiune \(n\times n\) într-o altă matrice \(A'\in\mathbb{R}^{(n-1)\times (n-1)}\) de dimensiune \((n-1)\times(n-1)\), unde \(n\in\mathbb{N}^*\).

Această tehnică este utilă în aflarea valorilor proprii unei matrice, în mod special dacă nu se caută valoarea dominantă, ci o altă valoare proprie.

Descoperirea algoritmului#

Fie \(A\in\mathbb{R}^{n\times n}\), \(n\in\mathbb{N}^*\), alături de valoarea sa proprie dominantă, \(\lambda_1\in\lambda(A)\) și de vectorul propriu asociat acesteia, \(\bm{x}\in\mathbb{R}^n\). Pentru început, vom crea un reflector Householder pornind de la un vector \(\bm{u}\in\mathbb{R}^n\):

\[ H=I_n-2\cdot\frac{\bm{u}\bm{u}^T}{||\bm{u}||^2} \]

Suplimentar, dacă \(||\bm{u}||=1\), formula se simplifică:

\[ H=I_n-2\bm{u}\bm{u}^T \]

Vrem să forțăm asupra acestui reflector Householder următoarea proprietate:

\[ H\bm{x}=\begin{bmatrix} \alpha\\0\\\vdots\\0 \end{bmatrix}=\alpha\bm{e_1} \]

unde \(\alpha\in\mathbb{R}\) este o constantă oarecare.

Pentru a realiza acest lucru, ne amintim că forma pe care trebuie să o ia vectorul \(\bm{u}\) este dictată de formula \(\bm{u}=\frac{\bm{x}-\bm{y}}{||\bm{x}-\bm{y}||}\), unde \(H\bm{x}=\bm{y}\) și se impune \(||\bm{x}||=||\bm{y}||\).

În situația noastră, vrem \(H\bm{x}=\alpha \bm{e_1}\), dar trebuie ca \(||\bm{x}||=||\alpha\bm{e_1}||\), așadar \(\alpha=||\bm{x}||\). Ajungem la următoarea relație:

\[ \bm{u}=\frac{\bm{x}-||\bm{x}||\bm{e_1}}{\Big|\Big|\bm{x}-||\bm{x}||\bm{e_1}\Big|\Big|} \]

Folosind această relație, este ușor de demonstrat (prin calcul brut) că \(H\bm{x}=\alpha\bm{e_1}\). Vom folosi acest reflector pentru a face următoarea transformare de asemănare:

\[ G=HAH^{-1}\Rightarrow G=HAH^T\Rightarrow G=HAH \]
Inversa lui \(H\)

Amintim că matricea \(H\) este ortogonală, deci \(H^{-1}=H^T\), respectiv că \(H^2=I_n\), așadar \(H=H^{-1}\).

Deoarece \(\bm{x}\) este un vector propriu al matricei \(A\), știm că acesta respectă relația \(A\bm{x}=\lambda_1\bm{x}\). Așadar:

\[\begin{split} A\bm{x}=\lambda_1\bm{x}&\Rightarrow HA\bm{x}=\lambda_1 H\bm{x}\\ &\Rightarrow HAI_n\bm{x}=\lambda_1 H\bm{x}\\ \text{dar } H^2=H^TH=I_n&\Rightarrow (HAH^T)(H\bm{x})=\lambda_1 H\bm{x}\\ \text{dar } G=HAH^T&\Rightarrow \boxed{G(H\bm{x})=\lambda_1 H\bm{x}} \end{split}\]

Cunoaștem însă o scriere alternativă pentru \(H\bm{x}\), anume \(H\bm{x}=\alpha\bm{e_1}\). Făcând înlocuirea, se obține ceva cu adevărat impresionant:

\[ \alpha G\bm{e_1}=\lambda_1\alpha\bm{e_1}\Rightarrow \boxed{G\bm{e_1}=\lambda_1\bm{e_1}=\begin{bmatrix} \lambda_1\\0\\\vdots\\0 \end{bmatrix}} \]

Motivul pentru care relația anterioară este atât de interesantă este simplu - aceasta implică faptul că matricea \(G\) va avea în colțul din stânga sus valoarea \(\lambda_1\), în timp ce toate celelalte valori de pe prima linie devin irelevante (când se calculează \(G\bm{e_1}\), acestea se înmulțesc cu \(0\)), respectiv toate celelalte elemente de pe prima coloană trebuie să fie \(0\) (operația \(G\bm{e_1}\) înseamnă că se păstrează exclusiv conținutul primei coloane). Cu alte cuvinte:

\[ G=\begin{bmatrix} \lambda_1 & \bm{\boxtimes}\in\mathbb{R}^{n-1}\\ \bm{0}\in\mathbb{R}^{n-1} & A'\in\mathbb{R}^{(n-1)\times(n-1)} \end{bmatrix} \]

Am marcat cu \(\boxtimes\in\mathbb{R}^{n-1}\) un vector linie ale cărui valori nu sunt utile.

După cum observăm, matricea \(G\) conține o submatrice \(A'\) pe care se va putea repeta calculul recursiv pentru a se extrage următoarele valori proprii.

Sinteza algoritmului#

În baza celor discutate anterior, vom sintetiza algoritmul - pentru fiecare valoare \(i\in\overline{1,n}\), vom afla a \(i\)-a valoare proprie (descrescător din punct de vedere al magnitudinii) astfel:

  • Se calculează următoarele aproximări: valoarea proprie dominantă și vectorul propriu asociat acesteia pentru matricea curentă utilizând orice tehnică obișnuită, de exemplu MPD, MPI etc.;

  • Se calculează vectorul \(\bm{u}\) în baza celor discutate anterior: \(\bm{u}=\bm{u}=\frac{\bm{x}-||\bm{x}||\bm{e_1}}{\big|\big|\bm{x}-||\bm{x}||\bm{e_1}\big|\big|}\);

  • Se calculează reflectorul Householder prin formula \(H=I_n-2\bm{u}\bm{u}^T\);

  • Se calculează matricea \(\bm{G}=HAH\);

  • Se extrage valoarea proprie curentă, \(\lambda_1=g_{11}\);

  • Se reia algoritmul pe submatricea \(A'\):

    \[A'=\begin{bmatrix} g_{22}&g_{23}&\dots&g_{2n}\\ g_{32}&g_{33}&\dots&g_{3n}\\ \vdots&\vdots&\ddots&\vdots\\ g_{n2}&g_{n3}&\dots&g_{nn}\\ \end{bmatrix} \]

Utilizarea în practică#

Această metodă este utilă exclusiv în situația în care sunt necesari doar primele câteva valori proprii (\(2\), poate maxim \(3\)). Pentru mai multe valori proprii, se poate utiliza metoda QR.

Licență#

The book "Metode Numerice", written by Valentin-Ioan Vintilă, is licensed under CC BY-NC-SA 4.0