Metoda konačnih diferencija

Neka je dan rubni problem

$$% latex2html id marker 40627 \begin{cases} u''(x) - q(x)\,u(... ...x) = 0,& \\ u(0) = 0, \hspace{1cm}u'(\ell) = 0. & \end{cases}$$

Podijelimo segment $[0,\ell{}]$ na $n$ jednakih podsegmenata.

$$0 = x_0 < x_1 < x_2 < \cdots < x_n = \ell{}.$$

Točke $x_0,x_1,x_2,\ldots,x_n$ se zovu čvorovi. Točke $x_1,x_2,\ldots,x_{n-1}$ se zovu unutrašnji čvorovi, a $x_0,x_n$ se zovu rubni čvorovi. Svaki segment ima duljinu $h=\frac{\ell{}}{n}.$ Broj $h$ se zove korak. On se izabire malen, da točnost bude veća, ali sa smanjivanjem koraka se javljaju drugi nepoželjni efekti, pa u odabiru koraka treba biti oprezan. Ideja metode se sastoji u tome da se u svakom čvoru diferencijalna jednadžba zamijeni odgovarajućom algebarskom jednadžbom i zatim riješi tako dobiveni sustav algebarskih jednadžbi. U tu svrhu treba derivacije zamijeniti odgovarajućim algebarskim aproksimacijama. Rješenje koje tako dobijemo predstavlja približne vrijednosti rješenja u čvorovima.

Uočimo jedan unutrašnji čvor $x_i.$ Koristit ćemo sljedeće oznake

$$u(x_i) = u_i,\hspace{1cm}u'(x_i) = u'_i,\hspace{1cm}u''(x_i) = u''_i.$$

Na isti način kao u 3.5.1 imamo

$$u_{i+1} \approx u_i + u'_i\,h.$$

Odatle

$$u'_i \approx \frac{u_{i+1}-u_i}{h}.$$

Ova aproksimacija derivacije se zove aproksimacija s desna.

Umjesto aproksimacije s desna ponekad se koriste aproksimacija s lijeva ili centralna aproksimacija. Polazeći od formule (3.24), u kojoj zamjenimo $h$ s $-h,$

$$u_{i-1} = u(x_i-h) = u_i - u'_i\,h + \frac{1}{2}\,u''(x_i-s\,h)\,h^2,$$

dobivamo, zanemarivanjem člana s $h^2,$ aproksimaciju s lijeva

$$u'_i \approx \frac{u_i-u_{i-1}}{h}.$$

Centralnu aproksimaciju dobivamo, ako Taylorove formule

$$u_{i+1} = u(x_i+h) = u_i + u'_i\,h + \frac{1}{2}\,u''_i\,h^2 + \frac{1}{3!}\,u'''(x_i+t\,h)\,h^3,$$

$$u_{i-1} = u(x_i-h) = u_i - u'_i\,h + \frac{1}{2}\,u''_i\,h^2 - \frac{1}{3!}\,u'''(x_i-s\,h)\,h^3.$$

za $h$ i $-h$ oduzmemo i podijelimo s $2\,h.$ Tada imamo

$$u'_i = \frac{u_{i+1}-u_{i-1}}{2\,h} + \frac{h^2}{3!}\,\left[u'''(x_i+t\,h) + u'''(x_i-s\,h)\right],$$

pa ako zanemarimo član s $h^2,$ dobivamo

$$u'_i \approx \frac{u_{i+1}-u_{i-1}}{2\,h}.$$

Za aproksimaciju druge derivacije, uzimamo

$$u_{i+1} = u(x_i+h) = u_i + u'_i\,h + \frac{1}{2}\,u''_i\,h^2 + \frac{1}{3!}\,u'''_i\,h^3 + \frac{1}{4!}\,u^{iv}(x_i+t\,h)\,h^4,$$

$$u_{i-1} = u(x_i-h) = u_i - u'_i\,h + \frac{1}{2}\,u''_i\,h^2 - \frac{1}{3!}\,u'''_i\,h^3 + \frac{1}{4!}\,u^{iv}(x_i-s\,h)\,h^4.$$

Zbrojimo ove jednakosti

$$u_{i+1} + u_{i-1} = 2\,u_i + u''_i\,h^2 + \frac{h^4}{4!}\,\left[u^{iv}(x_i+t\,h) + u^{iv}(x_i-s\,h)\right].$$

Podijelimo s $h^2,$ i izračunamo $u''_i$

$$u''_i = \frac{u_{i+1}-2\,u_i+u_{i-1}}{h^2} - \frac{h^2}{4!}\,\left[u^{iv}(x_i+t\,h) + u^{iv}(x_i-s\,h)\right].$$

Zanemarimo zadnji član s faktorom $h^2,$ i dobivamo

$$u''_i \approx{}\frac{u_{i+1}-2\,u_i+u_{i-1}}{h^2}.$$

U čvoru $x=x_i$ diferencijalna jednadžba iz rubnog problema glasi

$$u''_i - q(x_i)\,u_i + f(x_i) = 0.$$

