Priame metódy

Riešiť sústavu lineárnych rovníc $ {\bf Ax = b }$ s regulárnou maticou rádu $n$, kde $n$ je veľké číslo (rádovo tisíce aj viac), je v technickej praxi veľmi častá úloha. Táto úloha sa dnes už, samozrejme, realizuje na počítači. V posledných desaťročiach bolo vypracovaných niekoľko veľmi úspešných algoritmov na ich riešenie. Vzhľadom na dôležitosť akú majú lineárne systémy v praxi, tejto problematike sa v matematických aj programátorských kruhoch ešte stále venuje veľká pozornosť.
Metódu riešenia sústav $ {\bf Ax = b }$, ktorá vedie k riešeniu (až na zaokrúhľovacie chyby) po konečnom počte krokov, nazývame priama metóda. Základným princípom priamych metód je eliminácia neznámych, ktorú už poznáme. Priame metódy využívame hlavne pre plné matice. Sú to také matice, ktorých väčšina prvkov je nenulová.
Vzhľadom na to, že riešenia trojuholníkových sústav ako aj Gaussova eliminačná metóda boli vysvetlené už v predchádzajúcich častiach, tu si rozoberieme len ich algoritmické a numerické aspekty. Riadkové úpravy prevádzané GEM môžeme interpretovať ako postupné násobenie matice A transformačnými maticami ${\bf T}_1, {\bf T}_2,\dots, {\bf T}_{n-1}$, kde

\begin{displaymath}
{\bf T}_k =
{
\left(
\begin{array}{rrrrrrrrr}
1 & 0 & \...
..._{n,k} & 0 & \dots & \dots &0 & 1 \\
\end{array} \right)
};
\end{displaymath}

kde

\begin{displaymath}m_{ik}=- \frac{a_{ik}^{(k-1)}}{a_{kk}^{(k-1)}},
\ i=k+1, k+2, \dots,n\end{displaymath}

sú multiplikátory k-teho kroku.
Označujeme
$ { \bf A}\equiv { \bf A}^{(0)} = ( a_{ij}^{(0)}),\ { \bf T}_1 { \bf A}^{(0)} \e...
...\dots,
{ \bf T}_{n-1} { \bf A}^{(n-2)} \equiv { \bf A}^{(n-1)} \equiv { \bf U};$
${ \bf b}\equiv{ \bf b}^{(0)},\ { \bf b}^{(1)} \equiv { \bf T}_1 { \bf b}^{(0)},...
...\bf b}^{(k-1)},\ \dots
{ \bf b}^{(n-1)} \equiv { \bf T}_{n-1} { \bf b}^{(n-2)}.$ Po realizácii všetkých $n-1$ krokov eliminácie dostávame redukovanú maticu sústavy a stĺpec pravých strán (pozri časť Gaussova eliminačná metóda). Tú označíme ${ \bf U}$ a stĺpec pravých strán bude ${ \bf y}$. Máme:

\begin{displaymath}{ \bf U}= { \bf A}^{(n-1)} = { \bf T}{ \bf A}^{(0)}, \ \ { \bf y}= { \bf b}^{(n-1)} =
{ \bf T}{ \bf b}^{(0)}, \end{displaymath}

kde

\begin{displaymath}{ \bf T}= { \bf T}_{(n-1)} { \bf T}_{(n-2)} \dots { \bf T}_2 { \bf T}_1. \end{displaymath}

Matica ${ \bf T}$ je regulárna, a preto existuje jej inverzná matica ${ \bf T}^{-1} $, ktorú označíme ${ \bf L}$. Pre túto maticu platí:

