next up previous contents index
Next: Aproksimacija funkcije i numerička Up: Numerička matematika Previous: Rješavanje jednadžbi   Contents   Index

Subsections


Rješavanje sustava jednadžbi


Iterativne metode

Želimo riješiti sustav linearnih algebarskih jednadžbi

$\displaystyle a_{1\,1}\,x_1+a_{1\,2}\,x_2+\cdots+a_{1\,n}\,x_n$ $\displaystyle =$ $\displaystyle b_1$  
$\displaystyle a_{2\,1}\,x_1+a_{2\,2}\,x_2+\cdots +a_{2\,n}\,x_n$ $\displaystyle =$ $\displaystyle b_2$  
$\displaystyle \cdots$      
$\displaystyle a_{n\,1}\,x_1+a_{n\,2}\,x_2+\cdots +a_{n\,n}\,x_n$ $\displaystyle =$ $\displaystyle b_n,$  

koji matrično zapisan glasi

$\displaystyle A\,\boldsymbol{x} = \boldsymbol{b},$

gdje je

\begin{displaymath}
% latex2html id marker 38749
A = \left[
\begin{array}{cccc}
...
...n{array}{c}
x_1 \\  x_2 \\  \vdots \\  x_n
\end{array}\right].\end{displaymath}

To je matrična jednadžba, i mi ćemo često o sustavu jednadžbi govoriti kao o jednadžbi, misleći na ovu matričnu jednadžbu. Pretpostavimo da je $ A$ regularna matrica. Tada jednadžba ima rješenje, i označimo to rješenje sa $ \boldsymbol{s}.$

Osnovna ideja iterativnih metoda se sastoji u sljedećem. Stavimo

$\displaystyle A = B - C,$

gdje su $ B$ i $ C$ također kvadratne matrice $ n$-tog reda. Jednadžbu tada možemo prepisati kao

$\displaystyle B\,\boldsymbol{x} = C\,\boldsymbol{x} + \boldsymbol{b}.$

Ako s $ \boldsymbol{x}^{(k)}$ označimo $ k$-tu aproksimaciju rješenja, onda pomoću formule

$\displaystyle B\,\boldsymbol{x}^{(k+1)} = C\,\boldsymbol{x}^{(k)} + \boldsymbol{b},$

možemo naći $ k+1$-vu aproksimaciju rješenja.

Naravno, da bi postupak uopće krenuo, treba biti zadana početna aproksimacija $ \boldsymbol{x}_0.$ Nadalje, matrica $ B$ mora biti regularna i relativno jednostavna, da bismo imali rješenje i da bismo ga mogli relativno jednostavno izračunati. Također postupak nas mora približavati k rješenju, tj. postupak mora biti takav da

$\displaystyle \lim_{k\rightarrow{}\infty}\boldsymbol{x}^{(k)} = \boldsymbol{s}.$

Opišimo sada neke od iterativnih metoda.


Jacobijeva metoda

Pretpostavimo da su elementi na glavnoj dijagonali matrice $ A$ različiti od nule (ako je potrebno, premještanjem redaka u regularnoj matrici se to uvijek može postići). Rastavimo $ A$ kako slijedi

$\displaystyle A = L + D + R,$

gdje je

\begin{displaymath}
% latex2html id marker 38786
L = \left[
\begin{array}{ccccc}...
...{n\,1} & a_{n\,2} & \cdots & a_{n\,n-1} & 0
\end{array}\right],\end{displaymath}

\begin{displaymath}
% latex2html id marker 38788
D = \left[
\begin{array}{cccc}
...
...0 & a_{n-1\,n} \\
0 & 0 & \cdots & 0 & 0
\end{array}\right].\end{displaymath}

Tako je $ L$ donja, $ R$ gornja trokutasta, a $ D$ je dijagonalna matrica. Ako stavimo

$\displaystyle B = D,\hspace{1cm}C = -(L+R),$

imamo iterativni postupak oblika

$\displaystyle D\,\boldsymbol{x}^{(k+1)} = -(L+R)\,\boldsymbol{x}^{(k)} +
\boldsymbol{b}.$

Budući da su elementi na glavnoj dijagonali matrice $ A,$ a prema tome i matrice $ D$ različiti od nule, postoji $ D^{-1}.$ Inverz od dijagonalne matrice se vrlo lako računa. Iz $ D=[\delta_{ij}\,a_{ii}]$ slijedi $ D^{-1}=[\frac{\delta_{ij}}{a_{ii}}].$ Tako gornju jednadžbu možemo pomnožiti s lijeva s $ D^{-1},$ pa imamo sljedeći algoritam.

Algoritam 7   (Jacobijeva metoda) Proizvoljno izaberemo početnu aproksimaciju

% latex2html id marker 38813
$\displaystyle \boldsymbol{x}^{(0)}=\left[
\begin{array}{c}
x_1^{(0)} \\  x_2^{(0)} \\  \vdots \\  x_n^{(0)}
\end{array}\right],$

i zatim računamo sljedeće aproksimacije $ \boldsymbol{x}^{(k+1)}$ po formuli

$\displaystyle \boldsymbol{x}^{(k+1)} = - D^{-1}(L+R)\,\boldsymbol{x}^{(k)} +
D^{-1}\,\boldsymbol{b},\hspace{1cm}k=0,1,2,\ldots\ ,$

