Metoda puterii inverse#

Metoda puterii inverse reprezintă o adaptare a algoritmului practic discutat anterior, de aceea vă încurajăm să îl lecturați înainte de a trece mai departe.

Vom începe prin câteva considerente matematice, apoi vom descoperi împreună metoda (așa cum vom vedea, MPI nu reprezintă o singură metodă, ci o familie de metode asemănătoare).

MPI (versiunea nedeplasată)#

Fie matricea pătratică \(A\in\mathbb{R}^{n\times n}\), \(n\in\mathbb{N}^*\), aleasă astfel încât aceasta să fie nesingulară. Dacă \(\lambda_i\in\lambda(A)\) este o valoare proprie a lui \(A\) și \(\bm{x}\in\mathbb{R}^n\) este vectorul propriu asociat valorii \(\lambda_i\), atunci corespondentul perechii \((\lambda_i, \bm{x})\) pentru matricea \(A^{-1}\) este \(\left(\dfrac{1}{\lambda_i}, \bm{x}\right)\).

În acest algoritm (MPI), singura modificare adusă metodei MPD este înlocuirea calculului astfel încât să se utilizeze vectorul propriu al lui \(A^{-1}\). Deoarece acesta coincide și pentru \(A\), se poate utiliza această valoare în schimb.

Evident, pentru ca algoritmul să funcționeze, se presupune că \(\exists A^{-1}\), adică \(A\) este o matrice nesingulară.

Sinteza algoritmului#

Pe scurt, în baza celor discutate anterior, vom sintetiza direct algoritmul:

1. \(\bm{y'}^{(0)}\) primește o aproximare de început;

2. Pentru fiecare pas \(k\geq 1\), se va obține o aproximare tot mai bună astfel:

   • Se calculează un vector auxiliar: \(\bm{t}=A^{-1}\bm{y'}^{(k-1)}\);

   • Vectorul acesta normalizat va fi chiar \(\bm{y'}^{(k)}\), adică: \(\bm{y'}^{(k)}=\frac{\bm{t}}{||\bm{t}||}\);

   • Se calculează noua aproximare a lui \(\lambda_1\) folosind coeficientul Rayleigh, adică:

     \[\begin{split}\lambda_1^{(k)}&=r\left(\bm{y'}^{(k)}\right)=\frac{\left(\bm{y'}^{(k)}\right)^TA\left(\bm{y'}^{(k)}\right)}{\left|\left|\bm{y'}^{(k)}\right|\right|^2}\\&=\left(\bm{y'}^{(k)}\right)^TA\left(\bm{y'}^{(k)}\right)\end{split}\]

Pe parcursul algoritmului, se va verifica dacă eroarea relativă dintre \(\lambda_1^{(k)}\) și \(\lambda_1^{(k-1)}\) este suficient de mică, adică mai mică decât o valoare standard prestabilită, o toleranță aleasă.

MPI (versiunea deplasată)#

Fie matricea pătratică \(A\in\mathbb{R}^{n\times n}\), \(n\in\mathbb{N}^*\). Dacă \(\mu\in\mathbb{R}\), \(\lambda_i\in\lambda(A)\) este o valoare proprie a lui \(A\) și \(\bm{x}\in\mathbb{R}^n\) este vectorul propriu asociat valorii \(\lambda_i\), atunci corespondentul perechii \((\lambda_i, \bm{x})\) pentru matricea \((A-\mu I_n)\) este \((\lambda_i-\mu, \bm{x})\).

Matricea \((A-\mu I_n)\) este cunoscută drept deplasarea matricei \(A\) cu \(\mu\) unități. Totodată, dacă \(\mu\) nu este o valoare proprie a matricei \(A\), atunci matricea \((A-\mu I_n)\) va fi nesingulară, indiferent dacă \(A\) este sau nu singulară.

Așadar, punând aceste idei cap la cap, se poate adapta MPI (versiunea nedeplasată) astfel încât să utilizeze \((A-\mu I_n)\) în loc de \(A^{-1}\).

Sinteza algoritmului#

Scrisă ca un algoritm ce poate fi implementată ca atare, metoda poate fi sintetizată astfel:

1. \(\bm{y'}^{(0)}\) primește o aproximare de început;

2. Pentru fiecare pas \(k\geq 1\), se va obține o aproximare tot mai bună astfel:

   • Se calculează un vector auxiliar: \(\bm{t}=(A-\mu I_n)\bm{y'}^{(k-1)}\);

   • Vectorul acesta normalizat va fi chiar \(\bm{y'}^{(k)}\), adică: \(\bm{y'}^{(k)}=\frac{\bm{t}}{||\bm{t}||}\);

   • Se calculează noua aproximare a lui \(\lambda_1\) folosind coeficientul Rayleigh, adică:

     \[\begin{split}\lambda_1^{(k)}&=r\left(\bm{y'}^{(k)}\right)=\frac{\left(\bm{y'}^{(k)}\right)^TA\left(\bm{y'}^{(k)}\right)}{\left|\left|\bm{y'}^{(k)}\right|\right|^2}\\&=\left(\bm{y'}^{(k)}\right)^TA\left(\bm{y'}^{(k)}\right)\end{split}\]

Pe parcursul algoritmului, se va verifica dacă eroarea relativă dintre \(\lambda_1^{(k)}\) și \(\lambda_1^{(k-1)}\) este suficient de mică, adică mai mică decât o valoare standard prestabilită, o toleranță aleasă.

MPI (cu coeficient Rayleigh)#

Algoritmul anterior poate fi îmbunătățit prin găsirea celui mai potrivit coeficient \(\mu\in\mathbb{R}\). Deoarece se încearcă găsirea valorii proprii dominante, \(\mu\) poate fi ales în funcție de coeficientul Rayleigh (discutat anterior).

Deoarece calculăm la fiecare pas \(k\geq 1\) un coeficient Rayleigh și îl reținem în \(\lambda^{(k)}\), putem să utilizăm în cadrul oricărui pas \(k\) valoarea \(\lambda^{(k-1)}\) pe post de \(\mu\).

Sinteza algoritmului#

Algoritmul MPI cu deplasare va suferi o foarte mică modificare:

1. \(\bm{y'}^{(0)}\) primește o aproximare de început;

2. Pentru fiecare pas \(k\geq 1\), se va obține o aproximare tot mai bună astfel:

   • Se calculează un vector auxiliar, \(\bm{t}=(A-\lambda^{(k-1)} I_n)\bm{y'}^{(k-1)}\);

   • Vectorul acesta normalizat va fi chiar \(\bm{y'}^{(k)}\), adică: \(\bm{y'}^{(k)}=\frac{\bm{t}}{||\bm{t}||}\);

   • Se calculează noua aproximare a lui \(\lambda_1\) folosind coeficientul Rayleigh, adică:

     \[\begin{split}\lambda_1^{(k)}&=r\left(\bm{y'}^{(k)}\right)=\frac{\left(\bm{y'}^{(k)}\right)^TA\left(\bm{y'}^{(k)}\right)}{\left|\left|\bm{y'}^{(k)}\right|\right|^2}\\&=\left(\bm{y'}^{(k)}\right)^TA\left(\bm{y'}^{(k)}\right)\end{split}\]

Pe parcursul algoritmului, se va verifica dacă eroarea relativă dintre \(\lambda_1^{(k)}\) și \(\lambda_1^{(k-1)}\) este suficient de mică, adică mai mică decât o valoare standard prestabilită, o toleranță aleasă.

Analiza convergenței#

Motivul pentru care preferăm să utilizăm MPI (față de MPD) ține de minimizarea pașilor executați până la convergență. O evaluare amănunțită a acesteia poate fi găsită în cărțile lui Jun Lu[26][27] sau în subcapitolul dedicat.

Cod ilustrativ#

Pentru a implementa metoda puterii indirecte sub forma funcției mpi, vom utiliza notațiile:

• A - matricea pentru care se caută valoarea proprie dominantă;

• y0 - aproximația inițială pentru vectorul propriu dominant, echivalentul lui \(y_0\);

• tol - toleranța acceptată, adică diferența minimă relativă dintre soluțiile a două iterații consecutive pentru ca algoritmul să continue;

• max_iter - numărul maxim de iterații pe care îl poate executa algoritmul înainte să renunțe.

Funcția returnează tripletul [lambda, y, steps], alcătuit din:

• lambda - valoarea proprie predominantă, echivalentul lui \(\lambda_1\);

• y - vectorul asociat valorii proprii predominante;

• steps - numărul de pași până se ajunge la convergență.

function [lambda, y, steps] = MPI(A, y0, tol, max_iter)
    [m, n] = size(A);
    y = y0;
    lambda = y' * A * y;
    for steps = [1 : max_iter]
        t = (A - lambda) * y;
        y = t / norm(t);
        lambda_old = lambda;
        lambda = y' * A * y;
        err = abs((lambda - lambda_old) / lambda);
        if err < tol
            return
        end
    end
end

Licență#

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