\begin{displaymath}
{ \bf L}\equiv { \bf T}^{-1} =
{
\left(
\begin{array}{rr...
... -m_{n2} & \dots & -m_{n-1,n} & 1 \\
\end{array} \right)
}.
\end{displaymath}

To teda znamená, že matica ${\bf A}$ sa dá rozložiť na súčin dolnej ( L ) a hornej ( U ) trojuholníkovej matice. Platí:

\begin{displaymath}{ \bf A}= { \bf L}{ \bf U}.\end{displaymath}

Stanovenie matíc L a U nazývame trojuholníkovým rozkladom alebo LU rozkladom. K danej matici A môžeme určiť viacero trojuholníkových rozkladov podľa toho, ako volíme diagonálne prvky matice L.

Príklad 18. Urobme pre danú maticu LU rozklad:

\begin{displaymath}
{
\left(
\begin{array}{rr}
2&10 \\
7& 44 \\
\end{array} \right)
}
\end{displaymath}

Riešenie: Danú maticu vieme rozložiť dvojakým spôsobom:

\begin{displaymath}
{
\left(
\begin{array}{rr}
2&10 \\
7& 44 \\
\end{arr...
... \begin{array}{rr}
2&10 \\
0& 9 \\
\end{array} \right)
}
\end{displaymath}

alebo

\begin{displaymath}
{
\left(
\begin{array}{rr}
2&10 \\
7& 44 \\
\end{arr...
...
\begin{array}{rr}
1&5 \\
0& 3 \\
\end{array} \right)
}
\end{displaymath}

$\clubsuit$

Platí veta o jednoznačnosti rozkladu:

Veta 4.1   Trojuholníkový rozklad regulárnej matice A je určený jednoznačne diagonálnymi prvkami matice L a môže byť určený Gaussovou elimináciou.

V časti o GEM sme už spomínali triedy matíc, kedy môže byť GEM (a teda aj trojuholníkový rozklad) realizovaný. Existujú aj matice, pre ktoré sa GEM nedá realizovať. Aj tento fakt je dôvodom pre modifikácie GEM.
Ďalším dôvodom modifikácií GEM sú numerické aspekty GEM. Ich pôvod je hlavne v nepresnosti aritmetických operácií realizovaných na počítači (viď kapitola Reálne čísla, časť Chyby aritmetických operácií). Tieto nepresnosti spôsobia, že namiesto skutočných matíc rozkladu L a U vypočítame $\bar { \bf L}$ a $ \bar { \bf U}$, pričom $\bar { \bf L}\bar { \bf U}\neq { \bf A}$. Ak teraz označíme $\bar { \bf A}= \bar { \bf L}\bar { \bf U}$, zaujíma nás rozdiel ${ \bf A}- \bar { \bf A}$. Existujú matice chýb E a F, pre ktoré platí:

\begin{displaymath}\bar { \bf L}= { \bf L}+ {\bf E}, \ \ \ \bar { \bf U}= { \bf U}+{\bf F}.
\end{displaymath}

Potom

\begin{displaymath}{ \bf A}-\bar { \bf A}= { \bf L}{ \bf U}-\bar { \bf L}\bar { \bf U}= {\bf E} { \bf U}- { \bf L}{\bf F}
-{\bf EF}.\end{displaymath}

Odtiaľ plynie dôležitý záver: Ak pri realizácii GEM vychádzajú pre multilpikátory alebo prvky redukovanej matice sústavy veľké čísla, sú tiež prvky matice L veľké čísla a preto rozdiel ${ \bf A}- \bar { \bf A}$ (hlavne pre veľké $n$ ) môže byť rádovo mnohonásobne väčší ako je chyba E a F. Preto aj rozdiel presného a vypočítaného riešenia môže byť neprijateľne veľký. V takomto prípade hovoríme, že GEM je numericky nestabilná.
Tieto dva dôvody nás vedú k modifikovaným metódam GEM. Tou modifikáciou bude eliminácia s výberom hlavného prvku. Prvok, ktorého prostredníctvom určujeme multiplikátory k-teho kroku, budeme nazývať hlavným prvkom (pivotom) k-teho kroku. Pri transformačnej matici ${\bf T}_k$ sme ho označili $ a_{(k,k)}^{(k-1)}$. Rovnicu, z ktorej vyberáme hlavný prvok k-teho kroku, nazývame hlavnou rovnicou k-teho kroku. Aby sme minimalizovali zaokrúhľovacie chyby, je vhodné vyberať za hlavné prvky také prvky matice $ { \bf A}^{(k)}$, ktoré majú čo najväčšiu absolútnu hodnotu. Chyba sa v takomto prípade pri ďalších násobeniach u GEM ďalej nezväčšuje. Ak vyberáme v danej fáze hlavný prvok zo všetkých prvkov, ktoré prichádzajú do úvahy, hovoríme o algoritme GEM s úplným výberom hlavného prvku. Ak vyberáme hlavný prvok len z niektorých prvkov (napr. len z jedného riadku alebo stĺpca matice) hovoríme o algoritme GEM s čiastočným výberom hlavného prvku.

Príklad 19. Riešme sústavu

\begin{displaymath}
\begin{array}{rrrrr}
0,0001x_1 & + & 1,00x_2 & = & 1,00 \\
1,00x_1 & + & 1,00x_2 & = & 2,00 \\
\end{array}\end{displaymath}

za predpokladu, že všetky operácie a zaokrúhľovania robíme na tri platné miesta. Je to úloha, ktorá nám umožní ukázať citlivosť jednotlivých algoritmov na zaokrúhľovacie chyby. Najskôr aplikujeme algoritmus GEM, teda bez výberu hlavného prvku. Pre multiplikátor $m_{21}$ bude preto platiť: $m_{21} = -10^{4}$. Tento násobok prvej rovnice pripočítame k druhej. Dostávame (pracujeme na tri platné miesta):


\begin{displaymath}
\begin{array}{rrrrr}
0,0001x_1 & + & 1,00x_2 & = & 1,00 \\
& - & 10000x_2 & = & -10000 \\
\end{array}\end{displaymath}

Odtiaľ potom máme aproximáciu riešenia (označíme ju ${ \bf x}_c$).

\begin{displaymath}{ \bf x}_c =(0,000; 1,00)^T.\end{displaymath}

Ak za hlavný prvok zvolíme číslo v pozícii $(2,1)$ v matici (to jest hlavná rovnica prvého kroku bude druhá rovnica), potom $m_{21}= -\frac{0,0001}{1,00} = -0,0001$.
Dostávame redukovanú maticu v tvare:

\begin{displaymath}
\begin{array}{rrrrr}
\ & \ & 1,00x_2 & = & 1,00 \\
1,00x_1 & + & 1,00x_2 & = & 2,00. \\
\end{array}\end{displaymath}

Jej riešením je vektor ${ \bf x}_c = ( 1,00;1,00)^T$. Porovnanie s presnejšou aproximáciou riešenia pôvodnej sústavy (označíme ju ako ${ \bf x}_t$) ${ \bf x}_t = (1,00010;0,99990)^T$ ukazuje, ktorý z použitých algoritmov sa ukázal ako lepší.

Metóda LU rozkladu

Už vieme, za akých podmienok možno regulárnu maticu A rozložiť na súčin trojuholníkových matíc ${ \bf L}=(l_{ij})$ a $ { \bf U}=(u_{ij})$, to jest platí :

\begin{displaymath}{ \bf A}= { \bf L}{ \bf U}.\end{displaymath}

Metóda riešenia sústavy lineárnych rovníc LU- rozkladom spočíva v tom, že najskôr stanovíme matice L a U a potom riešime dve sústavy

\begin{displaymath}{ \bf L}{ \bf y}= { \bf b},\ \ { \bf U}{ \bf x}= { \bf y}.\end{displaymath}

V dolnej trojuholníkovej matici L volíme na diagonále jedničky, to jest $ l_{ii} =1, i=1,2,\dots,n$ a U je horná trojuholníková matica. Ako stanovíme prvky matíc L a U, si ukážeme na nasledujúcom príklade.

Príklad 20. Chceme vypočítať prvky matíc

\begin{displaymath}
{ \bf L}= {
\left(
\begin{array}{rrr}
1&0 & 0 \\
l_{21}...
...u_{22} & u_{23}\\
0 & 0 & u_{33} \\
\end{array} \right)
}
\end{displaymath}

tak, aby platila rovnosť

\begin{displaymath}
{
\left(
\begin{array}{rrr}
1&0 & 0 \\
l_{21}&1&0 \\
...
...6 \\
4& 13 & 19\\
6 & 27 & 50 \\
\end{array} \right)
}.
\end{displaymath}


Riešenie: Podľa definície súčinu matíc musí platiť:

\begin{displaymath}
\begin{array}{rrrr}
2 = & { \bf r}_1({ \bf L}) { \bf s}_1({ ...
...) = & l_{31}.u_{11},\hbox{ a preto } & l_{31}=3,\\
\end{array}\end{displaymath}

z toho dostávame

\begin{displaymath}
{
{ \bf s}_1({ \bf L}) = {
\left(
\begin{array}{r}
1 \\
...
...egin{array}{r}
2 \\
0 \\
0 \\
\end{array} \right)
},
}
\end{displaymath}

ďalej

\begin{displaymath}
\begin{array}{rrrr}
5 = &{ \bf r}_1({ \bf L}) { \bf s}_2({ \...
..._{12}+l_{32}.u_{22},
\hbox{ a preto } &l_{32}=4,\\
\end{array}\end{displaymath}


\begin{displaymath}
{
{ \bf s}_2({ \bf L}) = {
\left(
\begin{array}{r}
0 \\
...
...egin{array}{r}
5 \\
3 \\
0 \\
\end{array} \right)
},
}
\end{displaymath}


\begin{displaymath}
\begin{array}{rrrr}
6 = &{ \bf r}_1({ \bf L}) { \bf s}_3({ \...
..._{32}.u_{23}+u_{33}
,\hbox{ a preto} &u_{33}=4, \\
\end{array}\end{displaymath}


\begin{displaymath}
{
{ \bf s}_3({ \bf L}) = {
\left(
\begin{array}{r}
0 \\
...
...egin{array}{r}
6 \\
7 \\
4 \\
\end{array} \right)
},
}
\end{displaymath}

Takže máme

\begin{displaymath}
{
\left(
\begin{array}{rrr}
1&0 & 0 \\
2&1&0 \\
3 & ...
...6 \\
4& 13 & 19\\
6 & 27 & 50 \\
\end{array} \right)
}.
\end{displaymath}

Z tohoto príkladu vidíme, že matice L a U počítame po stĺpcoch. Vo všeobecnosti máme

\begin{displaymath}a_{ij} ={ \bf r}_i({ \bf L}) { \bf s}_j({ \bf U})\end{displaymath}

a odtiaľ pre $ i \leq j$ dostávame

\begin{displaymath}a_{ij}= l_{i1}u_{1j}+ l_{i2}u_{2j}+\dots + l_{i,i-1} u_{i-1,j}+ 1.u_{ij}\end{displaymath}

a pre $i>j$ máme

\begin{displaymath}a_{ij} =l_{i1}u_{1j} + l_{i2}u_{2j}+\dots + l_{i,j-1}u_{j-1,j}+l_{ij}u_{jj}.\end{displaymath}

Vo vzorcoch je $l_{ik}=0$, pokiaľ $i<k,$ a $ u_{kj} =0$ pre $k>j.$ Z prvého vzorca vypočítame $u_{ij}$ pre $ i=1,2,\dots,j$ a z druhého $ l_{ij}$ pre $ i=j+1,j+2,\dots,n.$
Výhoda metódy LU- rozkladu je značná hlavne v prípadoch, keď riešime viac sústav s rôznymi pravými stranami a rovnakou maticou. Spravíme najskôr rozklad matice A a v ďalšom už len počítame sústavy s trojuholníkovými maticami.