odnosno
$\displaystyle x_{1}^{(k+1)}$ $\displaystyle =$ $\displaystyle \alpha_{1\,2}\,x_{2}^{(k)} +
\alpha_{1\,3}\,x_{3}^{(k)} + \cdots + \alpha_{1\,n}\,x_{n}^{(k)} +
\beta_1$  
$\displaystyle x_{2}^{(k+1)}$ $\displaystyle =$ $\displaystyle \alpha_{2\,1}\,x_{1}^{(k)} +
\alpha_{2\,3}\,x_{3}^{(k)} + \cdots + \alpha_{2\,n}\,x_{n}^{(k)} +
\beta_2$  
$\displaystyle \cdots$      
$\displaystyle x_{n}^{(k+1)}$ $\displaystyle =$ $\displaystyle \alpha_{n\,1}\,x_{1}^{(k)} +
\alpha_{n\,2}\,x_{2}^{(k)} + \cdots +
\alpha_{n\,n-1}\,x_{{n-1}}^{(k)} + \beta_n,$  

gdje je $ \alpha{}_{ij}=-\frac{a_{ij}}{a_{ii}}, \beta{}_i=\frac{b_i}{a_{ii}}.$

Budući da dijelimo s $ a_{11},a_{22},\ldots,a_{nn},$ poželjno je, ako je to moguće, poredati retke u matrici $ A$ (jednadžbe u sustavu) tako da svaki element na glavnoj dijagonali bude po apsolutnoj vrijednosti veći od sume apsolutnih vrijednosti ostalih elemenata u retku u kojem se nalazi.


Gauss-Seidelova metoda

Gauss-Seidelova metoda je poboljšanje Jacobijeve metode u sljedećem smislu.

Pod pretpostavkom da je

% latex2html id marker 38847
$\displaystyle \boldsymbol{x}^{(k)}=\left[
\begin{array}{c}
x_1^{(k)} \\  x_2^{(k)} \\  \vdots \\  x_n^{(k)}
\end{array}\right]$

$ k$-ta aproksimacija rješenja, $ k+1$-vu nalazimo iz sustava
$\displaystyle a_{1\,1}\,x_{1}^{(k+1)} + a_{1\,2}\,x_{2}^{(k)} + \cdots +
a_{1\,n}\,x_{n}^{(k)}$ $\displaystyle =$ $\displaystyle b_1$  
$\displaystyle a_{2\,1}\,x_{1}^{(k+1)} + a_{2\,2}\,x_{2}^{(k+1)} + \cdots +
a_{2\,n}\,x_{n}^{(k)}$ $\displaystyle =$ $\displaystyle b_2$  
$\displaystyle \cdots$      
$\displaystyle a_{n\,1}\,x_{1}^{(k+1)} + a_{n\,2}\,x_{2}^{(k+1)} + \cdots +
a_{n\,n}\,x_{n}^{(k+1)}$ $\displaystyle =$ $\displaystyle b_n$  

Dakle $ x_1^{(k+1)},$ koji smo izračunali iz prve jednadžbe pomoću komponenti $ k$-te aproksimacije rješenja, koristimo odmah u drugoj, trećoj, ..., $ n$-toj jednadžbi; $ x_2^{(k+1)}$ koristimo u trećoj, četvrtoj, ..., $ n$-toj jednadžbi; itd. Taj sustav možemo zapisati u matričnom obliku

$\displaystyle (L + D)\,\boldsymbol{x}^{(k+1)} + R\,\boldsymbol{x}^{(k)} = \boldsymbol{b},$

odnosno

$\displaystyle (L + D)\,\boldsymbol{x}^{(k+1)} = -R\,\boldsymbol{x}^{(k)} + \boldsymbol{b}.$

Tako imamo formulu iterativnog postupka

$\displaystyle \boldsymbol{x}^{(k+1)} = -(L + D)^{-1}\,R\,\boldsymbol{x}^{(k)} + (L +
D)^{-1}\,\boldsymbol{b}.$

Nedostatak ove formule je u tome što treba naći inverz matrice $ L+D,$ što je teže nego naći inverz od $ D.$ Zato radimo malo drukčije. Pomnožimo ovu jednadžbu s $ D^{-1}.$

$\displaystyle (D^{-1}\,L + I)\,\boldsymbol{x}^{(k+1)} = -D^{-1}\,R\,\boldsymbol{x}^{(k)} +
D^{-1}\,\boldsymbol{b}.$

Odatle

$\displaystyle \boldsymbol{x}^{(k+1)} = -D^{-1}\,L\,\boldsymbol{x}^{(k+1)} -D^{-1}\,R\,\boldsymbol{x}^{(k)} +
D^{-1}\,\boldsymbol{b}.$

Tako dolazimo do sljedećeg algoritma.

Algoritam 8   (Gauss-Seidelova metoda) Izaberemo proizvoljnu početnu aproksimaciju

% latex2html id marker 38901
$\displaystyle \boldsymbol{x}^{(0)}=\left[
\begin{array}{c}
x_1^{(0)} \\  x_2^{(0)} \\  \vdots \\  x_n^{(0)}
\end{array}\right],$

i zatim računamo sljedeće aproksimacije $ \boldsymbol{x}^{(n)}$ po formuli

$\displaystyle \boldsymbol{x}^{(k+1)} = -D^{-1}\,L\,\boldsymbol{x}^{(k+1)} -D^{-1}\,R\,\boldsymbol{x}^{(k)} +
D^{-1}\,\boldsymbol{b},$

