Modern numerikus módszerek
Bene Gyula
Eötvös Loránd Tudományegyetem, Elméleti Fizikai Tanszék
1117 Budapest, Pázmány Péter sétány 1/A
3. hét
$QR$ dekompozíció
Legyen ${\bf A}$ tetszőleges valós elemű négyzetes mátrix. Ekkor mindig létezik az
$${\bf A}={\bf Q}{\bf R}$$
felbontás, ahol $Q$ ortogonális mátrix, ${\bf R}$ pedig felső
háromszögmátrix. Műveletigény: $\propto N^3$ (de néhányszor több, mint
LU-dekompozíció esetén).
Ennek segítségével lineáris egyenletrendszer megoldása:
$$\begin{eqnarray}
{\bf Q}{\bf R}{\bf x}&=&{\bf b}\\
{\bf R}{\bf x}&=&{\bf Q}^T{\bf b}\\
\end{eqnarray}$$
Az utóbbi egyenlet visszahelyettesítéssel az utolsó egyenlettel kezdve nehézség
nélkül megoldható. Műveletigény: $\propto N^2$.
$QR$ dekompozíció elvégzése Housholder-transzformációk sorozatával:
$${\bf R}={\bf Q}_n\dots {\bf Q}_2{\bf Q}_1{\bf A}$$
Legyen
$${\bf Q}_1={\bf 1}-2\frac{{\bf u}\cdot {\bf u}^T}{{\bf u}^T\cdot {\bf
u}}\;,$$
ahol
$${\bf u}={\bf x}-{\bf e}_1 |{\bf x}|\;,$$
itt ${\bf x}$ az ${\bf A}$ mátrix első oszlopa, ${\bf e}_1$ pedig $N$ elemű
egységvekter, amelynek az első komponense $1$, a többi $0$. Nyilván ${\bf Q}_1^T={\bf
Q}_1$ és ${\bf Q}_1^2={\bf 1}$, ezért ${\bf Q}_1$ ortogonális mátrix, másrészt
$${\bf Q}_1{\bf x}={\bf e}_1|{\bf x}|\;.$$
Következésképpen ${\bf Q}_1{\bf A}$ alakja
$${\bf Q}_1{\bf A}=\left(\begin{array}{ccccc}*& *& *& *& *\\ 0& {\color{red}*}& *& *& *\\0& {\color{red}*}& *& *& *\\0& {\color{red}*}& *& *& *\\0& {\color{red}*}& *& *& * \end{array}\right)$$
Legyen
$${\bf Q}_2=\left(\begin{array}{c@{}c@{}}\begin{array}{c}1\end{array}& \begin{array}{cccc}0& 0& 0&
0\end{array}\\ \begin{array}{c}0\\ 0\\ 0\\ 0\end{array}& {\bf 1}-2\frac{{\bf u}\cdot {\bf u}^T}{{\bf u}^T\cdot {\bf
u}}\end{array}\right)$$
A jobb alsó $(N-1)\times (N-1)$-es blokkban értelemszerűen az egységmátrix
$(N-1)\times (N-1)$-es, az ${\bf u}$ vektor $ (N-1)$-es méretű. Az utóbbit
ismét
$${\bf u}={\bf x}-{\bf e}_1 |{\bf x}|$$
definiálja, az itt szereplő vektorok $ (N-1)$ eleműek, az ${\bf x}$ vektort a
${\bf Q}_1{\bf A}$ mátrix pirossal jelölt elemei alkotják. Ezzel
$${\bf Q}_2{\bf Q}_1{\bf A}=\left(\begin{array}{ccccc}*& *& *& *& *\\ 0& *& *&
*& *\\0& 0& {\color{red}*}& *& *\\0& 0& {\color{red}*}& *& *\\0& 0&
{\color{red}*}& *& *\end{array}\right)$$
s.í.t.
Sajátértékproblémák
$${\bf A}{\bf v}=\lambda {\bf v}$$
Ismeretlenek: $\lambda$ sajátérték, ${\bf v}$ sajátvektor.
A sajátvektorra nézve az egyenletrendszer homogén lineáris:
$$\left({\bf A}-\lambda {\bf 1}\right){\bf v}=0$$
A megoldás feltétele (karakterisztikus egyenlet):
$${\rm det}\left({\bf A}-\lambda {\bf 1}\right)=0$$
Egy $n\times n$-es mátrixnak pontosan $n$ darab (nem feltétlenül különböző) sajátértéke van.
A mátrix hasonlósági transzformációval diagonalizálható, azaz
$${\bf T}^{-1}{\bf A}{\bf T}={\bf D}\equiv {\rm diag}(d_1\;,\;d_2\;,\;\dots\;d_n)$$
akkor és csak akkor, ha létezik $n$ darab lineárisan független sajátvektora, amelyek ${\bf T}$ oszlopait alkotják:
$${\bf A}{\bf T}={\bf T}{\bf D}$$
A ${\bf D}$ diagonálmátrix $d_j$ elemei a sajátértékek.
Normális mátrixok:
$${\bf A}{\bf A}^\dagger={\bf A}^\dagger{\bf A}$$
A normális mátrixok unitér transzformációval diagonalizálhatók, ami egyben azt is jelenti, hogy sajátvektoraik egymásra ortogonálisak.
$${\bf U}^{\dagger}{\bf A}{\bf U}={\bf D}$$
Néhány normális mátrixfajta:
- Hermitikus (önadjungált) mátrixok
$${\bf A}={\bf A}^\dagger$$
- Valós szimmetrikus mátrixok
$${\bf A}={\bf A}^T$$
- Unitér mátrixok
$${\bf A}^{-1}={\bf A}^\dagger$$
- (Valós) ortogonális mátrixok
$${\bf A}^{-1}={\bf A}^T$$
A hermitikus (és speciálisan a valós szimmetrikus) mátrixok minden sajátértéke valós. Az unitér (és speciálisan az ortogonális) mátrixok minden sajátértéke $1$ abszolút értékű (általában komplex).
Néhány általános tétel:
- Minden négyzetes mátrix unitér transzformációval háromszögmátrixra
transzformálható (Schur-dekompozíció)):
$${\bf U}^{\dagger}{\bf A}{\bf U}={\bf L}$$
ahol ${\bf L}$ alsó vagy felső háromszögmátrix.
- Minden négyzetes mátrix hasonlósági transzformációval Jordan-féle normálakra hozható (a színessel jelölt Jordan-blokkok mérete lehet más is, ez a mátrixtól függ):
$${\bf T}^{-1}{\bf A}{\bf T}=\left(\begin{array}{cccccccccccc}
\lambda_1&0&0&0&0&0&0&0&0&0&\dots &0\\
0&\lambda_2&0&0&0&0&0&0&0&0&\dots &0\\
0&0&{\color{red}{\lambda_3}}&{\color{red}{1}}&{\color{red}{0}}&0&0&0&0&0&\dots &0\\
0&0&{\color{red}{0}}&{\color{red}{\lambda_3}}&{\color{red}{1}}&0&0&0&0&0&\dots &0\\
0&0&{\color{red}{0}}&{\color{red}{0}}&{\color{red}{\lambda_3}}&0&0&0&0&0&\dots &0\\
0&0&0&0&0&{\color{blue}{\lambda_3}}&{\color{blue}{1}}&0&0&0&\dots &0\\
0&0&0&0&0&{\color{blue}{0}}&{\color{blue}{\lambda_3}}&0&0&0&\dots &0\\
0&0&0&0&0&0&0&{\color{green}{\lambda_4}}&{\color{green}{1}}&0&\dots &0\\
0&0&0&0&0&0&0&{\color{green}{0}}&{\color{green}{\lambda_4}}&0&\dots &0\\
\vdots&&&&&&&&&\ddots &&\vdots\\
\vdots&&&&&&&&&&\ddots &\vdots\\
0&\dots&&&&&&&&&&\lambda_k\\
\end{array}\right)$$
Ha fellépnek Jordan-blokkok, akkor a mátrixnak van(nak) elfajult
sajátértéke(i), és van olyan elfajult sajátértéke, amihez a
multiplicitásánal kevesebb számú sajátvektor tartozik.
- $QR$-dekompozíció ill. $QL$-dekompozíció: minden négyzetes mátrix felbontható
$${\bf A}={\bf Q}{\bf R}$$
ill.
$${\bf A}={\bf Q}{\bf L}$$
alakban, ahol ${\bf Q}$ unitér
mátrix (a kétféle felbontásnál persze különböző), ${\bf R}$ felső
háromszögmátrix, ${\bf L}$ alsó háromszögmátrix.
- $QR$-transzformáció ill. $QL$-transzformáció:
$${\rm Ha}\quad{\bf A}={\bf Q}{\bf R}\;,\quad{\rm akkor}\quad {\bf A'}={\bf R}{\bf Q}={\bf Q}^\dagger{\bf A}{\bf Q}$$
ill.
$${\rm Ha}\quad{\bf A}={\bf Q}{\bf L}\;,\quad{\rm akkor}\quad {\bf A'}={\bf L}{\bf Q}={\bf Q}^\dagger{\bf A}{\bf Q}$$
A $QR$-transzformáció ill. $QL$-transzformáció megőrzi a mátrix
- hermitikusságát (ha valós, akkor a szimmetrikusságát),
- Hessenberg-alakját (ha valós szimmetrikus, akkor a tridiagonális alakját),
ui.
$${\bf A'}={\bf R}{\bf Q}={\bf R}{\bf A}{\bf R}^{-1}\;,$$
belátható, hogy ${\bf R}^{-1}$ szintén felső háromszögmátrix. Ezenkívül
közvetlen számolással látszik, hogy
$$
\underbrace{\left(\begin{array}{cccccc}
\times& \times& \times& \times& \times& \times \\
0& \times& \times& \times& \times& \times\\
0& 0& \times& \times& \times& \times\\
0& 0& 0& \times& \times& \times\\
0& 0& 0& 0& \times& \times\\
0& 0& 0& 0& 0& \times\\
\end{array}\right)}_{\rm felső\; háromszögmátrix }\cdot
\underbrace{\left(\begin{array}{cccccc}
\times& \times& \times& \times& \times& \times \\
\times& \times& \times& \times& \times& \times\\
0& \times& \times& \times& \times& \times\\
0& 0& \times& \times& \times& \times\\
0& 0& 0& \times& \times& \times\\
0& 0& 0& 0& \times& \times\\
\end{array}\right)}_{\rm felső\; Hessenberg-mátrix }
=\underbrace{\left(\begin{array}{cccccc}
\times& \times& \times& \times& \times& \times \\
\times& \times& \times& \times& \times& \times\\
0& \times& \times& \times& \times& \times\\
0& 0& \times& \times& \times& \times\\
0& 0& 0& \times& \times& \times\\
0& 0& 0& 0& \times& \times\\
\end{array}\right)}_{\rm felső\; Hessenberg-mátrix }
$$
(itt az $\times$-ek a nullától különböző elemeket jelölik) és
$$
\underbrace{\left(\begin{array}{cccccc}
\times& \times& \times& \times& \times& \times \\
\times& \times& \times& \times& \times& \times\\
0& \times& \times& \times& \times& \times\\
0& 0& \times& \times& \times& \times\\
0& 0& 0& \times& \times& \times\\
0& 0& 0& 0& \times& \times\\
\end{array}\right)}_{\rm felső\; Hessenberg-mátrix }\cdot
\underbrace{\left(\begin{array}{cccccc}
\times& \times& \times& \times& \times& \times \\
0& \times& \times& \times& \times& \times\\
0& 0& \times& \times& \times& \times\\
0& 0& 0& \times& \times& \times\\
0& 0& 0& 0& \times& \times\\
0& 0& 0& 0& 0& \times\\
\end{array}\right)}_{\rm felső\; háromszögmátrix }
=\underbrace{\left(\begin{array}{cccccc}
\times& \times& \times& \times& \times& \times \\
\times& \times& \times& \times& \times& \times\\
0& \times& \times& \times& \times& \times\\
0& 0& \times& \times& \times& \times\\
0& 0& 0& \times& \times& \times\\
0& 0& 0& 0& \times& \times\\
\end{array}\right)}_{\rm felső\; Hessenberg-mátrix }
$$
Numerikus eljárások a sajátérték-probléma megoldására
- Jacobi-módszer (valós szimmetrikus mátrix esetén)
- Iteratív módszer
- Síkbeli forgatásokkal kinullázzuk az egyes nemdiagonális mátrixelemeket
- A különböző forgatások egymás hatását elronthatják, de a nemdiagonális
mátrixelemek négyzetösszegét minden alkalommal pontosan a nullává tett
mátrixelemek négyzetösszegével csökken.
- Givens-módszer (valós szimmetrikus mátrix transzformációja tridiagonális
alakra)
- Nem iteratív módszer
- Síkbeli forgatásokkal történik az egyes mátrixelemek kinullázása, de a
$p,q$ síkbeli forgatással valamelyik $a_{pr}$-t ($r\ne q$) tesszük nullává.
- Véges számú ($\propto N^2$) forgatással a mátrix tridiagonális alakra
hozható anélkül, hogy az egyes forgatások egymás hatását elrontanák.
- A Householder-módszernél nem előnyösebb.
- Householder-módszer (valós szimmetrikus mátrix transzformációja
tridiagonális alakra)
- Householder-transzformáció:
$${\bf P}={\bf 1}-2\frac{{\bf u}\cdot{\bf u}^T }{{\bf u}^T\cdot{\bf u} }$$
Itt
$$ {\bf u}={\bf x}\pm\left|{\bf x}\right|{\bf e}_1$$
ahol ${\bf x}$ a transzformálandó mátrix első oszlopa, $ {\bf e}_1$ pedig
egységvektor, melynek első eleme $1$.
- Ekkor ${\bf P}$ szimmetrikus, ${\bf P}^2=1$, tehát ${\bf P}$ ortogonális
mátrix, továbbá
$${\bf P}{\bf x}=\pm\left|{\bf x}\right|{\bf e}_1$$
- Egy Householder-transzformáció tehát adott elem alatt egy egész oszlopot
nullává tesz.
- Az eredeti mátrix egyre kisebb méretű blokkjaira alkalmazzuk
($(N-1)\times(N-1)$-es mérettel kezdve), a
transzformációt komplementer méretű egységmátrixszal kiegészítve. Így a már
elért tridiagonális szerkezetet a újabb transzformációk megőrzik.
- $n-2$ lépésben előáll a tridiagonális alak.
- Általános mátrix transzformációja Hessenberg-alakra
Ugyanaz az eljárás, mint szimmetrikus mátrixoknál: $n-2$
Householder-transzformációval kinullázzuk az alsó első mellékdiagonális alatti mátrixelemeket.
- Valós Hessenberg-mátrix ill. valós szimmetrikus tridiagonális mátrix
sajátérték-problémájának megoldása $QR$-transzformációval:
Egyszerűen iterálni kell a $QR$ (vagy a $QL$) algoritmust:
$$\begin{eqnarray}
{\bf A_1}& := &{\bf A}={\bf Q_1}{\bf R_1}\\
{\bf A_2}& := &{\bf R_1}{\bf Q_1}={\bf Q_2}{\bf R_2}\\
\vdots\\
{\bf A_{n+1}}& := &{\bf R_n}{\bf Q_n}={\bf Q_{n+1}}{\bf R_{n+1}}
\end{eqnarray}
$$
végül nemelfajult sajátértékek esetében
$$\lim_{n\rightarrow \infty} {\bf A_{n}}={\rm felső \; háromszögmátrix}\;.$$
Programcsomagok:
Numerical Recipes,
LAPACK,
GSL
bene@arpad.elte.hu