4 Πρόβλημα δύο σημείων: Πεπερασμένα Στοιχεία 4.3 Συνθήκες Neumann 4.5 Κυβικές splines

4.4 Υλοποίηση της μεθόδου Galerkin

Το βασικό πρόβλημα στην υλοποίηση της μεθόδου των πεπερασμένων στοιχείων είναι η κατασκευή του πίνακα (και, κατά δεύτερο λόγο του δεξιού μέλους) του συστήματος $Ac=F$ που προκύπτει από τη μέθοδο Galerkin. Δεδομένης της διαμέρισης $a=x_{0}<x_{1}<\cdots<x_{N}<x_{N+1}=b$ του $[a,b]$ , θα αναφερόμαστε σε κάθε διάστημα $I_{e}=[x_{e-1},x_{e}]$ , $e=1,\ldots,N+1$ ως πεπερασμένο στοιχείο. Το κάθε στοιχείο $I_{e}$ έχει δύο κόμβους, τους $x_{e-1}$ και $x_{e}$ , που αντιστοιχούν στους τοπικούς κομβικούς δείκτες $j=0$ και $j=1$ , αντίστοιχα. Ο δείκτης κάθε στοιχείου και οι τοπικοί κομβικοί δείκτες καθορίζουν τον λεγόμενο καθολικό κομβικό δείκτη $i$ κάθε κόμβου $x_{i}$ διαμέσου της απλής σχέσης

i=i(e,j)=e+j-1,\quad e=1,2,\ldots,N+1,\;j=0,1.

Οι συναρτήσεις $\{\phi_{i}\}_{i=0}^{N+1}$ , η βάση του χώρου πεπερασμένων στοιχείων της προηγούμενης παραγράφου, μπορούν να περιγραφούν με έναν εξαιρετικά απλό και αποτελεσματικό τρόπο, για την κατασκευή του πίνακα $A$ : συμβολίζουμε με $\phi_{j}^{e}$ , $j=0,1$ , τον περιορισμό στο διάστημα $[x_{e-1},x_{e}]$ των συναρτήσεων $\phi_{e-1},\phi_{e}$ , αντίστοιχα, και αναφερόμαστε σε αυτές ως τις τοπικές συναρτήσεις βάσης. Αν τώρα,

\Phi_{0}(y)=\begin{cases}1-y&\text{αν }\;0\leq y\leq 1\\ 0&\text{διαφορετικά}\end{cases},\qquad\Phi_{1}(y)=\begin{cases}y&\text{αν }\;0% \leq y\leq 1\\ 0&\text{διαφορετικά}\end{cases},

τότε

\phi_{j}^{e}(x)=\Phi_{j}\left(\frac{x-x_{e-1}}{h_{e}}\right),\quad j=0,1,

όπου $h_{e}=x_{e}-x_{e-1}$ . Ισοδύναμα, μπορούμε να αντιστοιχίσουμε σε κάθε διάστημα $I_{e}=[x_{e-1},x_{e}]$ τον ομοπαραλληλικό μετασχηματισμό $x=h_{e}y+x_{e-1}$ και να ορίσουμε τις τοπικές συναρτήσεις βάσης $\phi_{0}^{e},\phi_{1}^{e}$ ως

\phi_{0}^{e}(x)=\Phi_{0}(y),\;\phi_{1}^{e}(x)=\Phi_{1}(x),\;\text{όταν}\;x=h_{% e}y+x_{e-1}.

Συνεπώς, αν $v$ είναι μια συνάρτηση στον χώρο που παράγεται από τις συναρτήσεις $\{\phi_{i}\}_{i=0}^{N+1}$ , λαμβάνουμε την αναπαράσταση

v(x)=\sum_{i=0}^{N+1}\sum_{j=0}^{1}v_{i(e,j)}\phi_{j}^{e}(x),

όπου έχουμε θέσει $v_{i}=v(x_{i})$ . Ειδικώτερα, στο διάστημα $(x_{e-1},x_{e})$ έχουμε

v(x)=v_{e-1}\phi_{0}^{e}(x)+v_{e}\phi_{1}^{e}(x).

(4.42)

Ο συμβολισμός και οι παρατηρήσεις της προηγούμενης παραγράφου μας επιτρέπουν τώρα την εύκολη κατασκευή του πίνακα και του δεξιού μέλους των εξισώσεων της μεθόδου Galerkin. Στη γλώσσα των πεπερασμένων στοιχείων, η κατασκευή του πίνακα του συστήματος αναφέρεται ως συναρμολόγηση (assembly), ορολογία που δικαιολογείται από το γεγονός ότι