odnosno
$\displaystyle x_{1}^{(k+1)}$ $\displaystyle =$ $\displaystyle \alpha_{1\,2}\,x_{2}^{(k)} +
\alpha_{1\,3}\,x_{3}^{(k)} + \cdots + \alpha_{1\,n}\,x_{n}^{(k)} +
\beta_1$  
$\displaystyle x_{2}^{(k+1)}$ $\displaystyle =$ $\displaystyle \alpha_{2\,1}\,x_{1}^{(k+1)} +
\alpha_{2\,3}\,x_{3}^{(k)} + \cdots + \alpha_{2\,n}\,x_{n}^{(k)}+
\beta_2$  
$\displaystyle \cdots$      
$\displaystyle x_{n}^{(k+1)}$ $\displaystyle =$ $\displaystyle \alpha_{n\,1}\,x_{1}^{(k+1)} +
\alpha_{n\,2}\,x_{2}^{(k+1)} + \cdots +
\alpha_{n\,n-1}\,x_{{n-1}}^{(k+1)} + \beta_n,$  

gdje je $ \alpha{}_{ij}=-\frac{a_{ij}}{a_{ii}}, \beta{}_i=\frac{b_i}{a_{ii}}.$

Primjer 3.10   Riješiti sljedeći sustav jednadžbi
$\displaystyle 1.00\,x + 2.51\,y + 3.72\,z$ $\displaystyle =$ $\displaystyle 5.12$  
$\displaystyle 3.15\,x - 0.20\,y - 1.97\,z$ $\displaystyle =$ $\displaystyle 0.33$  
$\displaystyle 4.43\,x + 5.84\,y + 1.79\,z$ $\displaystyle =$ $\displaystyle 3.45$  

Rješenje. Radi jednostavnosti ćemo radije vektorstupac rješenja pisati u obliku $ \{x,y,z\}.$ Ako počnemo s $ \{1,1,1\},$ Jacobijeva metoda daje redom sljedeće aproksimacije.

$\displaystyle \{ 1,1,1\},\quad \{ -1.11,4.25,-3.81006\},\quad
\{ 8.62591,18.3966,-9.19145\},$

$\displaystyle \{ -6.86314,224.744,-79.4406\},\quad
\{ -263.468,672.745,-714.33\},\quad
\{ 973.837,2884.88,-1540.9\},$

$\displaystyle \{ -32614.,92744.1,-95831.1\},\quad
\{ 123709.,430265.,-221867.\},$

$\displaystyle \{ -254614.,4.13381\times {{10}^6},-1.70993\times {{10}^6}\},\quad
\{ -4.01492\times {{10}^6},1.28326\times {{10}^7},
-1.28567\times {{10}^7}\}.$

Gauss-Seidelova metoda daje

$\displaystyle \{ 1,1,1\},\quad \{ -1.11,-28.9825,99.2319\},\quad
\{ -291.277,-5566.69,18884.5\},$

$\displaystyle \{ -56272.9,-1.07231\times {{10}^6},3.63776\times {{10}^6}\},\quad
\{ -1.0841\times {{10}^7},-2.06577\times {{10}^8},
7.00802\times {{10}^8}\}.$

Međutim, ako se promijeni poredak jednadžbi, tako da prva jednadžba dođe na treće mjesto, druga na prvo, a treća na drugo, tj. ako sustav napišemo u obliku
$\displaystyle 3.15\,x - 0.20\,y - 1.97\,z$ $\displaystyle =$ $\displaystyle 0.33$  
$\displaystyle 4.43\,x + 5.84\,y + 1.79\,z$ $\displaystyle =$ $\displaystyle 3.45$  
$\displaystyle 1.00\,x + 2.51\,y + 3.72\,z$ $\displaystyle =$ $\displaystyle 5.12$  

onda Jacobijevom metodom dobivamo

$\displaystyle \{ 1,1,1\},\quad \{ 0.793651,-0.474315,0.432796\},\quad
\{ 0.345316,-0.143934,1.48303\},$

$\displaystyle \{ 1.02311,-0.125749,1.38063\},\quad
\{ 0.960222,-0.60851,1.18616\},\quad
\{ 0.807949,-0.501201,1.5288\},$

$\displaystyle \{ 1.02905,-0.490713,1.49733\},\quad
\{ 1.01003,-0.648784,1.43082\},\quad
\{ 0.958398,-0.613973,1.54259\},$

$\displaystyle \{ 1.03051,-0.609064,1.53298\},$

a Gauss-Seidelovom

$\displaystyle \{ 1,1,1\},\quad \{ 0.793651,-0.317786,1.37742\},\quad
\{ 0.946018,-0.549047,1.4925\},$

$\displaystyle \{ 1.0033,-0.627776,1.53022\},\quad
\{ 1.0219,-0.653441,1.54254\},\quad
\{ 1.02797,-0.661825,1.54656\},$

$\displaystyle \{ 1.02996,-0.664563,1.54788\},\quad
\{ 1.0306,-0.665458,1.54831\},\quad
\{ 1.03082,-0.66575,1.54845\},$

$\displaystyle \{ 1.03088,-0.665845,1.54849\}.$

Inače, rješenje na pet decimala točno, je

$\displaystyle \{ 1.03092,-0.665892,1.54851\}.$


OR (overrelaxation) metode

Navedene metode se mogu poboljšati, obzirom na konvergentnost i brzinu konvergencije, na sljedeći način.

Matrica $ A$ se rastavi kako slijedi

$\displaystyle A = B(\omega{}) - C(\omega{}),$

gdje su $ B(\omega{})$ i $ C(\omega{})$ matrice tipa $ (n,n)$ ovisne o parametru $ \omega{}.$ Osim toga $ B(\omega{})$ treba biti regularna. Parametar $ \omega{}$ se zove parametar relaksacije. Tada je

$\displaystyle B(\omega{})\,\boldsymbol{x} = C(\omega{})\,\boldsymbol{x} + \boldsymbol{b}.$ (3.11)

Zbog regularnosti $ B(\omega{}),$ postoji $ B(\omega{})^{-1}.$ Stavimo