Zamijenimo li drugu derivaciju s njezinom aproksimacijom, dobijemo algebarsku jednadžbu za $i$-ti čvor

$$\frac{u_{i+1}-2\,u_i+u_{i-1}}{h^2} - q_i\,u_i + f_i = 0.$$

Pomnožimo jednadžbu s $-h^2,$ pa imamo

$$-u_{i-1} + (2+h^2\,q_i)\,u_i - u_{i+1} = h^2\,f_i,\hspace{1cm} i=1,2,\ldots,n-1.$$

U rubnim čvorovima imamo zadan rubni uvjet. Tako je na lijevom rubu

$$u_0 = 0,$$ (3.25)
a na desnom

$$\frac{u_n-u_{n-1}}{h} = 0,$$
tj.

$$u_n = u_{n-1}.$$ (3.26)

Na taj način smo dobili sustav od $n-1$-ne linearne algebarske jednadžbe. Rubni uvjeti imaju utjecaja samo na prvu jednadžbu $(i=1),$ koja postaje, zbog (3.25),

$$(2+h^2\,q_1)\,u_1 - u_{2} = h^2\,f_1,$$

i na zadnju jednadžbu $(i=n-1),$ koja, zbog (3.27), postaje

$$-u_{n-2} + (1+h^2\,q_{n-1})\,u_{n-1} = h^2\,f_{n-1}.$$

Sustav možemo zapisati u matričnom obliku. Ako uzmemo prirodan poredak jednadžbi počevši od čvora $1$ preko $2$ sve do $n-1,$ imamo vektor nepoznanica

$$% latex2html id marker 40726 u = \left[ \begin{array}{cccc} u_1 & u_2 & \cdots{} & u_{n-1} \end{array} \right]^T,$$

matricu sustava

$$% latex2html id marker 40728 A(h) = \left[ \begin{array}{cc... ...0 & 0 & \cdots{} & -1 & (1+h^2\,q_{n-1}) \end{array} \right],$$

i vektor desne strane

$$% latex2html id marker 40730 b(h) = \left[ \begin{array}{cc... ...^2 & f_2\,h^2 & \cdots{} & f_{n-1}\,h^2 \end{array} \right]^T$$

Tada sustav možemo zapisati matrično

$$A(h)\,u = b(h).$$

Za kvadratnu matricu $A=[a_{ij}]$ kažemo da je striktno dijagonalno dominantna ako je

$$\sum_{\substack{i=0\\ i\neq k}}^n\vert a_{i\,j}\vert< \vert a_{i\,i}\vert,\qquad i=1,2,\ldots,n.$$

Može se pokazati da iz striktne dijagonalne dominantnosti matrice $A$ slijedi njezina regularnost. Ako je $q(x)>0,\forall x,$ onda je matrica $A(h)$ striktno dijagonalno dominantna, pa je regularna. Prema tome postoji jedinstveno rješenje.

Primjer 3.19   Riješiti metodom konačnih diferencija sljedeći rubni problem

$$% latex2html id marker 40745 \begin{cases} \ln (1 + x)\,y(x)... ...) + y''(x) = 0,& \\ y(0) = 1,\qquad y'(1) = 0. & \end{cases}$$

Rješenje. Podijelimo segment $[0,1]$ na deset jednakih dijelova, tako da je korak $h=0.1.$ Čvorovi su $x_i=i\,h,\;i=0,1,2,\ldots,10.$ Stavimo $y(x_i)=y_i.$ Za prvu derivaciju upotrebimo desnu aproksimaciju, osim na desnom rubu, gdje uzmemo lijevu. U svakom unutrašnjem čvoru zamijenimo derivacije odgovarajućim aproksimacijama. Dobivamo sustav jednadžbi

$$y_0 = 1$$

$$100\,y_0 - 200.905\,y_1 + 101\,y_2 = 0$$

$$100\,y_1 - 201.818\,y_2 + 102\,y_3 = 0$$

$$100\,y_2 - 202.738\,y_3 + 103\,y_4 = 0$$

$$100\,y_3 - 203.664\,y_4 + 104\,y_5 = 0$$

$$100\,y_4 - 204.595\,y_5 + 105\,y_6 = 0$$

$$100\,y_5 - 205.53\,y_6 + 106\,y_7 = 0$$

$$100\,y_6 - 206.469\,y_7 + 107\,y_8 = 0$$

$$100\,y_7 - 207.412\,y_8 + 108\,y_9 = 0$$

$$100\,y_8 - 208.358\,y_9 + 109\,y_{10} = 0$$

$$-10\,y_9 + 10\,y_{10} = 0$$

čije rješenje je

$${{y_0} = {1.}},\quad {{y_1} = {1.05352}},\quad {{y_2} = {1.10551}},\quad {{y_3} = {1.1545}},$$

$${{y_4} = {1.19913}},\quad {{y_5} = {1.23816}},\quad {{y_6} = {1.27055}},\quad {{y_7} = {1.29548}},$$

$${{y_8} = {1.31235}},\quad {{y_9} = {1.32083}},\quad {{y_{10}} = {1.32083}}.$$

Točke rješenja na sljedećoj slici spojene su ravnim linijama.