Kapitel 15 Die Methode der finiten Elemente Die Darstellung der Methode der FE stützt sich im wesentlichen auf die Kapitel II und III des Buches über FE von D. Braess. 15.1 Erläuterung der Grundideen anhand des MP’s Bevor wir die Methode der finiten Elemente systematisch einführen, wollen wir vorab schon die Grundideen dieses Ansatzes anhand des Modellproblems −∆u = f Ω ⊂ IR2 , u = 0 (∂Ω) (15.1) erläutern. 1.Grundidee: Die erste Grundidee besteht darin, der Randwertaufgabe eine Variationsaufgabe zuzuordnen und zwar in dem Sinne, dass die (noch zu spezifizierende) Lösung des Variationsproblems gleichzeitig auch Lösung der Randwertaufgabe ist. Für das Modellproblem betrachten wir hierzu das Variationsfunktional Z Z J [u] := u2x + u2y dΩ − 2 uf dΩ Ω Ω auf dem Raum D0 := u ∈ C2 (Ω) ∩ C Ω : u = 0 auf ∂Ω . Gesucht ist nun das Minimum u ∈ D0 , so dass J [u∗ ] ≤ J [u] für alle u ∈ D0 . Satz 15.1 (RWA ↔ Variationsproblem). Falls u∗ ∈ D0 obige Variationsaufgabe löst, dann ist u∗ auch eine klassische Lösung der Randwertaufgabe (15.1). 102 Die Methode der finiten Elemente 103 Beweis: u∗ löse die Variationsaufgabe. Betrachte u∗ +tv mit t ∈ IR, v ∈ D0 . Dann gilt: Z Z ∗ ∗ 2 2 2 J [u ] ≤ J [u + tv] = t vx + vy dΩ + t 2 u∗x vx + u∗y vy dΩ Ω Ω Z Z Z ∗2 ∗2 −2 f vdΩ + ux + uy dΩ − 2 f u∗ dΩ Ω Ω Ω Dabei sei v fest und t beliebig. Dann nimmt die Funktion Φv (t) = J [u∗ + tv] für t = 0 ihr Minimum an. Damit ergibt sich Z Z ∂Φv (0) ∗ ∗ ux vx + uy vy dΩ − f vdΩ = 0 = . (15.2) ∂t |Ω {z } | Ω {z } =:a(u∗ ,v) =:hf,vi Partielle Integration liefert dann Z Z Z ∗ ∗ − v uxx + uyy dΩ + v . . . − f vdΩ = 0 Ω ∂Ω Ω | {z } =0 Z ⇒ v (−∆u − f ) dΩ = 0 ∀v ∈ D0 . Ω Wäre jetzt −∆u − f 6= 0, also z.B. > 0 an einer Stelle (x∗ , y ∗ ), dann wäre −∆u − f > 0 auch in einer Umgebung von (x∗ , y ∗ ) aufgrund von Stetigkeitsargumenten. Man kann dann v so wählen, dass der Integrand in dieser Umgebung größer als Null und außerhalb der Umgebung gleich Null ist. Dies führt zu einem Widerspruch. Folglich ist −∆u(x, y) = f (x, y) ((x, y) ∈ Ω) . Die Randbedingungen werden offensichtlich wegen u ∈ D0 automatisch erfüllt. 2.Grundidee: Im zweiten Schritt drückt man die Variationsaufgabe mithilfe einer Bilinearform aus, siehe (15.2): Falls a (u∗ , v) = hf, vi ∀v ∈ D0 , (15.3) dann löst u∗ auch die RWA (15.1). Diese Formulierung der Aufgabe mithilfe einer geeigneten Bilinearform bezeichnet man als schwache Formulierung der RWA. Bezüglich der Formulierung (15.3) reicht offensichtlich, dass Die Methode der finiten Elemente 104 u∗ ∈ C1 (Ω) im Gegensatz zu einer klassischen Lösung u∗ ∈ C2 (Ω) der RWA (15.1). Eine derartige Abschwächung des Lösungsbegriffs wird systematisch in Abschnitt 15.3 behandelt. 3.Grundidee: Als nächstes ersetzt man den unendlichdimensionalen Raum D0 durch einen endlichdimensionalen Raum D0M = span {ψ1 , . . . , ψM }. Dies bedeutet also nichts anderes als eine P Diskretisierung der Variationsaufgabe, indem wir u∗ ∈ D0 durch u∗M = M i=1 ci ψi approximieren. Hierbei werden die Koeffizienten ci mithilfe der M Gleichungen a (u∗M , ψk ) = hf, ψk i ⇔ M X ci a (ψi , ψk ) = hf, ψk i (k = 1, . . . , M ) i=1 bestimmt. Dies liefert ein lineares Gleichungssystem c1 .. A . = b cM zur Berechnung der ci und zwar mit A = (aik )i,k=1,...,M = (a (ψi , ψk ))i,k=1,...,M und b = (hf, ψk i)k=1,...,M . Bemerkung 15.1 (Galerkin-, Ritz-Verfahren). Die unter 3. beschriebene Vorgehensweise nennt man Galerkin-Verfahren. Bei inhomogenen Randbedingungen wählt man zusätzlich eine Funktion ψ0 , die die inhomogenen Randbedingungen erfüllt: u∗M = ψ0 = M X ci ψi . i=0 Dies führt zu einer Modifikation der rechten Seite b. Falls die Bilinearform a spd ist, ergibt sich das sogenannte Ritz-Verfahren. Die Verwendung dieser Begriffe ist jedoch nicht einheitlich in der Literatur. Häufig spricht man auch allgemein vom Ritz-Galerkin-Verfahren unabhängig davon, ob a spd ist oder nicht. I Klassische Ansatzfunktionen zur Konstruktion von D0M : 1. Polynome 2. Produkte von Orthogonalpolynomen Die Methode der finiten Elemente 105 3. Produkte trigonometrischer Funktionen 4. Produkte von Spline-Funktionen 4.Grundidee. Bei der Methode der finiten Elemente wählt man die Basis/Ansatzfunktionen ψk derart, dass das resultierendes GS dünnbesetzt ist. Hierzu • unterteilt man Ω näherungsweise in Dreiecke und/oder Vierecke und • verwendet Ansatzfunktionen, die nur lokal (d.h. nur auf wenigen Elementen der Zerlegung) ungleich Null sind, aber global (d.h. über Elementgrenzen hinweg) wohldefinierte Glattheitsbedingungenen erfüllen. Insbesondere seien die ψk durch Polynome niedrigen Grades darstellbar. Beispiel 15.1 (Lineare FE, reguläre Triangulierung, MP). Als spezielles Beispiel betrachten wir das Modellproblem auf dem Einheitsquadrat mit einer gleichmäßigen regulären Triangulierung und linearen finiten Elementen: −∆u = f Ω = (0, 1)2 u=0 (∂Ω) D0M := v ∈ C Ω , v ist in jedem Element linear, und v = 0 auf ∂Ω . In jedem Dreieck hat v ∈ D0M die Form v (x, y) = a + bx + cy und ist durch die Werte in den 3 Ecken der Dreiecke eindeutig bestimmt. v ist somit global durch die Werte an den M Gitterpunkten festgelegt und es ergibt sich somit: M dim D0 = M . Eine geeignete Basis (d.h. eine Basis, die möglichst viele Nulleinträge in der resultierenden Matrix A nachsichzieht) stellt die sogenannte nodale Basis dar, für die gilt: ψk (x` , y` ) = δk` für k = 1, . . . , M . Die Methode der finiten Elemente 11111 00000 00000 11111 00000 11111 00000 11111 N NW 2. W 3. Z 1. O 4. 6. 5. SO S 106 Die Ansatzfunktion ψk ist offensichtlich nur im schraffierten Bereich der nebenstehenden Zeichnung ungleich Null. Hierbei sei Z = Zk = (xk , yk ) = h(ik , jk ). Derartige Ansatzfunktionen werden auch als Hutfunktionen bezeichnet. h Wir wollen nun die Darstellung von ψk in jedem der sechs schraffierten Dreiecke berechnen. Hierzu müssen wir die zugehörigen Koeffizienten a, b, c der allgemeinen Darstellung v(x, y) = a + bx + cy der linearen finiten Elemente in jedem der sechs Dreiecke bestimmen. Für die Berechnung der Matrixeinträge a(ψi , ψk ) benötigen wir die entsprechenden Gradienten von ψk , die wir deshalb direkt mit angeben. 1. Dreieck: a + bhik + chjk = 1 , a + bh (ik + 1) + chjk = 0 , a + bhik + ch (jk + 1) = 0 . x y ψk (x, y) = 1 + ik + jk − − . h h 1 1 T ⇒∇ψk (x, y) = − , − . h h ⇒ 2. Dreieck: a + bhik + chjk = 1 , a + bhik + ch (jk + 1) = 0 , a + bh (ik − 1) + ch (jk + 1) = 0 . ⇒ 3. Dreieck: . . . y . ψk (x, y) = 1 + jk − h 1 T ⇒∇ψk (x, y) = 0, − . h x ψk (x, y) = 1 − ik + h ⇒ ∇ψk = 1 ,0 h T . Die Methode der finiten Elemente 4. Dreieck: . . . 107 x y ψk (x, y) = 1 − ik − jk + + h h 6. Dreieck: . . . x ψk (x, y) = 1 + ik − h ∇ψk = y ψk (x, y) = 1 − jk + h 5. Dreieck: . . . ⇒ ⇒ ∇ψk = ⇒ ∇ψk = 1 0, h 1 1 , h h T T . T 1 − ,0 . h Schließlich bestimmen wir die Einträge von AM mit Z a (ψi , ψk ) = ∇ψi ∇ψk dx . Ω Z Z ∇ψz ∇ψz dx = 2 ∇ψz ∇ψz dx Ω 1.+2.+3. Z Z Z 2 1 1 =2 dx + dx + dx 2 2 2 1. h 2. h 3. h 2 h2 1 h2 1 h2 =2 2 + 2 + 2 =4 h 2 h 2 h 2 Z Z a(ψZ , ψN ) = ∇ψz ∇ψN dx = ∇ψz ∇ψN dx . a(ψZ , ψZ ) = Ω 1.+2. Z 1 1 T 1 T 1 T 1 1 T = − ,− 0, dx + 0, − , dx h h h h h h 1. 2. Z Z 1 1 h2 1 h2 1 − 2 dx = − 2 − 2 = − 2 dx + h h h 2 h 2 2. 1. = −1 = a (ψZ , ψS ) = a (ψZ , ψW ) = a (ψZ , ψO ) Z Z a(ψZ , ψN W ) = ∇ψZ ∇ψN W dx = ∇ψz ∇ψN w dx . Z Ω Z 2.+3. T 1 1 = 0, − − ,0 h h 2. = 0 = a (ψZ , ψSO ) . T Z dx + 3. 1 ,0 h T 1 0, h T dx Alle anderen Einträge mit nicht benachbarten ψk sind offensichtlich ebenfalls Null. In Sternschreibweise ergibt sich somit dasselbe lineare Gleichungssystem, wie wir es schon bei den FD und den FV kennen gelernt haben: −1 AM = b −1 4 −1 . −1 h . Die Methode der finiten Elemente 108 Dieser Zusammenhang zwischen den verschiedenen Diskretisierungsansätzen ist jedoch lediglich für dieses einfache Modellproblem charakteristisch und läßt sich auf allgemeinere Fälle nicht übertragen. J 15.2 Allgemeinere Aufgaben, weitere wichtige Beispiele 15.2.1 Inhomogene Dirichlet-Randbedingungen Häufig lassen sich RWA’en mit inhomogenen Dirichlet-RB’en auf RW’en mit homogenen Dirichlet-RB’en zurückführen. Man Betrachte den elliptischen Differentialoperator 2. Ordnung Lu = − N X ∂i (aik ∂k u) + a0 u i,k=1 mit a0 (x) ≥ 0 für x ∈ Ω. Ziel ist es das Dirichletproblem Lu = f u=g (x ∈ Ω) (x ∈ ∂Ω) auf ein solches mit homogenen RB’en zu transformieren. Hierzu sei eine zulässige Funktion u0 bekannt mit u0 = g auf ∂Ω. Daraus ergibt sich die modifizierte RWA Lw = f1 (x ∈ Ω) w=0 (x ∈ ∂Ω) mit w = u − u0 , f1 = f − Lu. Insbesondere für theoretische Zwecke wird daher oft angenommen, dass die Aufgabe mit homogenen Randbedingungen gegeben ist. 15.2.2 Poisson-Gleichung mit Neumann-Randbedingungen Betrachte die RWA −∆u = f un = 0 (Ω) (∂Ω) . Wir setzen voraus, dass ∂Ω hinreichend glatt ist (also z.B. ∈ C∞ ). J[u] sei b = C2 (Ω) ∩ wie in Abschnitt 15.1 definiert, wobei anstelle von D0 jetzt D Die Methode der finiten Elemente 109 b keine Randbedingungen C1 (Ω) betrachtet wird. Man beachte, dass in D auftauchen. Dafür wird jetzt einmal stetige Differenzierbarkeit bis zum Rand gefordert. Der Gedankengang ist wieder ähnlich wie in Abschnitt 15.1, mit Φv (t) = J[u∗ + tv] mit v ∈ C2 (Ω) . Damit erhält man dann Z Z ∗ ∗ f v d(x, y) = 0 , (ux vx + uy vy ) d(x, y) − Ω Ω und hieraus folgt durch partielle Integration Z Z ∗ (−∆u − f )v d(x, y) + u∗n v do(x, y) = 0 Ω für alle v ∈ C2 (Ω) , ∂Ω 0 (Ω). Mit dem gleichen Schluß wie oben folgt insbesondere für alle v ∈ C∞ −∆u = f (Ω) . Ein ähnlicher Schluß für das Randintegral (die Annahme un 6= 0 führt man durch geeignete Wahl von η zu einem Widerspruch) liefert un = 0 (∂Ω) . b erfüllt zu Offenbar brauchen hier die Randbedingungen nicht im Raum D sein; sie ergeben sich automatisch“. Es sind natürliche“ Randbedingungen. ” ” Die allgemeine Definition von wesentlichen und natürlichen Randbedingungen erfolgt später. 15.2.3 Minimalflächengleichung Betrachte die nichtlineare RWA (1 + (ux )2 )uxx − 2ux uy uxy + (1 + (uy )2 )uyy = 0 u=g (Ω) (∂Ω) und das zugehörige Variationsfunktional Z q J[u] = 1 + (ux )2 + (uy )2 d(x, y) . Ω Dieses Integral beschreibt den Flächeninhalt der durch u beschriebenen Fläche; er wird minimal für die Lösung der Minimalflächengleichung. Physikalisch beschreibt die RWA eine elastische Membran, die in eine Raumkurve Die Methode der finiten Elemente 110 g eingespannt wird. Sie stellt sich so ein, dass die eingespannte Fläche minimal wird. Bei einer eingespannten elastischen Membran unter einer Belastung f ist das Funktional Z q Z 2 2 J[u] = ( 1 + (ux ) + (uy ) − 1) d(x, y) − f u d(x, y) . Ω Ω Dabei beschreibt das zweite Integral die Arbeit der äußeren Kräfte, die −1“ ” im ersten Integral repräsentiert die Flächenänderungsenergie. Der Zusammenhang zur Poissongleichung wird deutlich, wenn man (für kleine Auslenkungen ux , uy ) eine Taylorentwicklung durchführt. In erster Näherung ist q 1 + (ux )2 + (uy )2 − 1 = ˙ (ux )2 + (uy )2 , d.h. die Poissongleichung stellt eine globale Linearisierung dar. 15.2.4 Biharmonische Gleichung Die DGL ∆∆u = uxxxx + 2uxxyy + uyyyy = f beschreibt die Auslenkung einer elastischen dünnen Platte unter einer Belastung f (bei kleinen Auslenkungen). Die Randbedingungen u = g1 (∂Ω) un = g2 (∂Ω) (15.4) beschreiben eine eingespannte Platte (u beschreibt die Lage, un die Richtung der Platte). Man benötigt hier zwei Randbedingungen, da die DGL von vierter Ordnung ist. Die Randbedingungen u = g1 (∂Ω) ∆u = 0 (∂Ω) (15.5) beschreiben eine frei aufliegende Platte (∆u beschreibt das Drehmoment). Das zugehörige Variationsfunktional lautet Z J[u] = ((∆u)2 − 2f u) d(x, y) . Ω Die Methode der finiten Elemente 111 Dabei sind im Falle von (15.4) beide Randbedingungen wesentlich und müssen im Raum Dg berücksichtigt werden. Im Falle von (15.5) ist nur die Bedingung u = g1 eine wesentliche Randbedingung, −∆u = 0 ist eine natürliche Randbedingung (ergibt sich automatisch). Allgemein haben bei einer DGL 2m-ter Ordnung die Randbedingungen eine Ordnung ≤ 2m − 1. Dabei sind die Ableitungen mit einer Ordnung ≤ m − 1 wesentliche Randbedingungen. 15.3 Variationsformulierung elliptischer RWA’en Eine systematische Verallgemeinerung der ersten beiden Grundideen aus Abschnitt 15.1 führt auf den Begriff der schwachen Lösung (Lösung der Variationsaufgabe bzw. der Formulierung mit Bilinearform). Hierbei wird eine Abschwächung der Glattheitsvoraussetzungen an die Lösung vorgenommen (vgl. Beispiel 15.1). Insgesamt ergeben sich somit deutlich schwächere Voraussetzungen für die Existenz von Lösungen als bei der Methode der FD (bzw. FV). 15.3.1 Sobolevräume Grundlegend für den Begriff der schwachen Lösung sind die sogenannten Sobolevräume. Sobolevräume basieren auf dem L2 (Ω), der mit dem Skalarprodukt Z hu, viL2 =: hu, vi0 = uv dΩ Ω und der zugehörigen Norm kuk0 = q hu, vi0 zu einem Hilbertraum wird. Sobolevräume fassen Funktionen mit sogenannten schwachen Ableitungen zusammen. Der Begriff der schwachen Ableitung wird durch die partielle Integration motiviert: Für u, φ ∈ C1 [a, b] gilt Z b Z b 0 uφ dx = − u0 φ dx + uφ|ba . (15.6) a a 0 [a, b] gilt, dass Falls für u, v ∈ L2 [a, b] für alle φ ∈ C∞ Z b Z b 0 uφ dx = − vφ dx a a
© Copyright 2024 ExpyDoc