$\displaystyle G(\omega{}) = B(\omega{})^{-1}\,C(\omega{}),\hspace{1cm}
\boldsymbol{g}(\omega{}) = B(\omega{})^{-1}\,\boldsymbol{b}.$

Uvrstimo u (3.11), dobivamo

$\displaystyle \boldsymbol{x} = G(\omega{})\,\boldsymbol{x} +
\boldsymbol{g}(\omega{}).$

Iterativni postupak je tada dan formulom

$\displaystyle \boldsymbol{x}^{(k+1)} = G(\omega{})\,\boldsymbol{x}^{(k)} + \boldsymbol{g}(\omega{}),\hspace{1cm}k=0,1,2,\ldots{}\ .$ (3.12)

Parametar relaksacije treba odabrati tako da postupak što brže konvergira.


Jacobijeva OR-metoda (JOR metoda)

Imamo

$\displaystyle A = L + D + R.$

Neka je $ D$ regularna matrica. Stavimo

$\displaystyle B(\omega{}) = \frac{1}{\omega{}}\,D.$

Tada je

$\displaystyle B(\omega{})^{-1} = \omega{}\,D^{-1},$

$\displaystyle C(\omega{}) = \left(\frac{1-\omega{}}{\omega{}}\right)\,D - (L + R),$

$\displaystyle G(\omega{}) = B(\omega{})^{-1}\,C(\omega{}) = (1 - \omega{})\,I -
\omega{}\,D^{-1}\,(L + R),$

$\displaystyle \boldsymbol{g}(\omega{}) = B(\omega{})^{-1}\,\boldsymbol{b} =
\omega{}\,D^{-1}\,\boldsymbol{b}.$

Tako imamo sljedeći algoritam

Algoritam 9   (Jacobijeva OR-metoda) Proizvoljno izaberemo početnu aproksimaciju

% latex2html id marker 39048
$\displaystyle \boldsymbol{x}^{(0)}=\left[
\begin{array}{c}
x_1^{(0)} \\  x_2^{(0)} \\  \vdots \\  x_n^{(0)}
\end{array}\right],$

i zatim računamo sljedeće aproksimacije $ \boldsymbol{x}^{(n)}$ po formuli

$\displaystyle \boldsymbol{x}^{(k+1)} = \left[(1 - \omega{})\,I -
\omega{}\,D^{-1}\,(L + R)\right]\,\boldsymbol{x}^{(k)} +
\omega{}\,D^{-1}\,\boldsymbol{b},$

odnosno
$\displaystyle x_{1}^{(k+1)}$ $\displaystyle =$ $\displaystyle (1-\omega)\,x_{1}^{(k)} +
\alpha_{1\,2}\,x_{2}^{(k)} +
\alpha_{1\,3}\,x_{3}^{(k)} + \cdots +
\alpha_{1\,n}\,x_{n}^{(k)} +
\beta_1$  
$\displaystyle x_{2}^{(k+1)}$ $\displaystyle =$ $\displaystyle (1-\omega)\,x_{2}^{(k)} +
\alpha_{2\,1}\,x_{1}^{(k)} +
\alpha_{2\,3}\,x_{3}^{(k)} + \cdots +
\alpha_{2\,n}\,x_{n}^{(k)} +
\beta_2$  
    $\displaystyle \vdots$  
$\displaystyle x_{n}^{(k+1)}$ $\displaystyle =$ $\displaystyle (1-\omega)\,x_{n}^{(k)} +
\alpha_{n\,1}\,x_{1}^{(k)} +
\alpha_{n\,2}\,x_{2}^{(k)} + \cdots +
\alpha_{n\,n-1}\,x_{{n-1}}^{(k)} +
\beta_n,$  

gdje je $ \alpha_{i\,j}=-\frac{\omega}{\alpha_{i\,i}}\,a_{i\,j},$ i $ \beta_i=\frac{\omega}{\alpha_{i\,i}}\,b_i.$

Primijetimo da za $ \omega{}=1$ Jacobijeva OR metoda postaje Jacobijeva metoda.


Gauss-Seidelova OR-metoda (SOR metoda)

Imamo

$\displaystyle A = L + D + R.$

Neka je $ D$ regularna matrica. Stavimo

$\displaystyle B(\omega{}) = \frac{1}{\omega{}}\,D + L = \frac{1}{\omega{}}\,(D +
\omega{}L).$

Tada je

$\displaystyle B(\omega{})^{-1} = \omega{}\,(D + \omega{}L)^{-1},$

$\displaystyle C(\omega{}) = \left(\frac{1-\omega{}}{\omega{}}\right)\,D - R,$

$\displaystyle G(\omega{}) = B(\omega{})^{-1}\,C(\omega{}) = (D +
\omega{}L)^{-1}\,\left[(1 - \omega{})\,D - \omega{}\,R\right],$

$\displaystyle \boldsymbol{g}(\omega{}) = B(\omega{})^{-1}\,\boldsymbol{b} =
\omega{}\,(D + \omega{}L)^{-1}\,\boldsymbol{b}.$

Tako je

$\displaystyle \boldsymbol{x}^{(k+1)} = (D + \omega{}L)^{-1}\,\left[(1 -
\omega{...
...right]\,\boldsymbol{x}^{(k)} +
\omega{}\,(D + \omega{}L)^{-1}\,\boldsymbol{b},$

Računanje $ (D + \omega{}L)^{-1}$ nije jednostavno, pa radimo nešto drukčije. Polazimo od prethodne jednadžbe

$\displaystyle B(\omega{})\,\boldsymbol{x}^{(k+1)} =
C(\omega{})\,\boldsymbol{x}^{(k)} + \boldsymbol{b},$

odnosno

$\displaystyle \frac{1}{\omega{}}\,(D +
\omega{}L)\,\boldsymbol{x}^{(k+1)} =
\le...
...omega{}}{\omega{}}\right)\,D -
R\right]\,\boldsymbol{x}^{(k)} + \boldsymbol{b}.$