a(v,w)=\sum_{e}a_{e}(v,w),

όπου

a_{e}(v,w)=\int_{I_{e}}v^{\prime}w^{\prime}+qvw=\int_{x_{e-1}}^{x_{e}}[v^{% \prime}(x)w^{\prime}(x)+q(x)v(x)w(x)]\,dx,

για οποιεσδήποτε συναρτήσεις $v, w$ του χώρου πεπερασμένων στοιχείων. Χρησιμοποιώντας την αναπαράσταση (4.42) και τον ομοπαραλληλικό μετασχηματισμό $x\mapsto(x-x_{e-1})/h_{e}$ , έχουμε διαδοχικά

\begin{split}\displaystyle a_{e}(v,w)&\displaystyle=\int_{x_{e-1}}^{x_{e}}% \left(\sum_{j}v_{i(e,j)}\phi_{j}^{e}(x)\right)^{\prime}\left(\sum_{j}w_{i(e,j)% }\phi_{j}^{e}(x)\right)^{\prime}\,dx\\ &\displaystyle\quad+\int_{x_{e-1}}^{x_{e}}q(x)\left(\sum_{j}v_{i(e,j)}\phi_{j}% ^{e}(x)\right)\left(\sum_{j}w_{i(e,j)}\phi_{j}^{e}(x)\right)\,dx\\ &\displaystyle=\frac{1}{h_{e}}\int_{0}^{1}\left(\sum_{j}v_{i(e,j)}\Phi_{j}(y)% \right)^{\prime}\left(\sum_{j}w_{i(e,j)}\Phi_{j}(y)\right)^{\prime}\,dy\\ &\displaystyle\quad+h_{e}\int_{0}^{1}q(x_{e-1}+h_{e}y)\left(\sum_{j}v_{i(e,j)}% \Phi_{j}(y)\right)\left(\sum_{j}w_{i(e,j)}\Phi_{j}(y)\right)\,dy.\end{split}

Μπορούμε, τέλος, να γράψουμε την παραπάνω σχέση σε μορφή πίνακα επί διάνυσμα ως

a_{e}(v,w)=\frac{1}{h_{e}}\begin{bmatrix}v_{i(e,0)}&v_{i(e,1)}\end{bmatrix}S_{% e}\begin{bmatrix}w_{i(e,0)}\\ w_{i(e,1)}\end{bmatrix}+h_{e}\begin{bmatrix}v_{i(e,0)}&v_{i(e,1)}\end{bmatrix}% M_{e}\begin{bmatrix}w_{i(e,0)}\\ w_{i(e,1)}\end{bmatrix},

όπου $S_{e}$ είναι ο $2\times 2$ τοπικός πίνακας ακαμψίας

S_{e}=\begin{bmatrix}\displaystyle\int_{0}^{1}(\Phi_{0}^{\prime})^{2}&% \displaystyle\int_{0}^{1}\Phi_{0}^{\prime}\Phi_{1}^{\prime}\\ \displaystyle\int_{0}^{1}\Phi_{1}^{\prime}\Phi_{0}^{\prime}&\displaystyle\int_% {0}^{1}(\Phi_{1}^{\prime})^{2}\end{bmatrix}=\begin{bmatrix}1&-1\\ -1&1\end{bmatrix}

και $M_{e}$ είναι ο $2\times 2$ τοπικός πίνακας μάζας

M_{e}=\begin{bmatrix}\displaystyle\int_{0}^{1}\tilde{q}_{e}(\Phi_{0})^{2}&% \displaystyle\int_{0}^{1}\tilde{q}_{e}\Phi_{0}\Phi_{1}\\ \displaystyle\int_{0}^{1}\tilde{q}_{e}\Phi_{1}\Phi_{0}&\displaystyle\int_{0}^{% 1}\tilde{q}_{e}(\Phi_{1})^{2}\end{bmatrix}.

Εδώ, έχουμε θέσει $\tilde{q}_{e}(y)=q(x_{e-1}+h_{e}y)$ . Η προσέγγιση των στοιχείων του πίνακα $M_{e}$ μπορεί να γίνει στην πράξη με τη χρήση ενός κανόνα αριθμητικής ολοκλήρωσης. Στην περίπτωση που εξετάζουμε στη συγκεκριμένη παράγραφο, μπορεί να αποδειχθεί ότι η χρήση του κανόνα του τραπεζίου είναι ικανή για να εξασφαλίσει το αποτέλεσμα του Θεωρήματος 4.2.