Metoda konačnih elemenata

Metoda konačnih elemenata je modifikacija Ritzove metode u tom smislu da što više elemenata matrice krutosti bude jednako nuli. Budući da su elementi matrice krutosti

$$K_{i\,j} = \int_0^{\ell}\,p\,v'_i\,v'_j\,dx,$$

funkcije $v_i$ treba birati tako da međusobni produkti njihovih derivacija budu nulfunkcije u što više slučajeva.

Pogledajmo na primjeru kako se to radi. Neka je dan rubni problem

$$% latex2html id marker 40967 \begin{cases} (p(x)\,u'(x))' +f(x) = 0,& \\ u(0) = 0,\qquad u'(\ell{}) = 0.& \end{cases}$$

Podijelimo segment $[0,\ell]$ na $n$ jednakih dijelova. Time smo dobili čvorove

$$x_i = i\,h,\qquad i=0,1,2,\ldots,n,$$

gdje je

$$h = \frac{\ell}{n}$$

korak diskretizacije. Rješenje tražimo u obliku

$$u_n = \sum_{i=0}^n c_i\,v_i,$$

gdje su $v_i$ koordinatne funkcije. Njih definiramo tako da je $v_i(x_j) = \delta_{ij},$ dakle da koordinatna funkcija $v_i$ ima u čvoru $x_i$ vrijednost $1,$ a u ostalim čvorovima vrijednost $0.$ Između čvorova graf funkcije $v_i$ dopunimo dijelovima pravaca. Odatle slijedi da su one definirane formulama

$$v_i(x) = \begin{cases}\frac{x}{h} - ... ...je $x_i\leqslant{}x\leqslant{}x_{i+1}$} \\ [1mm] 0,& \text{inače}, \end{cases}$$ (3.28)

za $i=1,2,\ldots,n-1.$ Koordinatne funkcije $v_0$ i $v_n$ su definirane ovako

$$v_0(x) = \begin{cases}-\frac{x}{h} +... ...o je $x_{0}\leqslant{}x\leqslant{}x_1$} \\ [1mm] 0,& \text{inače}, \end{cases}$$ (3.29)

$$v_n(x) = \begin{cases}\frac{x}{h} - ... ...je $x_{n-1}\leqslant{}x\leqslant{}x_n$} \\ [1mm] 0,& \text{inače}, \end{cases}$$ (3.30)

Imamo

$$v_i(x_j) = \delta_{i\,j},$$

pa je prema tome

$$u(x_j) = \sum_{i=0}^{n} c_i\,v_i(x_j) = \sum_{i=0}^{n} c_i\,\delta_{i\,j} = c_j.$$ (3.31)

Dakle, neodređeni koeficijenti $c_i$ su upravo vrijednosti približnog rješenja u odgovarajućim čvorovima. Budući da je uvjet na lijevom rubu homogen, imamo

$$u(0) = u(x_0) = c_0 = 0,$$

pa će, prema tome, rješenje biti oblika

$$u_n = \sum_{i=1}^{n} c_i\,v_i.$$

Budući da na desnom rubu imamo prirodni (Neumannov) uvjet, nema nikakvih daljnjih zahtjeva na oblik rješenja.

Iz formule (3.29) slijedi

$$v_i'(x) = \begin{cases} \frac{1}{h}... ... $x_i\leqslant{}x\leqslant{}x_{i+1}$} \\ [1mm] 0,& \text{inače}. \end{cases}$$

Odatle dobivamo

$$K_{i\,j} = \int_0^{\ell}\,p\,v'_i\,v... ...{i+1}} p(x)\,dx,& \text{ako je $j=i+1$} \\ [2mm] 0,& \text{inače}. \end{cases}$$ (3.32)

Dalje radimo isto kao kod Ritzove metode, da bismo na kraju dobili sustav linearnih algebarskih jednadžbi

$$\sum_{j=1}^n K_{i\,j}\,c_j = b_i,\qquad i=1,2,\ldots,n.$$

Elementi $K_{i\,j}$ za fiksni $i$ čine $i$-ti redak matrice krutosti $K.$ Formula (3.33) nam pokazuje da u jednom retku matrica $K$ ima najviše tri elementa različita od nule. To znatno pojednostavnjuje rješavanje sustava.

Primjer 3.21   Riješimo primjer 3.20 metodom konačnih elemenata.

Rješenje. Budući da je uvjet na desnom rubu homogen i prirodan, njega ne treba uzimati u obzir kod izbora koordinatnih funkcija. Uvjet na lijevom rubu nam kaže da treba isključiti funkciju $v_0.$ Dakle rješenje treba tražiti kao

$$u_{10} = {c_1}\,{v_1} + {c_2}\,{v_2} + {c_3}\,{v_3} + {c_4}\,{v_... ..._6}\,{v_6} + {c_7}\,{v_7} + {c_8}\,{v_8} + {c_9}\,{v_9} + {c_{10}}\,{v_{10}}.$$

Kad se izračunaju $K_{ij}$ i $b_i,$ dobije se sljedeći sustav jednadžbi

$$20.2667\,{c_1} - 10.2333\,{c_2} = 0.00883333,$$

$$-10.2333\,{c_1} + 20.8667\,{c_2} - 10.6333\,{c_3} = 0.0158333,$$

$$-10.6333\,{c_2} + 21.8667\,{c_3} - 11.2333\,{c_4} = 0.0208333,$$

$$-11.2333\,{c_3} + 23.2667\,{c_4} - 12.0333\,{c_5} = 0.0238333,$$

$$-12.0333\,{c_4} + 25.0667\,{c_5} - 13.0333\,{c_6} = 0.0248333,$$

$$-13.0333\,{c_5} + 27.2667\,{c_6} - 14.2333\,{c_7} = 0.0238333,$$

$$-14.2333\,{c_6} + 29.8667\,{c_7} - 15.6333\,{c_8} = 0.0208333,$$

$$-15.6333\,{c_7} + 32.8667\,{c_8} - 17.2333\,{c_9} = 0.0158333,$$

$$-17.2333\,{c_8} + 36.2667\,{c_9} - 19.0333\,{c_{10}} = 0.00883333,$$

$$-19.0333\,{c_9} + 19.0333\,{c_{10}} = 0.00158333$$

Iz ovog sustava se izračunaju koeficijenti

$${c_1} = 0.00883333, {c_2} = 0.0158333, {c_3} = 0.0208333, {c_4} = 0.0238333, {c_5} = 0.0248333,$$

$${c_6} = 0.0238333, {c_7} = 0.0208333, {c_8} = 0.0158333, {c_9} = 0.00883333, {c_{10}} = 0.00158333.$$

Budući da je rješenje između dva susjedna čvora linearna kombinacija polinoma prvog stupnja, rezultat će opet biti polinom prvog stupnja, tj. graf će biti dio pravca. Prema (3.32) rješenje je na rubu jednako odgovarajućim koeficijentima. Tako je za grafički prikaz rješenja dovoljno nanijeti točke $(x_i,c_i)$ kojima je dodana točka $(0,0)$ na lijevom rubu

$$(0,0), (0.1,0.00883333), (0.2,0.0158333), (0.3,0.0208333), (0.4,0.0238333),$$

$$(0.5,0.0248333), (0.6,0.0238333), (0.7,0.0208333),$$

$$(0.8,0.0158333), (0.9,0.00883333), (1.,0.00158333),$$

i zatim te točke spojiti ravnim linijama.