Pomnožimo ovu jednadžbu s $ \omega{}\,D^{-1}$

$\displaystyle (I + \omega{}D^{-1}L)\,\boldsymbol{x}^{(k+1)} =
\left[\left(1-\om...
...
\omega{}\,D^{-1}R\right]\,\boldsymbol{x}^{(k)} +
\omega{}D^{-1}\boldsymbol{b}.$

Tako imamo sljedeći algoritam

Algoritam 10   (Gauss-Seidelova OR-metoda) Proizvoljno izaberemo početnu aproksimaciju

% latex2html id marker 39109
$\displaystyle \boldsymbol{x}^{(0)}=\left[
\begin{array}{c}
x_1^{(0)} \\  x_2^{(0)} \\  \vdots \\  x_n^{(0)}
\end{array}\right],$

i zatim računamo sljedeće aproksimacije $ \boldsymbol{x}^{(n)}$ po formuli

$\displaystyle \boldsymbol{x}^{(k+1)} =
\left(1-\omega{}\right)\,\boldsymbol{x}^...
...{(k+1)} -
\omega{}D^{-1}R\,\boldsymbol{x}^{(k)} + \omega{}D^{-1}\boldsymbol{b},$

odnosno
$\displaystyle x_{1}^{(k+1)}$ $\displaystyle =$ $\displaystyle (1-\omega)\,x_{1}^{(k)} +
\alpha_{1\,2}\,x_{2}^{(k)} +
\alpha_{1\,3}\,x_{3}^{(k)} + \cdots +
\alpha_{1\,n}\,x_{n}^{(k)} +
\beta_1$  
$\displaystyle x_{2}^{(k+1)}$ $\displaystyle =$ $\displaystyle (1-\omega)\,x_{2}^{(k)} +
\alpha_{2\,1}\,x_{1}^{(k+1)} +
\alpha_{2\,3}\,x_{3}^{(k)} + \cdots +
\alpha_{2\,n}\,x_{n}^{(k)} +
\beta_2$  
    $\displaystyle \vdots$  
$\displaystyle x_{n}^{(k+1)}$ $\displaystyle =$ $\displaystyle (1-\omega)\,x_{n}^{(k)} +
\alpha_{n\,1}\,x_{1}^{(k+1)} +
\alpha_{n\,2}\,x_{2}^{(k+1)} + \cdots +
\alpha_{n\,n-1}\,x_{{n-1}}^{(k+1)} +
\beta_n.$  

gdje je $ \alpha_{i\,j}=-\frac{\omega}{\alpha_{i\,i}}\,a_{i\,j},$ i $ \beta_i=\frac{\omega}{\alpha_{i\,i}}\,b_i.$

Primijetimo da za $ \omega{}=1$ SOR metoda postaje Gauss-Seidelova metoda.

Broj $ \omega{}$ treba uzeti tako da je $ 0<\omega{}<2.$ Optimalna vrijednost ovisi o konkretnom problemu koji se rješava.

Primjer 3.11   Riješiti sustav u primjeru 3.10 OR metodama.

Rješenje. Jacobijeva OR metoda ne daje ništa bolje rezultate kad se uzme $ \omega\neq 1.$ Gauss-Seidelova OR metoda za $ \omega=1.075$ daje

$\displaystyle \{ 1,1,1\},\quad \{ 0.778175,-0.404,1.47273\},\quad
\{ 1.0168,-0.649051,1.54606\}$

$\displaystyle \{ 1.03148,-0.666804,1.5492\},\quad
\{ 1.03127,-0.666338,1.54868\},\quad
\{ 1.03098,-0.665961,1.54853\}$

$\displaystyle \{ 1.03092,-0.665897,1.54851\},\quad
\{ 1.03092,-0.665892,1.54851\},\quad
\{ 1.03092,-0.665892,1.54851\}$

$\displaystyle \{ 1.03092,-0.665892,1.54851\}.$


Konvergencija iterativne metode

Iterativni postupci, koje smo upravo definirali, opisani su formulom

$\displaystyle \boldsymbol{x}^{(k+1)} = G\,\boldsymbol{x}^{(k)} +
\boldsymbol{g},$

gdje je $ G$ matrica tipa $ (n,n),$ definirana od slučaja do slučaja na različite načine pomoću matrice $ A.$ Razmotrimo sada probleme konvergencije navedenih metoda.

Definicija 16   Neka je $ A$ simetrična matrica. Kažemo da je $ A$ pozitivno definitna, ako vrijedi

$\displaystyle A\,\boldsymbol{x} \cdot{} \boldsymbol{x} > 0,\hspace{1cm}
\forall{}\boldsymbol{x} \neq{} \textbf{0}.$

Lema 3   Neka je $ G$ simetrična pozitivno definitna matrica tipa $ (n,n).$ Tada su vlastite vrijednosti $ \lambda{}_i(G),\,i=1,2,\ldots{},n$ od matrice $ G$ pozitivne, i vrijedi

$\displaystyle \Vert G\,\boldsymbol{x}\Vert \leqslant{} \max{}_i \lambda{}_i(G)
\Vert\boldsymbol{x}\Vert,$

za svaki $ \boldsymbol{x}\in {\cal R}_{n}.$


Dokaz. Budući da je matrica $ G$ simetrična, postoji $ n$ međusobno okomitih i jediničnih vektora $ \boldsymbol{x}_1,\boldsymbol{x}_2,\ldots,\boldsymbol{x}_n$ u $ {\cal R}_n,$ i brojevi $ \lambda{}_1,\lambda{}_2,\ldots{},\lambda{}_n$ takvih, da je

$\displaystyle G\,\boldsymbol{x}_i = \lambda{}_i\,\boldsymbol{x}_i.$

Zbog pozitivne definitnosti

$\displaystyle G\,\boldsymbol{x}_i \cdot{} \boldsymbol{x}_i =
\lambda{}_i\,\boldsymbol{x}_i \cdot{} \boldsymbol{x}_i = \lambda{}_i >
0,\hspace{1cm}\forall{}i.$

Dokažimo sada drugu tvrdnju. Primijetimo najprije da vlastiti vektori $ \boldsymbol{x}_1,\boldsymbol{x}_2,\ldots,\boldsymbol{x}_n$ matrice $ G$ čine ortonormiranu bazu u $ {\cal R}_{n}.$ To znači da se bilo koji vektor $ \boldsymbol{x}\in{\cal R}_{n}$ može prikazati kao njihova linearna kombinacija.

$\displaystyle \boldsymbol{x} = \alpha{}_1\,\boldsymbol{x}_1 +
\alpha{}_2\,\boldsymbol{x}_2 + \cdots{} +
\alpha{}_n\,\boldsymbol{x}_n.$

Zbog linearnosti djelovanja matrice $ G$ na vektore

$\displaystyle G\,\boldsymbol{x} = \alpha{}_1\,G\,\boldsymbol{x}_1 +
\alpha{}_2\...
...ambda_2\,\boldsymbol{x}_2 + \cdots{} +
\alpha{}_n\,\lambda_n\,\boldsymbol{x}_n.$

Tako je

$\displaystyle \Vert G\,\boldsymbol{x}\Vert^2 = G\,\boldsymbol{x} \cdot{}
G\,\bo...
...,\lambda_1^2 +
\alpha{}_2^2\,\lambda_2^2 + \cdots{} +
\alpha{}_n^2\,\lambda_n^2$

$\displaystyle \leqslant{}
\left(\max{}_i\lambda{}_i\right)^2\, \left(\alpha{}_1...
...{}_n^2\right) =
\left(\max{}_i\lambda{}_i\right)^2\,\Vert\boldsymbol{x}\Vert^2.$

Odatle

$\displaystyle \Vert G\,\boldsymbol{x}\Vert \leqslant{}
\max{}_i\lambda{}_i\,\Vert\boldsymbol{x}\Vert.$

$ \heartsuit$

Teorem 25   Neka je $ G$ simetrična, pozitivno definitna matrica. Iterativni postupak

$\displaystyle \boldsymbol{x}^{(k+1)} = G\,\boldsymbol{x}^{(k)} +
\boldsymbol{g},\hspace{1cm}k=0,1,2,\ldots{}$

konvergira za bilo koju početnu aproksimaciju $ \boldsymbol{x}^{(0)},$ ako i samo ako je

$\displaystyle \max{}_i \lambda{}_i(G) = L < 1.$


Dokaz. 1. Dokažimo najprije da iz

$\displaystyle \max{}_i \lambda{}_i(G) = L < 1$

slijedi konvergencija iterativnog postupka. Doista

$\displaystyle \boldsymbol{x}^{(k+1)} - \boldsymbol{x}^{(k)} =
G\,(\boldsymbol{x}^{(k)} - \boldsymbol{x}^{(k-1)}),$

i odatle

$\displaystyle \Vert\boldsymbol{x}^{(k+1)} - \boldsymbol{x}^{(k)}\Vert =
\Vert G...
...l{x}^{(k-1)}\Vert = L\,\Vert\boldsymbol{x}^{(k)} -
\boldsymbol{x}^{(k-1)}\Vert.$

Odatle možemo, na isti način kao kod računanja apriorne ocjene greške kod metode iteracije (v. 3.2.2), zaključiti

$\displaystyle \Vert\boldsymbol{x}^{(m)}-\boldsymbol{x}^{(n)}\Vert \leqslant{}
\frac{L^n}{1-L}\,\Vert\boldsymbol{x}^{(1)}-\boldsymbol{x}^{(0)}\Vert.$

Budući da je $ 0<L<1,$

$\displaystyle \lim_{k\rightarrow{}\infty{}} L^k = 0,$

pa niz $ (\boldsymbol{x}^{(k)})$ konvergira. Neka je

$\displaystyle \lim_{k\rightarrow{}\infty}\boldsymbol{x}^{(k)} = \boldsymbol{s}.$

Tada je

$\displaystyle \lim_{k\rightarrow{}\infty{}} \boldsymbol{x}^{(k+1)} =
\lim_{k\rightarrow{}\infty{}} (G\,\boldsymbol{x}^{(k)} +
\boldsymbol{g}),$

$\displaystyle \boldsymbol{s} = G\,\boldsymbol{s} + \boldsymbol{g},$

i prema tome niz konvergira k rješenju $ \boldsymbol{s}.$
2. Dokažimo obrat, tj. da iz konvergencije postupka slijedi

$\displaystyle \max{}_i \lambda{}_i(G) = L < 1.$

Iz

$\displaystyle \boldsymbol{x}^{(k)} = G\,\boldsymbol{x}^{(k-1)} +
\boldsymbol{g}$

i

$\displaystyle \boldsymbol{s} = G\,\boldsymbol{s} +
\boldsymbol{g}$

slijedi

$\displaystyle \boldsymbol{x}^{(k)} - \boldsymbol{s} =
G\,(\boldsymbol{x}^{(k-1)...
...)} - \boldsymbol{s}) = \cdots{} =
G^k\,(\boldsymbol{x}^{(0)} - \boldsymbol{s}).$

Zbog pretpostavke

$\displaystyle \lim_{k\rightarrow{}\infty{}} (\boldsymbol{x}^{(k)} -
\boldsymbol{s}) = 0,$

slijedi

$\displaystyle \lim_{k\rightarrow{}\infty{}} G^k\,(\boldsymbol{x}^{(0)} -
\bolds...
...\rightarrow{}\infty{}}
G^k\right)\,(\boldsymbol{x}^{(0)} - \boldsymbol{s}) = 0.$

Ako se to događa za svaku početnu aproksimaciju $ \boldsymbol{x}^{(0)},$ onda mora biti

$\displaystyle \lim_{k\rightarrow{}\infty{}} G^k = O.$

Budući da je $ G$ simetrična matrica, ona se može dijagonalizirati, tj. postoji regularna matrica $ X$ takva, da vrijedi

% latex2html id marker 39275
$\displaystyle G^k = X\,\left[\begin{array}{cccc}
...
...ddots & \vdots \\
0 & 0 & \cdots & \lambda_n^k(G)
\end{array}\right]\,X^{-1},$   za $\displaystyle k=1,2,\ldots\,.$

Odatle

% latex2html id marker 39278
$\displaystyle \lim_{k\rightarrow{}\infty{}} G^k= \...
...dots & \vdots \\
0 & 0 & \cdots & \lambda_n^k(G)
\end{array}\right]\,X^{-1} =$

% latex2html id marker 39280
$\displaystyle = X\,\left(\lim_{k\rightarrow{}\inft...
...vdots \\
0 & 0 & \cdots & \lambda_n^k(G)
\end{array}\right]\right)\,X^{-1} = $

% latex2html id marker 39282
$\displaystyle = X\,
\left[\begin{array}{cccc}
\li...
...s & \lim_{k\rightarrow{}\infty{}}\lambda_n^k(G)
\end{array}\right]\,X^{-1} = O.$

Ako pomnožimo s lijeva s $ X^{-1},$ i s desna s $ X$, slijedi

% latex2html id marker 39288
$\displaystyle \left[\begin{array}{cccc}
\lim_{k\r...
... & \cdots & \lim_{k\rightarrow{}\infty{}}\lambda_n^k(G)
\end{array}\right] = O.$

Prema tome

$\displaystyle \lim_{k\rightarrow{}\infty{}}\lambda_i^k(G) = 0,\hspace{1cm}
i=1,2,\ldots{},n.$

Zbog $ \lambda_i(G)>0,$ to je moguće samo ako je

$\displaystyle \max{}_i \lambda{}_i(G) < 1.$

$ \heartsuit$

Na isti način kao kod metode iteracije prilikom traženja nultočke funkcije možemo naći apriornu ocjenu greške

$\displaystyle \Vert\boldsymbol{x}^{(k)} - \boldsymbol{s}\Vert \leqslant{}
\frac{L^k}{1-L}\,\Vert\boldsymbol{x}^{(1)} -
\boldsymbol{x}^{(0)}\Vert,$

i aposteriornu ocjenu greške

$\displaystyle \Vert\boldsymbol{x}^{(k)} - \boldsymbol{s}\Vert \leqslant{}
\frac{L}{1-L}\,\Vert\boldsymbol{x}^{(k)} -
\boldsymbol{x}^{(k-1)}\Vert.$

Ovi rezultati su dobiveni uz pretpostavku da je $ G$ simetrična pozitivno definitna matrica. To je važan slučaj, koji se pojavljuje prilikom približnog rješavanja rubnih i rubno-početnih problema varijacijskim metodama. Međutim, može se dogoditi da $ G$ ne ispunjava ovaj uvjet. Ipak, dobiveni rezultati vrijede i tada uz neke modifikacije.

Prije svega, budući da u općem slučaju vlastite vrijednosti mogu biti kompleksni brojevi, broj $ L=\max{}_i\lambda{}_i(G)$ više nema smisla, jer u skupu kompleksnih brojeva nemamo prirodni uređaj. Umjesto toga imamo ovu definiciju.

Definicija 17   Neka je $ G$ matrica tipa $ (n,n),$ i neka su $ \lambda{}_1,\lambda{}_2,\ldots{},\lambda{}_n$ njezine vlastite vrijednosti. Broj

$\displaystyle \nu(G) = \max{}_i\,\vert\lambda{}_i\vert$

se zove spektralni radius matrice $ G.$

Očito je u slučaju simetričnosti i pozitivne definitnosti matrice $ G$

$\displaystyle \nu(G) = \max{}_i \lambda{}_i(G).$

Može se dokazati sljedeći teorem

Teorem 26   Iterativni postupak

$\displaystyle \boldsymbol{x}^{(k+1)} = G\,\boldsymbol{x}^{(k)} +
\boldsymbol{g},\hspace{1cm}k=0,1,2,\ldots{}$

konvergira za bilo koju početnu aproksimaciju $ \boldsymbol{x}^{(0)},$ ako i samo ako je

$\displaystyle \nu(G) < 1.$

Tako se problem konvergencije iterativnog postupka svodi na problem ocjene spektralnog radiusa, odnosno najveće po modulu (apsolutnoj vrijednosti) vlastite vrijednosti matrice, i to ocjene odozgo. Ima raznih ocjena te vrste. Navedimo neke od njih. Neka je $ G=[g_{i\,j}].$ Tada vrijedi

1.
$ \nu(G)\leqslant
\left[\sum_{i=1}^n\sum_{j=1}^n\vert g_{ij}\vert^2\right]^{\frac{1}{2}},$
2.
$ \nu(G)\leqslant
\max{}\{\sum_{j=1}^n\vert g_{ij}\vert;\;i=1,2,\ldots,n\},$
3.
$ \nu(G)\leqslant \max{}\{\sum_{i=1}^n\vert g_{ij}\vert;\;j=1,2,\ldots,n\},$

Ako je bilo koji od ovih brojeva na desnoj strani manji od $ 1,$ spektralni radius je nužno manji od $ 1.$ Naglasimo da se ovdje ne radi o spektralnom radiusu matrice $ A$ već matrice $ G.$

Primjer 3.12   Ispitati konvergenciju Jacobijeve i Gauss-Seidelove metode u primjeru 3.10.

Rješenje. Kod Jacobijeve metode imamo

$\displaystyle G_J = - D^{-1}(L+R),$

a kod Gauss-Seidelove

$\displaystyle G_S = -(L + D)^{-1}\,R.$

Matrica sustava, kako je na početku napisan, je

% latex2html id marker 39354
$\displaystyle A = \left[
\begin{array}{rrr}
1.00 &...
...& 3.72 \\
3.15 & - 0.20 & - 1.97 \\
4.43 & 5.84 & 1.79
\end{array}\right].$

Ako primijenimo Jacobijevu metodu, imamo

% latex2html id marker 39356
$\displaystyle G_J = \left[\begin{array}{rrr}
0 & -...
...& -3.72 \\
15.75 & 0 & -9.85 \\
-2.47486 & -3.26257 & 0
\end{array}\right].$

Vlastite vrijednosti su $ 5.18447,-2.59224 + 4.28355\,i,
-2.59224 - 4.28355\,i,$ pa je spektralni radius $ \nu(G_J)=5.18447.$

Kod Gauss-Seidelove metode imamo još lošiji rezultat

% latex2html id marker 39362
$\displaystyle G_S = \left[\begin{array}{rrr}
0 & -...
...-3.72 \\
0 & -39.5325 & -68.44 \\
0 & 135.189 & 232.497
\end{array}\right].$

Vlastite vrijednosti su $ 192.647,0.317614,0,$ pa je spektralni radius $ \nu(G_S)=192.647.$

Kad prvu jednadžbu stavimo na treće mjesto, Jacobijeva metoda daje

% latex2html id marker 39368
$\displaystyle G_J = \left[
\begin{array}{rrr}
0 & ...
...
-0.758562 & 0 & -0.306507 \\
-0.268817 & -0.674731 & 0
\end{array}\right],$

vlastite vrijednosti su $ -0.341587 + 0.599596\,i,
-0.341587 - 0.599596\,i,0.683174,$ pa je spektralni radius $ \nu(G_J)=0.69007.$

% latex2html id marker 39374
$\displaystyle G_S = \left[
\begin{array}{rrr}
0 & ...
...
0 & -0.0481626 & -0.780909 \\
0 & 0.0154291 & 0.358786
\end{array}\right],$

vlastite vrijednosti su $ 0.326639,-0.0160158,0,$ pa je spektralni radius $ \nu(G_S)=0.326639.$ Ako se radi Gauss-Seidelovom OR metodom, spektralni radius postaje $ \nu(G_S)=0.126947$ za $ \omega=1.075.$

Primjer 3.13   Neka je matrica sustava

% latex2html id marker 39385
$\displaystyle A = \left[
\begin{array}{rrrr}
2 & -...
... & 2 & -1 & 0 \\
0 & -1 & 2 & -1 \\
0 & 0 & -1 & 2 \\
\end{array}\right].$

Treba ispitati da li Jacobijev i Gauss-Seidelov postupak konvergiraju.

Rješenje. U Jacobijevom postupku je

% latex2html id marker 39387
$\displaystyle G = - D^{-1}(L+R) = \left[
\begin{ar...
... & 1/2 & 0 \\
0 & 1/2 & 0 & 1/2 \\
0 & 0 & 1/2 & 0 \\
\end{array}\right].$

Njezine vlastite vrijednosti su

$\displaystyle -0.809017,\;\;-0.309017,\;\;0.309017,\;\;0.809017,$

spektralni radius je

$\displaystyle \nu(G) = 0.809017 < 1,$

i prema tome Jacobijeva metoda konvergira za svaku početnu aproksimaciju.

Kod Gauss-Seidelove metode je

% latex2html id marker 39393
$\displaystyle G = -(L + D)^{-1}\,R = \left[
\begin...
...& 0 \\
0 & 1/8 & 1/4 & 1/2 \\
0 & 1/16 & 1/8 & 1/4 \\
\end{array}\right].$

Njezine vlastite vrijednosti su

$\displaystyle 0,\hspace{1cm}0,\hspace{1cm}0.0954915,\hspace{1cm}0.654508,$

spektralni radius je sada

$\displaystyle \nu(G) = 0.654508 < 1,$

pa prema tome Gauss-Seidelova metoda također konvergira za svaku početnu aproksimaciju, i to brže nego Jacobijeva.

S porastom reda matrice, uz isti tridijagonalni oblik, vlastite vrijednosti matrice $ G$ rastu prema $ 1.$ U slučajevima kad su vrlo blizu broju $ 1,$ OR metode postaju vrlo efikasne u smislu brzine konvergencije.


next up previous contents index
Next: Aproksimacija funkcije i numerička Up: Numerička matematika Previous: Rješavanje jednadžbi   Contents   Index
Salih Suljagic
1999-12-17