<?page no="1"?> Jürgen Ulm Numerische Lösung gewöhnlicher und partieller Differenzialgleichungen <?page no="3"?> Numerische Lösung gewöhnlicher und partieller Differenzialgleichungen Finite-Elemente-Methode (FEM) - Finite-Differenzen-Methode (FDM) - Aufgaben mit Lösungen Prof. Dr.-Ing. Jürgen Ulm Mit 46 Bildern und 20 Tabellen Kontakt & Studium Band 707 Herausgeber: Prof. Dr.-Ing. Dr. h.c. Wilfried J. Bartz Dipl.-Ing. Hans-Joachim Mesenholl Dipl.-Ing. Elmar Wippler TAE <?page no="4"?> Bei der Erstellung des Buches wurde mit großer Sorgfalt vorgegangen; trotzdem lassen sich Fehler nie vollständig ausschließen. Verlag und Autoren können für fehlerhafte Angaben und deren Folgen weder eine juristische Verantwortung noch irgendeine Haftung übernehmen. Für Verbesserungsvorschläge und Hinweise auf Fehler sind Verlag und Autoren dankbar. © 2018 by expert verlag, Wankelstr. 13, D -71272 Renningen Tel.: + 49 (0) 71 59 - 92 65 - 0, Fax: + 49 (0) 71 59 - 92 65 - 20 E-Mail: expert@expertverlag.de, Internet: www.expertverlag.de Alle Rechte vorbehalten Printed in Germany Das Werk einschließlich aller seiner Teile ist urheberrechtlich geschützt. Jede Verwertung außerhalb der engen Grenzen des Urheberrechtsgesetzes ist ohne Zustimmung des Verlags unzulässig und strafbar. Dies gilt insbesondere für Vervielfältigungen, Übersetzungen, Mikroverfilmungen und die Einspeicherung und Verarbeitung in elektronischen Systemen. ISBN 978-3-8169-3420-2 Bibliografische Information Der Deutschen Bibliothek Die Deutsche Bibliothek verzeichnet diese Publikation in der Deutschen Nationalbibliografie; detaillierte bibliografische Daten sind im Internet über http: / / www.dnb.de abrufbar. Bibliographic Information published by Die Deutsche Bibliothek Die Deutsche Bibliothek lists this publication in the Deutsche Nationalbibliografie; detailed bibliographic data are available on the internet at http: / / www.dnb.de <?page no="5"?> Vorwort Wenn Sie wissen m¨ochten, wie Sie mit einem Taschenrechner (Abb. 1) gew¨ohnliche oder partielle Differenzialgleichungen l¨osen, Temperaturverl¨aufe oder elektrische und magnetische Felder berechnen k¨onnen, dann treffen wir uns auf den kommenden Seiten wieder. Abbildung 1: Ihr Taschenrechner kann mehr als Sie ihm zutrauen! Es sei noch ein Hinweis zur erweiterten Nutzung des Skriptes gestattet: Neue ¨ Ubungsaufgaben lassen sich durch einfaches Ab¨andern der gestellten und bereits gel¨osten Originalaufgabe generieren. Die Ab¨anderung der Originalaufgabe soll in der Weise vorgenommen werden, dass deren L¨osung bereits im Voraus bekannt ist. Damit besteht die M¨oglichkeit, die Ergebnisse zu vergleichen und die Einarbeitung weiter zu vertiefen. i <?page no="6"?> Denn ”Ungl¨uck kommt von mangelhaften Berechnungen“ (aus Bertolt Brecht, ”Das Leben des Galilei.“) Mein Dank gilt meinen Institutsmitarbeitern Herrn Wilhelm Feucht, Herrn Reinhardt Erli, Herrn Oliver Vogel und Herrn Dimitri Delkov f¨ ur die Bereitstellung der Bilder, Simulationsergebnisse und Grafiken. F¨ ur die kritische Durchsicht dieses Skriptes und fruchtbaren Diskussionen gilt nochmal Herrn Delkov mein Dank. Last but not least ein herzliches Dankesch¨on an meine Frau Monika und Tochter Natascha f¨ ur die unendliche Geduld mit meiner Person und meiner fortdauernden Ausrede ”hab gerade keine Zeit...“. Mit freundlichen Gr¨ ußen der Autor im Herbst 2017 ii <?page no="7"?> Inhaltsverzeichnis 1 Erforderliche mathematische Grundlagen 1 1.1 Matrizen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1.1 Rechenoperationen mit Matrizen . . . . . . . . . . . . . . . . . 1 1.1.2 Determinante . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.1.3 Inverse Matrix . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.2 Differenzialgleichungen . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.2.1 Definitionen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2.2 Partielle Differenzialgleichungen . . . . . . . . . . . . . . . . . . 4 1.2.3 Partielle Integration . . . . . . . . . . . . . . . . . . . . . . . . 5 1.2.4 Klassifikation von Differenzialgleichungen . . . . . . . . . . . . . 6 1.2.5 Anfangswertaufgabe . . . . . . . . . . . . . . . . . . . . . . . . 7 1.2.6 Randwertaufgabe . . . . . . . . . . . . . . . . . . . . . . . . . . 7 1.2.7 Inneres Produkt . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 1.2.8 Starke Form/ Formulierung einer Differenzialgleichung . . . . . . 11 1.2.9 Schwache Form/ Formulierung einer Differenzialgleichung . . . . 11 1.3 Vektoroperatoren . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.3.1 Nabla- und Laplaceoperator . . . . . . . . . . . . . . . . . . . . 12 1.3.2 Vektoroperator Gradient . . . . . . . . . . . . . . . . . . . . . . 12 1.3.3 Vektoroperator Divergenz . . . . . . . . . . . . . . . . . . . . . 13 1.3.4 Vektoroperator Rotation . . . . . . . . . . . . . . . . . . . . . . 13 1.3.5 Gegen¨ uberstellung der Vektoroperatoren . . . . . . . . . . . . . 13 1.3.6 N¨ utzliche Normen . . . . . . . . . . . . . . . . . . . . . . . . . . 14 2 Differenzialgleichungen und Finite Elemente 16 2.1 Physik-Beispiele f¨ ur Differenzialgleichungen erster Ordnung . . . . . . . 16 iii <?page no="8"?> 2.2 Physik-Beispiele f¨ ur Differenzialgleichungen zweiter Ordnung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 2.3 Finite Elemente . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3 Von der Momentenmethode zur Galerkin-Methode 21 3.1 Grundprinzip der Momentenmethode . . . . . . . . . . . . . . . . . . . 21 3.2 Anmerkungen zur Momentenmethode . . . . . . . . . . . . . . . . . . . 23 3.2.1 Matrix [l mn ] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 3.2.2 Wahl der Basis- und Wichtungsfunktionen f n und w n . . . . . . 23 3.3 Galerkins Idee . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 3.4 Traditionelle Galerkin-Methode . . . . . . . . . . . . . . . . . . . . . . 24 3.5 Galerkin-FEM-Methode . . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.6 Vorgehen zur L¨osung mit der Galerkin-Methode . . . . . . . . . . . . . 27 4 L¨osung der Gleichung dy dx − y = 0 mit der Galerkin-Methode 29 4.1 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . . . . . 29 4.2 Formulierung der schwachen Form mit Basis- und Wichtungsfunktion . 30 4.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . . . . . 31 4.4 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . . . . . 31 5 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode 33 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion . . . . . . . . . . . . 33 5.1.1 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 34 5.1.2 Diskretisierung des zu l¨osenden Gebiets Ω . . . . . . . . . . . . 35 5.1.3 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 36 5.1.4 Formulierung der schwache Form mit Dreiecksfunktionen φ(x) . 37 5.1.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 38 5.1.6 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 41 5.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion . . . . . . . . . 44 5.2.1 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 44 5.2.2 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 45 5.2.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 45 5.2.4 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 45 6 L¨osung der Gleichung u(x) = − 1 2 x 2 + 2x mit der Galerkin-Methode 47 6.1 L¨osung mit linearer Basis- und Wichtungsfunktion . . . . . . . . . . . . 48 iv <?page no="9"?> 6.1.1 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 49 6.1.2 Diskretisierung des zu l¨osenden Gebiets Ω . . . . . . . . . . . . 49 6.1.3 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 49 6.1.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) 50 6.1.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 50 6.1.6 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 50 6.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion . . . . . . . . . 51 6.2.1 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 52 6.2.2 Formulierung der schwachen Form mit Basis- und Wichtungsfunktion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 6.2.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 53 6.2.4 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 54 7 L¨osung physik. Bsp. DGL 1’ter Ordnung mit Galerkin-Methode 56 7.1 Durchflutungsgesetz gel¨ost mit nichtlinearer Basis- und Wichtungsfunktion 56 7.1.1 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 57 7.1.2 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 59 7.1.3 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 59 7.2 Gegen¨ uberstellung FEMmit Galerkin-Ergebnis . . . . . . . . . . . . . 60 8 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode 63 8.1 Elektrostatische Feldberechnung . . . . . . . . . . . . . . . . . . . . . . 63 8.1.1 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 63 8.1.2 Diskretisierung des zu l¨osenden Gebiets Ω . . . . . . . . . . . . 64 8.1.3 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 64 8.1.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) 64 8.1.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 66 8.1.6 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 68 8.2 Ortsabh¨angige Temperaturberechnung . . . . . . . . . . . . . . . . . . 70 8.2.1 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 71 8.2.2 Diskretisierung des zu l¨osenden Gebiets Ω . . . . . . . . . . . . 72 8.2.3 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 72 8.2.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) 72 8.2.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 73 8.2.6 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 75 v <?page no="10"?> 8.2.7 Diffusionsvorgang vollendet . . . . . . . . . . . . . . . . . . . . 78 8.3 Ortsabh¨angige Magnetfeldberechnung . . . . . . . . . . . . . . . . . . . 79 8.3.1 Schwache Formulierung der Differenzialgleichung . . . . . . . . . 80 8.3.2 Diskretisierung des zu l¨osenden Gebiets Ω . . . . . . . . . . . . 80 8.3.3 Wahl der Basis- und Wichtungsfunktion . . . . . . . . . . . . . 81 8.3.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) 81 8.3.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung . 81 8.3.6 L¨osung des linearen Gleichungssystems . . . . . . . . . . . . . . 83 9 Einf¨ uhrung in die Finite-Differenzen-Methode 88 9.1 Numerische Notation der linearen Felddiffusionsgleichung . . . . . . . . 88 9.2 L¨osung mit impliziter Methode nach Crank-Nicolson . . . . . . . . . . 89 9.2.1 ¨ Uberf¨ uhrung der Diffusionsgleichung in eine Matrizengleichung . 89 9.2.2 L¨osung der Matrizengleichung . . . . . . . . . . . . . . . . . . . 91 9.2.3 Anwendungsbeispiel . . . . . . . . . . . . . . . . . . . . . . . . 94 9.3 L¨osung mit expliziter Methode . . . . . . . . . . . . . . . . . . . . . . . 96 9.3.1 ¨ Uberf¨ uhrung der Diffusionsgleichung in eine Matrizengleichung . 97 9.3.2 L¨osung der Matrizengleichung . . . . . . . . . . . . . . . . . . . 98 9.3.3 Anwendungsbeispiel . . . . . . . . . . . . . . . . . . . . . . . . 99 10 Anwendungen der FEM zur Produktentwicklung 104 10.1 Analyse eines Proportionalmagnetaktors . . . . . . . . . . . . . . . . . 104 10.1.1 Preprocessing . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105 10.1.2 Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 10.1.3 Postprocessing . . . . . . . . . . . . . . . . . . . . . . . . . . . 107 10.2 Synthese eines planaren Asynchron-Scheibenl¨aufermotors . . . . . . . . 108 10.2.1 Preprocessing . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 10.2.2 Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 10.2.3 Postprocessing . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 10.2.4 Musterbau des planaren Asynchronmotors . . . . . . . . . . . . 109 11 Anwendung der FEM zur Produktoptimierung 111 A 114 A.1 MATLAB-Code - W¨armediffusionsskript . . . . . . . . . . . . . . . . . 114 A.2 MATLAB-Code - Magnetfelddiffusionsskript . . . . . . . . . . . . . . . 118 vi <?page no="11"?> A.3 MATLAB-Ergebnisse vs. COMSOL Multiphysics-Ergebnisse . . . . . . 124 Literaturverzeichnis 127 vii <?page no="12"?> Symbole und Abk¨ urzungen Symbol Bedeutung Einheit A Koeffizient, Matrix, Fl¨ache B Koeffizient, Matrix B Magnetische Flussdichte V s/ m 2 B Vektor der magnetischen Flussdichte T B h Interpolations-, Ansatzfunktion C Koeffizient, Matrix C spezifische W¨armekapazit¨at J/ (kgK) D Koeffizient, Matrix E Koeffizient, Matrix E Vektor der elektrischen Feldst¨arke V / m F Koeffizient, Funktion, Kraft G Koeffizient H Vektor des magnetischen Feldes A/ m H Φ Interpolations-, Ansatzfunktion J Vektor der elektrischen Stromdichte A/ m 2 K Konstante L Induktivit¨at V s/ A L Linearoperator M Matrix N Anzahl Knoten P Leistung W R Residuum R Radius m R Widerstand Ω S P Scheitelpunkt viii <?page no="13"?> ix Symbol Bedeutung Einheit a Koeffizient a 0 Beschleunigung m/ s 2 c Konstante f Hilfsvariable, Funktion f n Entwicklungs-, Basisfunktion g Hilfsvariable, Funktion, Matrix h Elementl¨ange m i Laufvariable j Laufvariable k Konstante l Matrix m Laufvariable n Normale, Anzahl Teilintervalle r Radius m s Konstante t Zeit s u Funktion u h Interpolations-, Ansatzfunktion v Funktion w Gewichts-, Wichtungs-, Test-, Formfunktion x Koordinate, Weg m y Koordinate, Weg, Funktion m z Koordinate m Symbol Bedeutung Einheit Γ Rand des FEM-Gebietes Γ l linker Rand des FEM-Gebietes Γ r rechter Rand des FEM-Gebietes Δ Delta, differenziell Ω Gebiet Ω n Element, Teilgebiete <?page no="14"?> x Symbol Bedeutung Einheit α Koeffizient β Koeffizient γ Randwert δ Abklingkonstante ε Permittivit¨at As/ (V m) ε Permittivit¨at des Vakuums As/ (V m) υ Temperatur ◦ C κ spezifische elektrische Leitf¨ahigkeit m/ (Ωmm 2 ) λ W¨armeleitf¨ahigkeit W/ (mK) μ Permeabilit¨at V s/ (Am) μ 0 Permeabilit¨at des Vakuums [4π10 −7 V s/ (Am)] V s/ (Am) ρ Dichte kg/ m 3 υ h Ansatz-, Testfunktion ϕ Potenzial eines L¨osungsraums V ϕ e Potenzial eines Elements V ϕ ei Potenzial eines Elementknotens V ϕ h Interpolations-, Ansatzfunktion φ Entwicklungs-, Basis-, Dreiecksfunktion ω Winkelgeschwindigkeit 1/ s ∇ Nabla-Operator <?page no="15"?> Kapitel 1 Erforderliche mathematische Grundlagen Die zur numerischen L¨osung von Differenzialgleichungen erforderlichen Grundlagen sind in diesem Kapitel zusammengestellt worden. Diese beinhalten im Wesentlichen Matrizen, Definitionen und Klassifikationen von Differenzialgleichungen sowie Anfangs- und Randwertaufgaben und Vektoroperatoren. 1.1 Matrizen Hier werden die wesentlichen Matrizenoperationen vorgestellt. Diese beinhalten die erforderlichen Matrizen-Rechenregeln, die Determinantenberechnung sowie die Invertierung einer Matrix. 1.1.1 Rechenoperationen mit Matrizen In Tab. 1.1 werden die wichtigsten algebraischen Axiome zusammengefasst. Tabelle 1.1: Zusammenfassung der wichtigsten Rechenregeln Assoziativgesetz A (BC) = (AB) C Distributivgesetz A (B+C) = AB + AC (A+B) C = AC + BC Transponieren (AB) T = B T A T 1 <?page no="16"?> 2 Erforderliche mathematische Grundlagen Man beachte, dass die Matrizenmultiplikation nicht kommutativ ist. 1.1.2 Determinante Unter der Determinante einer quadratischen 3x3-Matrix A = (a ik ) versteht man die reelle Zahl det A = ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ a 11 a 12 a 13 a 21 a 22 a 23 a 31 a 32 a 33 ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ . Die 3-reihige Determinante wird nach der Regel von Sarrus detA = a 11 a 22 a 33 + a 12 a 23 a 31 + a 13 a 21 a 32 − a 13 a 22 a 31 − a 11 a 23 a 32 − a 12 a 21 a 33 berechnet. 1.1.3 Inverse Matrix Die Berechnung der inversen Matrix A −1 A −1 = 1 det A ⎛ ⎜ ⎜ ⎜ ⎝ A 11 A 12 ... A 1 m A 21 A 22 ... ... ... ... ... ... A n 1 A n 2 ... A nm ⎞ ⎟ ⎟ ⎟ ⎠ erfolgt unter Verwendung von Unterdeterminanten mit dem Gauß’schen Algorithmus. Weiterhin gilt A · A −1 = A −1 · A = E. 1.2 Differenzialgleichungen Viele Vorg¨ange in Naturwissenschaft und Technik werden mittels Differenzialgleichungen beschrieben. Um den Zugang zu den Differenzialgleichungen zu erleichtern, werden diese hier vorgestellt. Nach anf¨anglichen Begriffsdefinitionen wird eine Klassifizierung <?page no="17"?> 1.2 Differenzialgleichungen 3 von Differenzialgleichungen vorgenommen. Des Weiteren werden Anfangswertaufgaben und Randwertaufgaben vorgestellt. Bei der nachfolgenden Zusammenfassung wurde sich insbesondere der Literatur [3], [22], [25] bedient. 1.2.1 Definitionen • Wenn x und y zwei variable Gr¨oßen sind und wenn sich einem gegebenen x-Wert genau ein y-Wert zuordnen l¨asst, dann nennt man y eine Funktion von x und schreibt y=f(x). • Die ver¨anderliche Gr¨oße x heißt unabh¨angige Variable oder Argument der Funktion y. Die ver¨anderliche Gr¨oße y heißt abh¨angige Variable. • Differenzialgleichung wird eine Gleichung genannt, in der neben einer oder mehreren unabh¨angigen Ver¨anderlichen und einer oder mehreren Funktionen dieser Ver¨anderlichen auch noch die Ableitung dieser Funktion nach den unabh¨angigen Ver¨anderlichen auftreten. Die Ordnung einer Differenzialgleichung ist gleich der Ordnung der h¨ochsten in ihr auftretenden Ableitung. • Eine Gleichung, in der Ableitungen einer Funktion y = y(x) bis zur n-ten Ordnung auftreten, heißt gew¨ohnliche Differenzialgleichung n-ter Ordnung. • Partielle Differenzialgleichungen enthalten partielle Ableitungen einer Funktion von mehreren Variablen. • Eine Differenzialgleichung heißt linear, wenn die Funktion und ihre Ableitungen nur linear, d. h. in der ersten Potenz, auftreten. • Eine Differenzialgleichung heißt homogen, wenn in allen Termen x oder die Ableitung von x vorkommt. Andernfalls heißt sie inhomogen. • Eine Funktion, deren Gleichung nach der abh¨angigen Ver¨anderlichen aufgel¨ost ist, heißt explizit (explicitus (lat.) auseinandergewickelt). Die allgemeine Form der expliziten Funktion ist y = f(x). Bei der expliziten Form einer mathematischen Funktion lassen sich deren Werte ohne Umformen der Funktion unmittelbar berechnen (Bsp.: y = √ 1 − x 2 ). • Eine Funktion, deren Gleichung nicht nach der abh¨angigen Ver¨anderlichen aufgel¨ost ist, heißt implizit (implicitus (lat.) hineingewickelt). Die allgemeine Form <?page no="18"?> 4 Erforderliche mathematische Grundlagen einer impliziten Funktion ist f(x,y) = 0. Die implizite Form einer Gleichung f(x, y) = 0 erh¨alt man, wenn sich diese Gleichung eindeutig nach y aufl¨osen l¨asst (Bsp.: x 2 − y 2 − 1 = 0). In Tab. 1.2 sind Beispiele von Differenzialgleichungen aufgef¨ uhrt. Tabelle 1.2: Beispiele zur Darstellung und Benennung von Differenzialgleichungen y ′ = 2 x Explizite Dgl 1. Ordnung x + yy ′ = 0 Implizite Dgl 1. Ordnung y ′ + yy ′′ = 0 Implizite Dgl 2. Ordnung ¨ s = − g Explizite Dgl 2. Ordnung y ′′′ + 2 y ′ = cos(x) Implizite Dgl 3. Ordnung y (6) − y (4) + y ′′ = e x Implizite Dgl 6. Ordnung 1.2.2 Partielle Differenzialgleichungen Eine allgemeine Differenzialgleichung mit n-unabh¨angigen Variablen m-ter Ordnung in impliziter Form nennt man die Gleichung F ( x 1 , x 2 , ..., x n , u, ∂u ∂x 1 , ..., ∂u ∂x n , ∂ 2 u ∂x 2 1 , ∂ 2 u ∂x 1 ∂x 2 ..., ) = 0. Diese heißt partielle Differenzialgleichung. Ist m die h¨ochste Ordnung der darin vorkommenden partiellen Ableitung, so heißt die Gleichung partielle Differenzialgleichung m-ter Ordnung. Wird die Gleichung nach u m (x) aufgel¨ost, so erh¨alt man die explizite Form der gew¨ohnlichen Differenzialgleichung m-ter Ordnung. Empfehlenswerte Literatur: [3], S. 504 ff.; [25], S. 549 ff. Der Wert einer gemischten Ableitung ∂ 2 u ∂x 1 ∂x 2 ist f¨ ur gegebene Werte von x 1 und x 2 unabh¨angig von der Reihenfolge der Ableitungsbildung ∂ 2 u ∂x 1 ∂x 2 = ∂ 2 u ∂x 2 ∂x 1 <?page no="19"?> 1.2 Differenzialgleichungen 5 (Schwarz’scher Vertauschungssatz). Partielle Ableitungen h¨oherer Ordnung sind analog definiert ([3], S. 410). Im Allgemeinen ist ∂ 2 u ∂x 1 ∂x 2 = ∂u ∂x 1 · ∂u ∂x 2 . Hier seien noch Hinweise f¨ ur die Schreibweise von Ableitungen gegeben. Es ist f ′ = ∂f ∂x = df dx . F¨ ur die zweite Ableitung folgt hieraus f ′′ = d 2 f dx 2 , was wiederum nicht mit (f ′ ) 2 = df dx · df dx verwechselt werden darf ([19], S. 327). 1.2.3 Partielle Integration Die Gleichung zur partiellen Integration ist ∫ u(x) v ′ (x) dx = u(x) v(x) − ∫ u ′ (x) v(x) dx. Diese Gleichung gilt f¨ ur bestimmte Integrale ∫ b a u(x) v ′ (x) dx = [u(x) v(x)] b a − ∫ b a u ′ (x) v(x) dx. In einigen F¨allen kann eine mehrmalige partielle Integration erforderlich werden. <?page no="20"?> 6 Erforderliche mathematische Grundlagen 1.2.4 Klassifikation von Differenzialgleichungen Zur weiteren Betrachtung werden die gebr¨auchlichen Abk¨ urzungen f¨ ur partielle Ableitungen nach Tab. 1.3 verwendet. Die allgemeine lineare partielle Differenzialgleichung zweiter Ordnung f¨ ur die Funktion f(x, y) der zwei unabh¨angigen Variablen x und y hat die Form A f xx + B f xy + C f yy + D f x + E f y + F f = G, (1.1) wobei die Koeffizienten A, B, C, D, E, F, G im Allgemeinen Funktionen von x und y sind. Diese partiellen Differenzialgleichungen werden je nach den Werten der Koeffizienten A, B, C in drei Gruppen eingeteilt: Tabelle 1.3: Abk¨ urzungen f¨ ur partielle Ableitungen Partielle Ableitung Abk¨ urzung ∂f ∂x f x ∂f ∂y f y ∂ 2 f ∂x 2 f xx ∂ 2 f ∂y 2 f yy • elliptische Differenzialgleichungen: Diese beschreiben Zust¨ande, d. h. abgeschlossene Prozesse, die nicht von der Zeit t abh¨angen (station¨are Vorg¨ange). F¨ ur die L¨osung elliptischer Differenzialgleichungen gilt das Extremum-Prinzip. Das Maximum oder Minimum der L¨osung wird an den R¨andern angenommen und nicht im Inneren des Gebiets. • parabolische Differenzialgleichungen: Sie beschreiben Ausgleichsprozesse, die von der Zeit t abh¨angen (instation¨are Prozesse). Parabolische Probleme sind typischerweise Anfangsbzw. Randwertprobleme. Sie beschreiben thermische oder magnetische Diffusionsprozesse. F¨ ur parabolische Differenzialgleichungen gilt ebenfalls das Extremum-Prinzip. • hyperbolische Differenzialgleichungen: Dieser Typ beschreibt Wellenausbreitungen und Transportvorg¨ange, die von der Zeit t abh¨angen. Hyperbolische <?page no="21"?> 1.2 Differenzialgleichungen 7 Probleme sind reine Anfangswertprobleme. In einem endlichen Gebiet k¨onnen damit Randwerte nicht beliebig vorgegeben werden, sondern werden durch Vertr¨aglichkeitsbedingungen ersetzt. Mit Hilfe von Tab. 1.4 kann die Klassifizierung von Differenzialgleichungen vorgenommen werden. Da die Koeffizienten A, B und C Funktionen von x und y sind, kann die partielle Differenzialgleichung Gl. (1.1) in einem bestimmten Teilgebiet elliptisch, in einem anderen Teilgebiet parabolisch oder hyperbolisch sein [25], S. 463. Tabelle 1.4: Klassifizierung von Differenzialgleichungen elliptisch B 2 − 4AC < 0 parabolisch B 2 − 4AC = 0 hyperbolisch B 2 − 4AC > 0 1.2.5 Anfangswertaufgabe Werden der L¨osung y = y(x) einer gew¨ohnlichen Differenzialgleichung n-ter Ordnung an einer Stelle x 0 die n-Werte y(x 0 ), y ′ (x 0 ), y ′′ (x 0 ), ..., y n −1 (x 0 ) vorgegeben, spricht man von einer Anfangswertaufgabe. Die vorgegebenen Werte werden als Anfangswerte oder Anfangsbedingungen bezeichnet [3], S. 504. 1.2.6 Randwertaufgabe Die allgemeine L¨osung einer gew¨ohnlichen Differenzialgleichung n-ter Ordnung enth¨alt n freie Integrationskonstanten, die in einer speziellen L¨osung durch Anfangs- oder Randbedingungen festgelegt werden. Werden an die L¨osung einer gew¨ohnlichen Differenzialgleichung Bedingungen an mehreren Stellen des Definitionsbereichs gestellt, so nennt man diese Bedingungen Randbedingungen. Eine Differentialgleichung mit Randbedingungen heißt Randwertaufgabe [3], S. 504. Dies sind zus¨atzliche Bedingungen, die die spezielle L¨osung in einem oder mehreren Punkten erf¨ ullen m¨ ussen. Die allgemeine L¨osung einer partiellen Differenzialgleichung enth¨alt dagegen <?page no="22"?> 8 Erforderliche mathematische Grundlagen im allgemeinen willk¨ urliche Funktionen als Integrationskonstanten. Um diese festzulegen, m¨ ussen Randbedingungen l¨angs der Berandung Γ vorgegeben werden. Es sind die Randbedingungen • Randbedingung 1. Art: u = γ auf Γ (Dirichlet-Bedingung) • Randbedingung 2. Art: ∂u/ ∂n = γ auf Γ (Neumann-Bedingung) • Randbedingung 3. Art: ∂u/ ∂n + βu = γ auf Γ (Cauchy-Bedingung) zu unterscheiden [22], S. 18 f; [25], S. 464 f. Die Richtungsableitung in der Randbedingung ∂u/ ∂n ist definiert durch ∂u ∂n = n grad u = n 1 u x + n 2 u y . Abbildung 1.1: Beispiel zur Normalenableitung eines rechteckigen Gebiets In Abb. 1.1 ist ein rechteckiges Gebiet mit der Berandung Γ 1 bis Γ 4 zu sehen. F¨ ur n = n(x, y) auf dem Rand Γ 1 gilt ∂u ∂n = − ∂u ∂x und damit n = (-1, 0). F¨ ur n auf dem Rand Γ 2 gilt <?page no="23"?> 1.2 Differenzialgleichungen 9 ∂u ∂n = ∂u ∂y und damit n = (0, 1). F¨ ur n auf dem Rand Γ 3 gilt ∂u ∂n = ∂u ∂x und damit n = (1, 0). F¨ ur n auf dem Rand Γ 4 gilt ∂u ∂n = − ∂u ∂y und damit n = (0, -1). 1.2.7 Inneres Produkt Das innere Produkt wird als die Verallgemeinerung des Punktprodukts eingef¨ uhrt. Das innere Produkt entspricht 〈 f, g 〉 = ∫ ∫ Ω f g dx dy der Multiplikation zweier Vektoren oder Funktionen mit anschließender Integration ¨ uber das Gebiet Ω. Das Ergebnis ist immer ein Skalar. Aus typografischen Gr¨ unden erfolgt die Darstellung des inneren Produkts mit den beiden Klammern 〈 〉 . • Inneres Produkt von Vektoren: Das innere Produkt, auch Punktprodukt oder Skalarprodukt genannt, beschreibt die skalare Multiplikation von Vektoren. Das Ergebnis ist ein Skalar. Hierbei ist a = (a 1 , a 2 , ..., a n ) und b = (b 1 , b 2 , ..., b n ). Das innere Produkt beider Vektoren erfolgt mit der Schreibweise a · b = 〈 a, b 〉 = a 1 b 1 + a 2 b 2 + a 3 b 3 ([20], S. 24, S. 65). Gezeigt wird der besondere Fall, bei dem das Skalarprodukt gleich dem Integral <?page no="24"?> 10 Erforderliche mathematische Grundlagen 〈 a, b 〉 = ∫ ∫ Ω a b dx ist. Als Beispiel wird das innere Produkt der ortsabh¨angigen Vektoren a(s) und b(s) 〈 a, b 〉 = 〈 ⎛ ⎜ ⎜ ⎜ ⎝ 1 2, 5 2 3 ⎞ ⎟ ⎟ ⎟ ⎠ , ⎛ ⎜ ⎜ ⎜ ⎝ 4 1, 5 1, 5 1 ⎞ ⎟ ⎟ ⎟ ⎠ 〉 = 1 · 4 · m 2 + 2, 5 · 1, 5 · m 2 + 2 · 1, 5 · m 2 + 3 · 1 · m 2 = 4 m 2 + 3, 75 m 2 + 3 m 2 + 3 m 2 = 13, 75 m 2 gebildet. Die Integration erfolgt in diesem Beispiel entlang der Strecke s. Die Annahme besteht darin, dass die ¨ortliche ¨ Anderung der Vektorkomponenten pro Meter stattfinden. Im Beispiel entspricht die integrale Wegstrecke 4 m. Die geometrische Interpretation entspricht der Fl¨ache unter dem ”Kurvenverlauf“. • Inneres Produkt von Funktionen: Als Beispiel soll das innere Produkt zweier Funktionen dienen, die aus der Spannung u(t) und dem Strom i(t) 〈 u, i 〉 = ∫ Ω u(t) i(t) dt = ∫ Ω P (t) dt = W el berechnet werden, mit Ω ∈ [0, T ]. Die Multiplikation beider Funktionen mit anschließender Integration ¨ uber die Zeit t f¨ uhrt zur elektrischen Energie W el . F¨ ur eine runde Einf¨ uhrung in die Thematik des inneren Produktes siehe auch [4], S. 162 und S. 166. <?page no="25"?> 1.3 Vektoroperatoren 11 1.2.8 Starke Form/ Formulierung einer Differenzialgleichung Die gew¨ohnliche Differenzialgleichung L(u) beispielsweise in der Form L(u) = 0, x ∈ Ω d 2 u(x) dx 2 + 1 = 0 wird als starke Form (die f¨ ur uns bekannte Form) bezeichnet. 1.2.9 Schwache Form/ Formulierung einer Differenzialgleichung Die gew¨ohnliche Differenzialgleichung L(u) wird mit einer Funktion v(x) multipliziert und ¨ uber das Intervall Ω integriert 〈 L, v 〉 = ∫ Ω L(u) v(x) dx = ∫ Ω ( d 2 u(x) dx 2 + 1 ) v(x) dx = 0 und damit in ihre schwache Form ¨ uberf¨ uhrt. Dies kann mit der Schreibweise des inneren Produkts vereinfacht dargestellt werden. Die schwache Form erfordert nur eine einfache Stetigkeit der Ableitung (schwache Stetigkeitsforderung). 1.3 Vektoroperatoren Vorgestellt werden die drei Vektoroperatoren • Nabla- und Laplaceoperator: erm¨oglicht vereinfachte Schreibweise f¨ ur folgende Vektoroperatoren • Gradient: Richtungsableitung einer skalaren Funktion • Divergenz: Untersucht den Fluss des Vektorfeldes (Fluss pro Volumeneinheit) hinsichtlich Flussquellen und Flusssenken • Rotation: Untersucht ein Vektorfeld nach Wirbeln in kartesischen Koordinaten. <?page no="26"?> 12 Erforderliche mathematische Grundlagen 1.3.1 Nabla- und Laplaceoperator Der Nabla-Operator (eng.: del-operator) wird beschrieben mit ∇ = ∂ ∂x 1 e 1 + ∂ ∂x 2 e 2 + ∂ ∂x 3 e 3 = ⎛ ⎜ ⎝ ∂ ∂x 1 ∂ ∂x 2 ∂ ∂x 3 ⎞ ⎟ ⎠ . Der Operator hat gleichzeitig die Eigenschaft eines Vektors und eines mathematischen Operators. Der Laplace-Operator Δ wird als Delta-Operator bezeichnet. Diese Bezeichnung darf nicht mit der englischen Bezeichnung f¨ ur Nabla-Operator (del-operator) verwechselt werden. Der Laplace-Operator Δ ist das skalare Produkt des Nabla-Operators Δ = ∇ · ∇ = ∇ 2 = ⎛ ⎜ ⎝ ∂ x 1 ∂ x 2 ∂ x 3 ⎞ ⎟ ⎠ · ⎛ ⎜ ⎝ ∂ x 1 ∂ x 2 ∂ x 3 ⎞ ⎟ ⎠ = ∂ 2 ∂x 2 1 + ∂ 2 ∂x 2 2 + ∂ 2 ∂x 2 3 mit sich selbst. 1.3.2 Vektoroperator Gradient Der Gradient der skalaren Ortsfunktion (Potenzialfunktion) ϕ ist in kartesischen Koordinaten grad ϕ = ∂ϕ ∂x 1 e 1 + ∂ϕ ∂x 2 e 2 + ∂ϕ ∂x 3 e 3 = ∇ ϕ = ⎛ ⎜ ⎝ ∂ ∂x 1 ∂ ∂x 2 ∂ ∂x 3 ⎞ ⎟ ⎠ · ϕ = ⎛ ⎜ ⎝ ∂ϕ ∂x 1 ∂ϕ ∂x 2 ∂ϕ ∂x 3 ⎞ ⎟ ⎠ . <?page no="27"?> 1.3 Vektoroperatoren 13 Der Gradient ist ein Vektorfeld. Er ist ein Maß f¨ ur die ¨ Anderung der skalaren Funktion ϕ(x 1 , x 2 , x 3 ) in Richtung der Koordinaten x 1 , x 2 , x 3 im betrachteten Punkt P. Der Gradient grad ϕ(x 1 , x 2 , x 3 ) zeigt stets in Richtung des gr¨oßten Anstiegs seiner Potenzialfunktion. Der Gradient steht immer senkrecht auf den Fl¨achen ϕ(x 1 , x 2 , x 3 ). 1.3.3 Vektoroperator Divergenz In kartesischen Koordinaten wird die Divergenz mit div B = ∂B 1 ∂x 1 + ∂B 2 ∂x 2 + ∂B 3 ∂x 3 = ∂ B 1 ∂x 1 e 1 + ∂ B 2 ∂x 2 e 2 + ∂ B 3 ∂x 3 e 3 = ∇ B berechnet. Die Divergenz ist ein Skalar. Durch die skalare Multiplikation des Operators mit einem Vektor wird das Skalarprodukt und damit die Divergenz des Vektors B gebildet. 1.3.4 Vektoroperator Rotation Die Rotation eines Vektorfeldes wird in kartesischen Koordinaten mit rot B = ( ∂B 3 ∂x 2 − ∂B 2 ∂x 3 ) e 1 + ( ∂B 1 ∂x 3 − ∂B 3 ∂x 1 ) e 2 + ( ∂B 2 ∂x 1 − ∂B 1 ∂x 2 ) e 3 berechnet. Die Berechnungsformel l¨asst sich in kartesischen Koordinaten durch eine Determinante rot B = ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ e 1 e 2 e 3 ∂ ∂x 1 ∂ ∂x 2 ∂ ∂x 3 B 1 B 2 B 3 ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ = ∇ × B darstellen. Die Rotation ist ein Vektor. 1.3.5 Gegen¨ uberstellung der Vektoroperatoren In der Tabelle 1.5 sind die Vektoroperatoren gegen¨ ubergestellt. Diese beinhaltet die Argumente sowie die Ergebnisse. <?page no="28"?> 14 Erforderliche mathematische Grundlagen Tabelle 1.5: Gegen¨ uberstellung der Vektoroperatoren Operator: grad ϕ div A rot A Argument Skalar Vektor Vektor Ergebnis Vektor Skalar Vektor Nablaoperator ∇ · ϕ ∇ · A ∇ × A 1.3.6 N¨ utzliche Normen Im Fortgang werden empfohlene und n¨ utzliche Normen zur Erstellung von Dokumentationen, Thesen, wissenschaftlicher Berichte und dergleichen aufgelistet. • DIN 1301-Teil 1: ”Einheiten - Einheitennamen, Einheitenzeichen“. In der Norm sind Einheiten des Internationalen Einheitensystems (SI) sowie weitere empfohlene Einheiten mit Gr¨oße, Einheitenname, Einheitenzeichen und Definitionen aufgef¨ uhrt [6]. • DIN 1301-Beiblatt: Das Beiblatt enth¨alt keine weiteren Normen, sondern zus¨atzliche Informationen zum Teil 1 [5]. • DIN 1302: ”Allgemeine mathematische Zeichen und Begriffe“. Diese Norm legt mathematische Zeichen und Begriffe sowie ihre Benennungen fest. [7]. • DIN 1303: ”Vektoren, Matrizen, Tensoren - Zeichen und Begriffe“. In der Norm werden Zeichen und Begriffe behandelt, die Vektoren, Matrizen und Tensoren betreffen. Dabei wird die algebraische Struktur dargestellt [8]. • DIN 1304-Teil 1: ”Formelzeichen - Allgemeine Formelzeichen“. In dieser Norm werden Formelzeichen f¨ ur physikalische Gr¨oßen festgelegt. Zudem sind allgemeine Formelzeichen aufgef¨ uhrt, die in Physik und Technik Anwendung finden [9]. • DIN 1338: ”Formelschreibweise und Formelsatz“. Diese Norm gilt f¨ ur die Schreibweise und den Satz mathematischer, physikalischer und chemischer Formeln. Sie dient dem Verfasser (Studierende, Korrektoren, ...) zur Erstellung guter Formels¨atze [10]. • DIN 4895 Teil 1: ”Orthogonale Koordinatensysteme - Allgemeine Begriffe“. Diese Norm befasst sich mit Koordinatensystemen in dreidimensionalen euklidischen <?page no="29"?> 1.3 Vektoroperatoren 15 R¨aumen und mit der Darstellung physikalischer Gr¨oßen in solchen Koordinatensystemen [12]. • DIN 4895 Teil 2: ”Orthogonale Koordinatensysteme - Differentialoperatoren der Vektoranalysis“. In dieser Norm werden Differenzialoperatoren der Vektoranalysis in ihren orthogonalen Koordinaten dargestellt [11]. <?page no="30"?> Kapitel 2 Differenzialgleichungen und Finite Elemente Vorgestellt werden h¨aufig verwendete Differenzialgleichungen zweiter Ordnung und deren L¨osung mit der Galerkin-Methode als meist verwendeter Methode innerhalb der Finite-Element-Methode (FEM). Zudem werden die im FEM-Einsatz verwendeten Elemente gezeigt. Die FEM z¨ahlt zu den dominierenden numerischen Methoden zur Berechnung von Gleichungssystemen. Sie hat ihren Ursprung in der Strukturanalysis. Eine mathematische Beschreibung der Methode erfolgte bereits 1943. Erst 1968 erfolgte die Anwendung der Methode auf elektromagnetische Feldprobleme. Die FEM hat sich als leistungsf¨ahiger und flexibler in ihrer Anwendung auf komplexe Geometrien und inhomogene Stoffen im Vergleich zu anderen Methoden erwiesen und damit durchgesetzt. Aus diesem Grund konnten Computerprogramme entwickelt werden, welche eine große Vielfalt von Anwendungsm¨oglichkeiten schufen. Viele Vorg¨ange in Naturwissenschaft und Technik werden mittels partieller Differenzialgleichungen erster und zweiter Ordnung beschrieben. Die folgenden Abhandlungen des Manuskripts wurden im Wesentlichen mit der Literatur [1], [2], [13], [14] und [21] erstellt. 2.1 Physik-Beispiele f¨ ur Differenzialgleichungen erster Ordnung Differenzialgleichungen erster Ordnung enthalten eine Ableitung der unabh¨angigen Variable. Bei Ableitungen nach der Zeit verk¨orpern diese oft einen Energiespeicher. Beispiele von Differenzialgleichungen erster Ordnung sind in Abb. 2.1 zusammengefasst. 16 <?page no="31"?> 2.2 Physik-Beispiele f¨ ur Differenzialgleichungen zweiter Ordnung 17 Abbildung 2.1: Beispiele von Differenzialgleichungen 1’ter Ordnung 2.2 Physik-Beispiele f¨ ur Differenzialgleichungen zweiter Ordnung Beispiele von Differenzialgleichungen zweiter Ordnung sind in Abb. 2.2 zusammengefasst. In Abb. 2.2 a) sind die Schwingungsgleichungen des elektrischen und mechanischen Schwingkreises zu entnehmen. Schwingungsf¨ahige Systeme sind durch zwei unabh¨angige Energieformen gekennzeichnet, in welche die Energie beim Austausch jeweils ¨ uberf¨ uhrt wird. Die zeitlichen Ableitungen in den Gleichungen verk¨orpern jeweils einen Energiespeicher. Bei den Differenzialgleichungen handelt es sich jeweils um einen ged¨ampften Reihenschwingkreis mit Erregung. Gel¨ost wird nach dem Strom i bzw. nach dem Weg x. Die Diffusionsgleichung ist Abb. 2.2 b) zu entnehmen. Gezeigt werden die Felddiffusions- und W¨armediffusionsgleichung. Die Gleichungen sind durch zwei Orts- und eine Zeitableitung charakterisiert. Gel¨ost wird nach der Flussdichte B bzw. nach der Temperatur υ. Die Poisson’sche Differenzialgleichung der Elektrostatik zur L¨osung <?page no="32"?> 18 Differenzialgleichungen und Finite Elemente des Potenzials ϕ ist der Abb. 2.2 c) zu entnehmen. Die Wellengleichung in Abb. 2.2 d) ist mit zwei Ortsableitungen und zwei Zeitableitungen charakterisiert. Gel¨ost wird jeweils nach der elektrischen Feldst¨arke E bzw. nach der magnetischen Feldst¨arke H. Abbildung 2.2: Beispiele von Differenzialgleichungen 2’ter Ordnung <?page no="33"?> 2.3 Finite Elemente 19 2.3 Finite Elemente Die Diskretisierung des L¨osungsraums oder L¨osungsgebietes erfordert die Aufteilung in einzelne R¨aume bzw. einzelne Gebiete, die finite Elemente. In Abb. 2.3 sind typische finite Elemente nach Raumdimensionen und Ordnung der Ansatzfunktion geordnet dargestellt. Abbildung 2.3: Klassifizierung gew¨ahlter finiter Elemente Zur Diskretisierung der Randwertprobleme im zweidimensionalem Gebiet finden vorzugsweise Dreiecks- oder Viereckselemente Anwendung. Bei Dreieckselementen werden h¨aufig elementweise lineare, quadratische oder kubische Ansatzfunktionen verwendet. <?page no="34"?> 20 Differenzialgleichungen und Finite Elemente Diese sind genau in einem Knoten gleich Eins und in allen verbleibenden Knoten gleich Null. Bei elementweise quadratischen Ansatzfunktionen werden zus¨atzlich noch Kantenmittelknoten verwendet. Hier werden quadratische Ansatzfunktionen definiert, welche in einem Dreiecksknoten gleich Eins und in den verbleibenden Eckknoten gleich Null sind. Bei den Viereckselementen werden bilineare, biquadratische oder quadratische Funktionen verwendet ([16], S. 227 ff.). <?page no="35"?> Kapitel 3 Von der Momentenmethode zur Galerkin-Methode Die Momentenmethode oder Method of Moments (MOM) ist eine numerische Methode zur L¨osung von Randwertproblemen vorzugsweise im Gebiet Elektromagnetismus. Gleich wie die Finite-Elemente-Methode ¨ uberf¨ uhrt die Momentenmethode Gleichungen eines vorliegenden Randwertproblems in eine Matrixgleichung zur L¨osung mit einem Computer. Eine geschlossene Formulierung der Methode wird in [15] vorgestellt. Sie wurde zu einer vorherrschenden Methode f¨ ur Berechnungen im Gebiet des Elektromagnetismus entwickelt. In [23] wird das Phasen-Diskretisierungsraster-Verfahren (PDRV) als eine Erweiterung der MOM zur Berechnung des Feld-Streuverhaltens elektrisch leitf¨ahiger K¨orper vorgestellt. 3.1 Grundprinzip der Momentenmethode Eine gelungene Einf¨ uhrung in die Momentenmethode ist in [15], S.5 ff. gegeben. Das Grundprinzip der MOM beruht auf der Umwandlung einer Funktionsgleichung (Randwertproblem) mittels numerischer Approximation in eine Matrixgleichung zur L¨osung mit bekannten Vorgehensweisen. Zur Veranschaulichung der Methode wird die inhomogene Gleichung L(f) = g betrachtet. Wobei L ein linearer Operator, f die unbekannte, zu bestimmende Funktion und g die bekannte Funktion ist, welche den Quellterm repr¨asentiert. Der L¨osungsbe- 21 <?page no="36"?> 22 Von der Momentenmethode zur Galerkin-Methode reich ist durch Ω gekennzeichnet. Zur L¨osungsfindung wird f in eine Reihe von Funktionen f 1 , f 2 , f 3 , ... im Berich von Ω f = N ∑ n =1 α n f n (3.1) entwickelt. Dabei sind α n die unbekannten Entwicklungskoeffizienten und f n die Entwicklungs- oder Basisfunktionen. Die dabei entstehende endliche Summe muss entsprechend abgebrochen werden. Damit folgt die inhomogene Gleichung N ∑ n =1 α n L(f n ) = g. Folgend werden Wichtungs- oder Testfunktionen w 1 , w 2 , w 3 definiert und mit diesen das innere Produkt N ∑ n =1 α n 〈 w m , L(f n ) 〉 = 〈 w m , g 〉 mit m = 1, 2, 3, ... gebildet. Die Gleichung wird in Matrixschreibweise [l mn ] [α n ] = [g m ] notiert. Dabei sind [l mn ] = ⎛ ⎜ ⎝ 〈 w 1 , L(f 1 ) 〉 〈 w 1 , L(f 2 ) 〉 ... 〈 w 2 , L(f 1 ) 〉 〈 w 2 , L(f 2 ) 〉 ... ... ... ... ⎞ ⎟ ⎠ [α n ] = ⎛ ⎜ ⎜ ⎜ ⎝ α 1 α 2 . . ⎞ ⎟ ⎟ ⎟ ⎠ [g m ] = ⎛ ⎜ ⎜ ⎜ ⎝ 〈 w 1 , g 〉 〈 w 2 , g 〉 . . ⎞ ⎟ ⎟ ⎟ ⎠ . <?page no="37"?> 3.2 Anmerkungen zur Momentenmethode 23 Unter der Voraussetzung der Nichtsingularit¨at der Matrix l ist die Bestimmung der Koeffizienten α mit [α n ] = [l −1 nm ] [g m ] m¨oglich. Dies entspricht der L¨osung von f. Die Ausdr¨ ucke zur L¨osung von f werden mit [f n ] = [f 1 f 2 f 3 ...] und f = [f n ] [α n ] = [f n ] [l −1 nm ] [g m ] gek¨ urzt wiedergegeben. Ob die L¨osung exakt oder approximiert ist, h¨angt von der Wahl von f n und w n ab. 3.2 Anmerkungen zur Momentenmethode Im Folgenden werden dem Anwender der MOM Hinweise und Empfehlungen hinsichtlich der Matrix, der Wahl der Basis- und Wichtungsfunktion gegeben. 3.2.1 Matrix [l mn ] Ist die Matrix [l mn ] unendlicher Ordnung, so kann ihre Invertierung nur f¨ ur einige spezielle F¨alle wie beispielsweise eine Diagonalmatrix durchgef¨ uhrt werden. Die klassische Eigenfunktionsmethode f¨ uhrt auf eine Diagonalmatrix und kann als spezieller Fall der MOM aufgefasst werden. Sind dagegen f n und w n endlich, so nimmt die Matrix eine endliche Ordnung an und kann mittels bereits bekannten Methoden invertiert werden. 3.2.2 Wahl der Basis- und Wichtungsfunktionen f n und w n Eine der Hauptaufgaben aller spezieller Probleme ist die Wahl der Basisf n und Wichtungsfunktion w n . Die Basisfunktion sollte linear unabh¨angig sein und dabei so gew¨ahlt werden, dass die Superponierung von Gl. (3.1) die Funktion f angemessen und gut approximiert. Die Wichtungsfunktion w n sollte ebenfalls linear unabh¨angig und so <?page no="38"?> 24 Von der Momentenmethode zur Galerkin-Methode gew¨ahlt werden, dass das innere Produkt 〈 w n , g 〉 relativ unabh¨angig von der Eigenschaft von g ist. Einige weitere Faktoren, welche die Wahl von f n und w n beeinflussen, sind • die gew¨ unschte Genauigkeit der L¨osung, • der Aufwand zur Entwicklung der Matrixelemente, • die zu invertierende Matrixgr¨oße, • die Realisierung einer wohl konditionierten Matrix. Dieses Kapitel informiert ¨ uber die Idee des Mathematikers Galerkin und f¨ uhrt in dessen Methode ein. Die Galerkin-Methode findet eine sehr breite Anwendung zur L¨osung von Differenzialgleichungen beispielsweise in der Strukturmechanik, der Str¨omungsmechanik, dem W¨arme- und Massetransport, der Akustik sowie in Mikrowellenanwendungen. Probleme, beschrieben durch gew¨ohnliche Differenzialgleichungen, partielle Differenzialgleichungen und Integralgleichungen k¨onnen mittels des Galerkin-Formalismus untersucht werden. ”Any problem for which governing equation can be written down is a candidate for a Galerkin method.“ [14], S. 1. Der Ursprung der Methode wird auf eine Ver¨offentlichung von Galerkin (1915) zur¨ uckgef¨ uhrt. Galerkin war ein Maschinenbauingenieur, der 1899 in in St. Petersburg seinen Abschluss machte [14]. 3.3 Galerkins Idee Die spezielle Wahl, als Wichtungsfunktion die Basisfunktion w n = f n zu nehmen, wird als Galerkin-Methode bezeichnet. 3.4 Traditionelle Galerkin-Methode Ein 2D-Problem wird durch die lineare Differenzialgleichung (DGL) L(u) = 0 (3.2) <?page no="39"?> 3.4 Traditionelle Galerkin-Methode 25 im Gebiet D(x,y) mit den Randbedingungen S(u) = 0, ∂D beschrieben. Die Galerkin-Methode nimmt an, dass u durch die N¨aherungsl¨osung u a = u 0 (x, y) + N ∑ j =1 a j φ j (x, y) (3.3) genau dargestellt wird. Dabei gilt: • φ j ist eine bekannte analytische Funktion (Basisfunktion) • u 0 wird eingef¨ uhrt, um die Randbedingungen zu erf¨ ullen • a j sind zu bestimmende Koeffizienten. Mit dem Einsetzen von Gl. (3.3) in Gl. (3.2) folgt L(u a ) = L(u 0 ) + N ∑ j =1 a j L(φ j ) = R(a 0 , a 1 , ..., a N , x, y) (3.4) = 0. Mit der Galerkin-Methode werden die unbekannte Koeffizienten a j durch die Entwicklung des inneren Produkts 〈 φ k , R 〉 = 0, k = 1, ..., N (3.5) im Gebiet D ∈ [0, 1], mit dem Residuum R und der aus Gl. (3.3) bekannten analytischen Funktion φ k (Wichtungsfunktion) bestimmt. Da im vorliegenden Beispiel eine lineare Differenzialgleichung zugrunde liegt, kann die Gl. (3.5) direkt als Matrizengleichung f¨ ur den Koeffizienten a j als N ∑ j =1 a j 〈 φ k , L(φ j ) 〉 = −〈 φ k , L(u 0 ) 〉 (3.6) geschrieben werden. Das L¨osen und Einsetzen von a j in die Gl. (3.3) f¨ uhrt zur N¨aherungsl¨osung von u a . Die Wirkung der Galerkin-Methode wird in den nachfolgend bearbeiteten Beispielen transparent. <?page no="40"?> 26 Von der Momentenmethode zur Galerkin-Methode 3.5 Galerkin-FEM-Methode ”The Galerkin finite-element method has been the most popular method of weighted residuals, used with piecewise polynomials of low degree, since the early 1970s.“ [14], S. 86. Die Galerkin-FEM dient zur L¨osung von Differenzialgleichungen ≥ 2’ter Ordnung. In der Galerkin-FEM-Methode findet eine abschnittsweise definierte lineare Wichtungs- oder Testfunktion Anwendung. In der Literatur wird diese oft als Form- oder Interpolations- oder Dreiecksfunktion nach Abb. 3.1 b) bezeichnet, welche zudem die einfachste Form darstellt. Abbildung 3.1: Interpolation finiter Elemente mittels Dreiecksfunktion Die Testl¨osung in einem eindimensionalem Gebiet x 1 ≤ x ≤ x N ist beispielsweise mit der globalen, f¨ ur das gesamte Gebiet geltenden Gleichung <?page no="41"?> 3.6 Vorgehen zur L¨osung mit der Galerkin-Methode 27 u h = N ∑ i =1 u i φ i (x) (3.7) gegeben, wobei φ i (x) die Dreiecksfunktion und u i die zu l¨osenden nodalen Werte (Koeffizienten) darstellt. In Abb. 3.1 a) ist ersichtlich, dass u h die Funktion u linear zwischen den unbekannten Knotenwerten interpoliert und das f¨ ur jedes Element. Entnommen werden kann in Abb. 3.1 b) der lineare Abfall vom Wert Eins an einem jeweiligen Knoten auf den Wert Null bei den beiden benachbarten Knoten und dar¨ uber hinaus im verbleibendem Gebiet. Nur zwei Formfunktionen und zwei unbekannte Knotenwerte liefern einen von Null verschiedenen Beitrag zur Gl. (3.7). Beispielsweise liefern am Element 2 nur die Formfunktion φ 2 und φ 3 einen Beitrag zu u h der Gl. (3.7). Zudem kann der Abb. 3.1 entnommen werden, dass u a entlang der Elemente kontinuierlich und die Ableitung du h / dx an den Elementr¨andern diskontinuierlich verl¨auft. Weiterhin ist erkennbar, dass nur die Knotenwerte mit u ¨ ubereinstimmen. Mit der Einf¨ uhrung von u h ist demzufolge ein Interpolationsfehler verbunden. 3.6 Vorgehen zur L¨osung mit der Galerkin-Methode Der Leser wird in die L¨osung einer Differenzialgleichung mittels der Galerkin-Methode anhand eines 1D-Beispiels eingef¨ uhrt. Die Methode nach Galerkin ist in die Klasse der Methoden der gewichteten Residuen einzuordnen [14], S. 24. Die erforderlichen Bedingungen der Galerkin-Methode sind nach [14], S. 30: • Die Wichtungsfunktion w k ist von derselben Klasse wie die Basisfunktion φ. • Die Wichtungs- und Basisfunktionen sind bei der Galerkin-FEM linear unabh¨angig. • Die Basisfunktion sollte den Anfangs-, wie auch den Randbedingungen exakt gen¨ ugen. Die allgemeine Vorgehensweise zur Anwendung der Galerkin-Methode zur L¨osung einer partiellen Differenzialgleichung wird in die folgenden Schritte unterteilt: 1. ¨ Uberf¨ uhrung der starken Form der zu l¨osenden partiellen Differenzialgleichung in die schwache Form (schwache Formulierung) <?page no="42"?> 28 Von der Momentenmethode zur Galerkin-Methode 2. Diskretisierung des zu l¨osenden Gebietes Ω in eine endliche Anzahl n von Teilgebieten Ω n mit N Knoten bei der Galerkin-FEM 3. Wahl der Basis- und Wichtungsfunktionen φ 4. Formulierung der schwachen Form der erhaltenen Differenzialgleichung mittels gew¨ahlter Basis- und Wichtungsfunktion 5. ¨ Uberf¨ uhrung der Gleichung in eine Matrizengleichung 6. L¨osung des erhaltenen lineare Gleichungssystems <?page no="43"?> Kapitel 4 L¨osung der Gleichung dy dx − y = 0 mit der Galerkin-Methode Es sei die gew¨ohnlich DGL zweiter Ordnung dy dx − y = 0 (4.1) im Gebiet D mit 0 ≤ x ≤ 1 und der Randbedingung y(x=0) = 1 gegeben. 4.1 Wahl der Basis- und Wichtungsfunktion Im Fortgang wird die N¨aherungs- oder Ansatzfunktion y a = N ∑ j =0 a j x j = a 0 + N ∑ j =1 a j x j , (4.2) und die Wichtungsfunktion W = x k −1 eingef¨ uhrt. Wobei a 0 = 1 und die Basisfunktion x j die Randbedingung erf¨ ullt. Vergleiche hierzu Gl. (3.3). Mit dem Einsetzen von Gl. (4.2) in Gl. (4.1) folgt d dx ( a 0 + N ∑ j =1 a j x j ) − ( a 0 + N ∑ j =1 a j x j ) = 0 = R. 29 <?page no="44"?> 30 L¨osung der Gleichung dy dx − y = 0 mit der Galerkin-Methode Nach dem Ableiten verbleibt N ∑ j =1 a j jx j −1 − 1 − N ∑ j =1 a j x j = 0 − 1 + N ∑ j =1 a j ( jx j −1 − x j ) = 0 = R. 4.2 Formulierung der schwachen Form mit Basis- und Wichtungsfunktion Im Fortgang erfolgt die Entwicklung des inneren Produktes 〈 R, W 〉 = 0 ∫ 1 0 R x k −1 dx = 0, k = 1, ..., N ∫ 1 0 [ − 1 + N ∑ j =1 a j ( jx j −1 − x j ) ] x k −1 dx = 0 ∫ 1 0 [ − x k −1 + N ∑ j =1 a j ( jx j −1 − x j ) x k −1 ] dx = 0 ∫ 1 0 N ∑ j =1 a j ( jx j −1 − x j ) x k −1 dx = ∫ 1 0 x k −1 dx. Durch Einsetzen der Integrationsgrenzen j j + k − 1 x j + k −1 ∣ ∣ ∣ 1 0 − 1 j + k x j + k ∣ ∣ ∣ 1 0 = 1 k x k ∣ ∣ ∣ 1 0 j j + k − 1 − 1 j + k = 1 k folgt die schwache Formulierung von Gl. (4.1). <?page no="45"?> 4.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung 31 4.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Die Gleichungen werden in ein System linearer Gleichungen M A = D mit den Elementen von M und D m kj = 〈 jx j −1 − x j , x k −1 〉 = j j + k − 1 − 1 j + k d k = 〈 1, x k −1 〉 = 1 k ¨ uberf¨ uhrt. A ist der Vektor der unbekannten Koeffizienten a j . Die einzelnen Elemente der Matrix M werden nun mit N = 3 berechnet. Beispielhaft erfolgt die Berechnung der Elemente m(3, j) f¨ ur k = 3: m 31 = 1 1 + 3 − 1 − 1 1 + 3 = 1 12 m 32 = 2 2 + 3 − 1 − 1 2 + 3 = 3 10 m 33 = 3 3 + 3 − 1 − 1 3 + 3 = 13 30 . 4.4 L¨osung des linearen Gleichungssystems Das so erhaltene lineare Gleichungssystem ⎛ ⎜ ⎝ 1 2 2 3 3 4 1 6 5 12 11 20 1 12 3 10 13 30 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ M · ⎛ ⎜ ⎝ a 1 a 2 a 3 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ A = ⎛ ⎜ ⎝ 1 1 2 1 3 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ D wird mit A = M −1 D <?page no="46"?> 32 L¨osung der Gleichung dy dx − y = 0 mit der Galerkin-Methode gel¨ost. Die Galerkin-Methode ¨ uberf¨ uhrt die gew¨ohnliche DGL in ein algebraisches Gleichungssystem. Die L¨osung der Koeffizienten ist A = ⎛ ⎜ ⎝ a 1 a 2 a 3 ⎞ ⎟ ⎠ = ⎛ ⎜ ⎝ 1, 014 0, 423 0, 828 ⎞ ⎟ ⎠ . Aus der anschließenden Substitution von Gl. (4.2) folgt die N¨aherungsl¨osung y a = 1 x 0 + 1, 014 x 1 + 0, 423 x 2 + 0, 282 x 3 . Die exakte L¨osung ist y = e x . In Abb. 4.1 wurde die Gegen¨ uberstellung zwischen der exakten und der N¨aherungsl¨osung vorgenommen. 0 0.2 0.4 0.6 0.8 1 x 1 1.2 1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 y(x) analytische Lösung y = 1 x 0 + 1.014 x 1 + 0.423 x 2 + 0.282 x 3 y = 1 x 0 + 1.014 x 1 + 0.423 x 2 y = 1 x 0 + 1.014 x 1 Abbildung 4.1: Gegen¨ uberstellung der analytischen und der numerischen Funktionsverl¨aufe <?page no="47"?> Kapitel 5 L¨osung der Gleichung u ( x ) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode Im folgenden Beispiel wird eine Differenzialgleichung zweiter Ordnung mittels Galerkin- Methode auf ein System algebraischer Gleichungen reduziert. Beginne mit dem Beispiel einer quadratischen Gleichung (nach unten ge¨offnete Parabel) u(x) = − 1 2 x 2 + 1 2 x = 1 2 ( x − x 2 ) , (5.1) bei welcher sich die Koeffizienten nur durch das Vorzeichen unterscheiden, der Scheitelpunkt S P (1/ 2 | 1/ 8) ist und die Nullstellen a = 0 und b = 1 aufweist. In Abb. 5.1 ist der exakte Verlauf der Gleichung ersichtlich. Die L¨osung erfolgt mittels • linearer Basis- und Wichtungsfunktion (Dreiecksfunktion) • nichtlinearer Basis- und Wichtungsfunktion. Beide Ergebnisse werden im Anschluss gegen¨ ubergestellt und diskutiert. 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion Die Gl. (5.1) soll mit Hilfe linearer Wichtungsfunktionen in eine Matrixgleichung ¨ uberf¨ uhrt und gel¨ost werden. 33 <?page no="48"?> 34 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode 0 0.2 0.4 0.6 0.8 1 x 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 u(x) analytische Lösung Abbildung 5.1: Verlauf der Gleichung 2’ter Ordnung u(x) 5.1.1 Schwache Formulierung der Differenzialgleichung Durch zweimaliges Differenzieren erh¨alt man die partielle Differenzialgleichung 2’ter Ordnung (Poisson’sche Differenzialgleichung) in der starken Formulierung d 2 u(x) dx 2 + 1 = 0, x ∈ Ω (5.2) u(x) = 0, x ∂ Ω. (5.3) Gesucht ist die Funktion u(x) auf dem Intervall Ω = [a,b] = [0,1]. Die starke Form der Gl. (5.2) wird folgend in die schwache Form zur L¨osung mit der Galerkin-Methode ¨ uberf¨ uhrt ∫ Ω R W dx = 〈 R, W 〉 = 0. <?page no="49"?> 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion 35 Es sind W = Wichtungs- oder Testfunktion und R = Residuum. Mit R = d 2 u(x) dx 2 + 1 W = W (x) folgt ∫ Ω ( d 2 u(x) dx 2 + 1 ) W (x) dx = 0 ∫ Ω d 2 u(x) dx 2 W (x) dx + ∫ Ω W (x) dx = 0. Mit partieller Integration W (x) du(x) dx − ∫ Ω du(x) dx dW (x) dx dx + ∫ Ω W (x) dx = 0 und Einbezug der Dirichlet-Randbedingungen an den ¨außeren Knoten (R¨ander) u(x) = W (x) = 0; x ∈ ∂Ω wird der erste Term der Gleichung gem. Kap. 1.2.3 W (x) du(x) dx ∣ ∣ ∣ b a = 0. Dies ist der Fall, da die Wichtungsfunktionen nur f¨ ur die inneren Knoten Anwendung finden und an den ¨außeren Knoten den Wert Null annehmen. Daraus folgt die schwache Form der Differenzialgleichung (5.1) ∫ Ω du(x) dx dW (x) dx dx − ∫ Ω W (x) dx = 0. (5.4) 5.1.2 Diskretisierung des zu l¨osenden Gebiets Ω Das Intervall [a,b] wird in n Teilintervalle Ω n mit N Knoten unterteilt. F¨ ur die grafische Darstellung gem. Abb. 5.2 sind dies n = 5 Teilintervalle (Elemente, Teilgebiete) mit N = 6 Knoten. Die Intervallgrenzen werden mit x 0 = a und x 6 = b gesetzt und als ¨außere Knoten bezeichnet. <?page no="50"?> 36 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode Abbildung 5.2: 1D-Diskretisierung des Gebiets Ω in n Teilgebiete Ω i 5.1.3 Wahl der Basis- und Wichtungsfunktion Gew¨ahlt werden lineare Funktionen vom Typ Geradengleichungen, welche gem. Abb. 5.3 den einzelnen Knoten zugeordnet werden. Aus Gr¨ unden der Anschauung wird diese Funktion in der Literatur als Dreiecksfunktion bezeichnet. Die Dreiecksfunktion wird f¨ ur einen Knoten definiert, an welchem sie den Wert Eins annimmt. Die Dreiecksfunktion besitzt an dieser Stelle eine Unstetigkeit, was eine abschnittsbzw. elementweise Funktionsdefinition erfordert. Den Funktionswerten aller benachbarten Knoten wird der Wert Null zugeordnet. Der Vorteil dieses Typs besteht in deren Ableitungen, welche eine Konstante ergeben. Jedem inneren Knoten x i , i = 1, ... , N-1 wird eine Dreiecksfunktion φ i (x) mit φ i (x) = ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ 1 x i − x i−1 x − x i−1 x i − x i−1 = x − x i−1 h , x i −1 < x ≤ x i −1 x i+1 − x i x + x i+1 x i+1 − x i = x i+1 − x h , x i < x ≤ x i +1 0, sonst. zugeordnet. An den R¨andern des Gebiets Ω sind keine Basisfunktionen definiert. Dabei ist h = Ω n der Abstand zwischen zwei Knoten (Elementl¨ange). In Abb. 5.3 sind die Dreiecksfunktionen im Gebiet Ω andeutungsweise eingezeichnet. Mit Blick auf Gl. (5.4) werden Ableitungen der Funktion erforderlich. Somit schließen sich die abschnittsweisen Ableitungen der Dreiecksfunktionen mit <?page no="51"?> 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion 37 Abbildung 5.3: Nodale Zuordnung der Dreiecksfunktionen φ(i) mit deren Ableitungen im Gebiet Ω dφ i (x) dx = ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ 1 h , x i −1 < x ≤ x i −1 h , x i < x ≤ x i +1 0, sonst. an. Dreiecksfunktionen und deren Ableitungen nehmen mit dieser Definition außerhalb ihrer nodalen Zuordnungen den Wert Null an. 5.1.4 Formulierung der schwache Form mit Dreiecksfunktionen φ(x) Die Galerkin-Methode beinhaltet, dass die Basisfunktion gleich der Wichtungsfunktion ist. Im Fortgang wird die Dreiecksfunktion φ(x) als Basis- und Wichtungsfunktion verwendet. Die zu l¨osende Funktion u(x) wird durch die gen¨aherte Ansatzfunktion u h (x) f¨ ur ein 1D-Element mit den zwei Knoten x i und x i +1 u h (x) = 2 ∑ i =1 u i φ i (x) = u 1 φ 1 (x) + u 2 φ 2 (x) W (x) = φ(x) ersetzt. Dabei ist u i die gesuchte Variable, nach der zu l¨osen ist und φ i (x) die Basis- und φ(x) die Wichtungsfunktion. Eingesetzt in Gl. (5.4) folgt <?page no="52"?> 38 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode n ∑ i =1 [ ∫ Ω d(u i φ i (x)) dx dφ(x) dx dx ] − ∫ Ω φ(x) dx = 0 (5.5) mit anschließendem Ausmultiplizieren ∫ Ω [ u 1 dφ 1 (x) dx dφ(x) dx + u 2 dφ 2 (x) dx dφ(x) dx ] dx = ∫ Ω φ(x) dx. Die Methode nach Galerkin sieht vor, dass die Basisfunktion gleich der Wichtungsfunktion ist, und diese ein Produkt bilden. Insofern m¨ ussen im Fortgang die beiden Ableitungen der identischen und abschnittsweise definierten Basis- und Wichtungsfunktionen miteinander multipliziert werden. Dies betrifft jeweils die Multiplikation der Ableitungen der steigenden und fallenden Geraden der Dreiecksfunktionen. Hieraus entstehen zwei Gleichungen. Zur Bestimmung von u 1 wird der verbleibenden Funktion φ(x) die Funktion φ 1 (x) und zur Bestimmung von u 2 die Funktion φ 2 (x) zugeordnet. Damit folgt f¨ ur u 1 am Knoten x i u 1 ∫ Ω dφ 1 dx dφ 1 dx dx + u 2 ∫ Ω dφ 2 dx dφ 1 dx dx = ∫ Ω φ 1 (x) dx und f¨ ur u 2 am Knoten x i +1 u 1 ∫ Ω dφ 1 dx dφ 2 dx dx + u 2 ∫ Ω dφ 2 dx dφ 2 dx dx = ∫ Ω φ 2 (x) dx, oder in Matrizenschreibweise zusammengefasst ∫ Ω ( dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx ) dx ( u 1 u 2 ) = ∫ Ω φ(x) dx ( 1 1 ) . (5.6) Bei dem rechten Term der Integration der Funktion ¨ uber Ω ist eine Unterscheidung mittels Indizes 1, 2 nicht erforderlich, da eine derartige Unterscheidung das Integral nicht beeinflussen wird. 5.1.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Jedes Element wird mit zwei Knoten begrenzt. Im Fortgang werden mit Hilfe der Gl. (5.6) die beiden Knotenmatrizen hergeleitet, zu einer Elementmatrix zusammengef¨ uhrt <?page no="53"?> 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion 39 und diese f¨ ur alle Elemente in eine globale Matrix, die Koeffizientenmatrix ¨ uberf¨ uhrt. Dem Gedanken von Galerkin folgend ist ¨ uber die Ableitung der Basis- und Wichtungsfunktion zu integrieren. Beide sind gem. Abb. 5.3 abschnittsweise definiert und erfordern eine abschnittsweise Integration der Ableitung der steigenden und fallenden Flanken beider Dreiecksfunktionen. • Knotenmatrix des inneren Knotens x i mit Wichtungsfunktion φ = φ 1 : Der innere Knoten x i entspricht dem inneren Knoten x 1 der Abb. 5.3. Damit folgt das Integrationsintervall [x 0 , x 2 ]. Mit Hilfe von Gl. (5.6) folgt u 1 ∫ x 1 x 0 φ 1 dx φ 1 dx dx + u 2 ∫ x 1 x 0 φ 2 dx φ 1 dx dx+u 1 ∫ x 2 x 1 φ 1 dx φ 1 dx dx + u 2 ∫ x 2 x 1 φ 2 dx φ 1 dx dx = ∫ Ω φ 1 dx jeweils f¨ ur die Ableitungen der aufsteigenden und fallenden Flanken der Dreiecksfunktionen φ 1 und φ 2 . Das Integrationsintervall entspricht der Elementl¨ange h, dem Teilintervall Ω i . Die Ableitung der Basis- und Wichtungsfunktionen erfolgt mit Betrachtung von Abb. 5.3, womit u 1 ∫ x 1 x 0 1 h 1 h dx + u 2 ∫ x 1 x 0 0 1 h dx+u 1 ∫ x 2 x 1 − 1 h − 1 h dx + u 2 ∫ x 2 x 1 1 h − 1 h dx = ∫ Ω φ 1 dx erreicht wird. Im zweiten Term ist die Ableitung f¨ ur φ 2 im Intervall [x 0 , x 1 ] nicht definiert und nimmt daher den Wert Null an. Eine weitere Zusammenfassung f¨ uhrt zu u 1 ∫ x 1 x 0 1 h 2 dx + 0 + u 1 ∫ x 2 x 1 1 h 2 dx + u 2 ∫ x 2 x 1 − 1 h 2 dx = ∫ Ω φ 1 dx. Die Integration erfolgt gliedweise ¨ uber die Teilintervalle Ω i (Elementl¨ange h), bzw. bei der Funktion φ 1 ¨ uber das Intervall Ω mit ∫ Ω i 1 h 2 dx = ∫ h 0 1 h 2 dx = 1 h 2 ∫ h 0 dx = 1 h 2 x ∣ ∣ ∣ h 0 = 1 h ∫ Ω i − 1 h 2 dx = ∫ h 0 − 1 h 2 dx = − 1 h 2 ∫ h 0 dx = − 1 h 2 x ∣ ∣ ∣ h 0 = − 1 h ∫ Ω φ 1 dx = h. <?page no="54"?> 40 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode Vereinfachend wird die untere Integrationsgrenze des Intervalls Ω i mit Null und die obere Integrationsgrenze mit h angenommen. Durch Einsetzen der Integrationsergebnisse folgt die Knotenmatrix u 1 1 h + 0 + u 1 1 h − u 2 1 h = h 1 h ( 2 − 1 ) ( u 1 u 2 ) = h (5.7) f¨ ur den Knoten x 1 . • Knotenmatrix des inneren Knotens x i +1 mit Wichtungsfunktion φ = φ 2 : Der innere Knoten x i +1 entspricht dem inneren Knoten x 2 der Abb. 5.3. Daraus folgt das Integrationsintervall [x 1 , x 3 ]. Mit Hilfe von Gl. 5.6 folgt u 1 ∫ x 2 x 1 φ 1 dx φ 2 dx dx + u 2 ∫ x 2 x 1 φ 2 dx φ 2 dx dx+u 1 ∫ x 3 x 2 φ 1 dx φ 2 dx dx + u 2 ∫ x 3 x 2 φ 2 dx φ 2 dx dx = ∫ Ω φ 2 dx das identische Vorgehen wie beim Knoten x 1 . Die Ableitung der Basis- und Wichtungsfunktionen folgt aus der Betrachtung von Abb. 5.3, womit u 1 ∫ x 2 x 1 − 1 h 1 h dx + u 2 ∫ x 2 x 1 1 h 1 h dx+u 1 ∫ x 3 x 2 0 − 1 h dx + u 2 ∫ x 3 x 2 − 1 h − 1 h dx = ∫ Ω φ 2 dx erreicht wird. Im dritten Term ist die Ableitung f¨ ur φ 1 im Intervall [x 2 , x 3 ] nicht definiert und nimmt daher den Wert Null an. Eine weitere Zusammenfassung f¨ uhrt zu u 1 ∫ x 2 x 1 − 1 h 2 dx + u 2 ∫ x 2 x 1 1 h 2 dx + 0 + u 2 ∫ x 3 x 2 1 h 2 dx = ∫ Ω φ 2 dx. Die Integration erfolgt gliedweise ¨ uber die Teilintervalle wie beim Knoten x 1 . Durch Einsetzen der Integrationsergebnisse folgt die Knotenmatrizengleichung u 1 − 1 h + u 2 1 h + 0 + u 2 1 h = h 1 h ( − 1 2 ) ( u 1 u 2 ) = h (5.8) f¨ ur den Knoten x 2 . <?page no="55"?> 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion 41 • Elementmatrix: Stellvertretend wird f¨ ur das Element, begrenzt durch die inneren Knoten x 1 und x 2 unter Einbezug der Knotenmatrizengleichungen (5.7) und (5.8) die Elementmatrix und damit die Elementmatrizengleichung 1 h ( 2 − 1 − 1 2 ) · ( u 1 u 2 ) = h ( 1 1 ) erstellt. • Koeffizientenmatrix: Die einzelnen Elementgleichungen der inneren Knoten x 1 bis x 4 werden in der globalen Matrix (Koeffizientenmatrix) S durch Erweiterung 1 h ⎛ ⎜ ⎜ ⎜ ⎝ 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ S · ⎛ ⎜ ⎜ ⎜ ⎝ u 1 u 2 u 3 u 4 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ u h = h ⎛ ⎜ ⎜ ⎜ ⎝ 1 1 1 1 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ f in der Koeffizientenmatrizengleichung zusammengefasst. Bleiben noch die Dirichlet-Randbedingungen der ¨außeren Knoten zu ber¨ ucksichtigen. In der angenommenen Beispielgleichung wurde die starke Form der Gl. (5.2) durch zweimaliges Differenzieren von Gl. (5.1) erreicht. Leicht ersichtlich ist, dass damit auch Informationen verloren gingen, welche im Fortgang als Randbedingungen wieder hinzugef¨ ugt werden m¨ ussen, um den Verlauf der Gl. (5.1) in Abb. 5.1 anzun¨ahern. 5.1.6 L¨osung des linearen Gleichungssystems Das so erhaltene lineare Gleichungssystem 1 h S u h = h f wird nach S −1 (S u h ) = h 2 S −1 f <?page no="56"?> 42 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode Abbildung 5.4: Ergebnisdarstellung mittels Galerkin-Methode ( S −1 S ) ︸ ︷︷ ︸ E u h = h 2 S −1 f u h = h 2 S −1 f = h 2 ⎛ ⎜ ⎜ ⎜ ⎝ 2 3 3 2 ⎞ ⎟ ⎟ ⎟ ⎠ gel¨ost. Das Gleichungssystem ist unter Einbezug der Dirichlet-Randbedingungen nach Gl. (5.3) zu erweitern ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0 0 − 1 2 -1 0 0 0 0 -1 2 -1 0 0 0 0 -1 2 -1 0 0 0 0 -1 2 − 1 0 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ S · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ u 0 u 1 u 2 u 3 u 4 u 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ u h = h 2 ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 0 1 1 1 1 0 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ f . (5.9) <?page no="57"?> 5.1 L¨osung mit linearer Basis- und Wichtungsfunktion 43 0 0.2 0.4 0.6 0.8 1 x 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 u(x) analytische Lösung numerische Lösung Abbildung 5.5: Gegen¨ uberstellung des analytischen mit numerischen Funktionsverlaufs Mit diesem Vorgehen k¨onnen im Spaltenvektor f beliebige Randbedingungen vorgegeben werden, wie zum Beispiel der Funktionswert Null an den beiden ¨außeren Knoten. Damit folgen die Werte der inneren Knoten an den Stellen x i u h = 4 ∑ i =1 u i (x) φ i (x) = h 2 (2 φ 1 + 3 φ 2 + 3 φ 3 + 2 φ 4 ) . Die Basisfunktionen nehmen an den Stellen x i den Wert Eins an. Elementl¨ange, Randbedingungen und Ergebnisse der numerischen N¨aherung an den Verlauf von Gl. (5.1) sind in Tab. 5.1 zusammengefasst und in Abb. 5.5 grafisch dargestellt. Der Abb. 5.4 ist die grafische Ergebnisinterpretation der Galerkin-Methode zu entnehmen. Die Dirichlet-Randbedingungen u(x 0 ) = u(x 5 ) = 0 sind im Funktionsverlauf nicht eingezeichnet. <?page no="58"?> 44 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode Tabelle 5.1: Numerische Ergebnisse d. Gl. u(x)= 1 2 (x − x 2 ) Elementl¨ange h: 0,02 Dirichlet Bedingung linker Rand 0 Dirichlet Bedingung rechter Rand 0 Knotennummer x i , i = 0 1 2 3 4 5 Position 0 0,2 0,4 0,6 0,8 1,0 u(x i ) 0 0,08 0,12 0,12 0,08 0 5.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion Zur L¨osung der Differenzialgleichung wird eine nichtlineare Wichtungsfunktion gew¨ahlt. Beziehend auf die Gl. (5.2), welcher die Bedingung nach Gl. (5.3) auferlegt wird, folgt erneut 〈 W, R 〉 = ∫ Ω W R dx = ∫ Ω W ( d 2 u(x) dx 2 + 1 ) dx = 0, mit dem Residuum R und Wichtungs- oder Testfunktion W. 5.2.1 Wahl der Basis- und Wichtungsfunktion Im Fortgang wird u(x) = N ∑ n =1 α n ( x − x n +1 ) (5.10) W (x) = x − x m +1 nach Galerkin die zu l¨osende Funktion u(x) als Polynom mit der gleichen Funktionsklasse wie die Wichtungsfunktion angenommen. <?page no="59"?> 5.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion 45 5.2.2 Schwache Formulierung der Differenzialgleichung Daraus folgt die schwache Formulierung der Gl. (5.2) mit N ∑ n =1 α n ∫ Ω ( x − x m +1 ) d 2 dx 2 ( x − x n +1 ) dx = − ∫ Ω ( x − x m +1 ) dx. (5.11) 5.2.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Der linke Term der Gl. (5.11) wird mit Hilfe der partiellen Integration und der rechte Term durch Integration in die Form N ∑ n =1 α n ( n 2 + n n x n + m +1 − n + mn + n 2 m + n + 1 x n +1 ) ∣ ∣ ∣ ∣ 1 0 ︸ ︷︷ ︸ [ lmn ] = − ( 1 2 x 2 − 1 m + 2 x m +2 ) ∣ ∣ ∣ ∣ 1 0 ︸ ︷︷ ︸ [ g m ] ¨ uberf¨ uhrt. Das Gebiet wird mit Ω = [0; 1] gem. Abb. 5.5 eingegrenzt. 5.2.4 L¨osung des linearen Gleichungssystems Das so entstandene lineare Gleichungssystem [α] [l mn ] = [g m ] wird nach [α] gel¨ost. Es ist ⎛ ⎜ ⎝ α 1 α 2 α 3 ⎞ ⎟ ⎠ = ⎛ ⎜ ⎝ 80/ 3 256 8448/ 5 112 5504/ 5 7424 1968/ 5 3968 191232/ 7 ⎞ ⎟ ⎠ −1 ︸ ︷︷ ︸ [ l −1 mn ] · ⎛ ⎜ ⎝ 40/ 3 56 984/ 5 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ [ g m ] = ⎛ ⎜ ⎝ 1/ 2 0 0 ⎞ ⎟ ⎠ . Damit folgt die L¨osung der Differenzialgleichung von Gl. (5.2) in der Darstellung von Gl. (5.10) <?page no="60"?> 46 L¨osung der Gleichung u(x) = − 1 2 x 2 + 1 2 x mit der Galerkin-Methode u(x) = α 1 ( x − x 1+1 ) + α 2 ( x − x 2+1 ) + α 3 ( x − x 3+1 ) = 1 2 ( x − x 2 ) = − 1 2 x 2 + 1 2 x. <?page no="61"?> Kapitel 6 L¨osung der Gleichung u ( x ) = − 1 2 x 2 + 2 x mit der Galerkin-Methode In diesem Kapitel wird eine Differenzialgleichung zweiter Ordnung mittels Galerkin- Methode auf ein System algebraischer Gleichungen reduziert. Hierzu wird die quadratische Gleichung (nach unten ge¨offnete Parabel) mit Scheitelpunkt S P (2 | 2) u(x) = − 1 2 x 2 + 2 x (6.1) mit Nullstellen a = 0 und b = 4. verwendet. Die Gleichung ist durch zwei, auch im Vorzeichen unterscheidbare Koeffizienten gekennzeichnet. In Abb. 6.1 ist der exakte Verlauf der Gleichung ersichtlich. Die L¨osung erfolgt mittels • linearer Basis- und Wichtungsfunktion (Dreiecksfunktion) • nichtlinearer Basis- und Wichtungsfunktion. Das Vorgehen zur L¨osung erfolgt nach den in Kap. 3.6 definierten Schritten und ist identisch mit dem Vorgehen in Kap. 5. 47 <?page no="62"?> 48 L¨osung der Gleichung u(x) = − 1 2 x 2 + 2x mit der Galerkin-Methode 0 0.5 1 1.5 2 2.5 3 3.5 4 x 0 0.5 1 1.5 2 u(x) analytische Lösung Abbildung 6.1: Verlauf der Gleichung 2’ter Ordnung u(x) 6.1 L¨osung mit linearer Basis- und Wichtungsfunktion Durch zweimaliges Differenzieren der Gl. (6.1) folgt die partielle Differenzialgleichung 2’ter Ordnung (Poisson’sche Differenzialgleichung) in der starken Formulierung d 2 u(x) dx 2 + 1 = 0, x ∈ Ω (6.2) u(x) = 0, x ∂ Ω. (6.3) Es ist leicht zu erkennen, dass Gl. (6.2) damit wieder der Gl. (5.2) entspricht, und nur durch die Randbedingungen exakt bestimmbar ist. Gesucht ist im Fortgang die Funktion u(x) auf dem Intervall Ω = [a,b] = [0,4]. <?page no="63"?> 6.1 L¨osung mit linearer Basis- und Wichtungsfunktion 49 6.1.1 Schwache Formulierung der Differenzialgleichung Die starke Form der Gl. (6.2) wird folgend in die schwache Form zur L¨osung mit der Galerkin-Methode ¨ uberf¨ uhrt ∫ Ω W R dx = 〈 W, R 〉 = 0. Es sind W = Wichtungs- oder Testfunktion und R = Residuum. Aus R = d 2 u(x) dx 2 + 1 W = W (x) folgt ∫ Ω ( d 2 u(x) dx 2 + 1 ) W (x) dx = 0 ∫ Ω d 2 u(x) dx 2 W (x) dx + ∫ Ω W (x) dx = 0. Das Vorgehen nach Kap. 5.1 f¨ uhrt zu Gl. (5.4) ∫ Ω du(x) dx dW (x) dx dx − ∫ Ω W (x) dx = 0, der schwachen Form von Gl. (6.2). 6.1.2 Diskretisierung des zu l¨osenden Gebiets Ω Die Diskretisierung erfolgt nach Abb. 5.2. Das Gebiet Ω oder Intervall [a,b] wird erneut in f¨ unf Teilintervalle Ω n mit sechs Knoten diskretisiert, wobei x 0 und x 5 die ¨außeren Knoten bilden. 6.1.3 Wahl der Basis- und Wichtungsfunktion Die Wahl der Basis- und Wichtungsfunktion ist gleich der in Kap. 5.1.3 gew¨ahlten Dreiecksfunktionen. Die Funktionen und deren Ableitungen sind abschnittsweise definiert. Siehe hierzu auch Abb. 5.3. <?page no="64"?> 50 L¨osung der Gleichung u(x) = − 1 2 x 2 + 2x mit der Galerkin-Methode Tabelle 6.1: Numerische Ergebnisse d. Gl. u(x)= − 1 2 x 2 + 2x Elementl¨ange h: 0,8 Dirichlet Bedingung linker Rand 0 Dirichlet Bedingung rechter Rand 0 Knotennummer x i , i = 0 1 2 3 4 5 Position 0 0,8 1,6 2,4 3,2 4,0 u(x i ) 0 1,28 1,92 1,92 1,28 0 6.1.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) Unter Anwendung des Vorgehens nach Kap. 5.1.4 folgt die schwache Form der Differenzialgleichung in Matrizenschreibweise nach Gl. (5.6). 6.1.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Hier gen¨ ugt die Anwendung des Kapitels 5.1.5. Es werden die beiden Knotenmatrizen erstellt, in einer Elementmatrix zusammengef¨ uhrt und anschließend die Koeffizientenmatrix erstellt. 6.1.6 L¨osung des linearen Gleichungssystems Aus dem in Kap. 5.1.6 beschriebenen Vorgehen folgt erneut die Gl. (5.9). Die erforderlichen Werte und Ergebnisse sind in Tab. 6.1 zusammengefasst. Die ¨außeren Knoten erhalten gem. den Randbedingungen den Wert Null. Es verbleiben die inneren Knoten mit u h = 4 ∑ i =1 u i (x) φ i (x) = h 2 (2 φ 1 + 3 φ 2 + 3 φ 3 + 2 φ 4 ) . In Abb. 6.2 ist das Ergebnis grafisch mittels Dreiecksfunktionen (Basis- und Wichtungsfunktionen) dargestellt. In Abb. 6.3 ist die Gegen¨ uberstellung zwischen analytisch und <?page no="65"?> 6.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion 51 numerisch errechnetem Ergebnis ersichtlich. Abbildung 6.2: Ergebnisdarstellung mittels Galerkin-Methode 6.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion Es gilt die Gleichung (6.2) d 2 u(x) dx 2 + 1 = 0 = R durch Entwicklung des inneren Produkts 〈 W, R 〉 = 0 ∫ Ω W (x) ( d 2 u(x) dx 2 + 1 ) dx = 0 mit der Galerkin-Methode zu l¨osen. <?page no="66"?> 52 L¨osung der Gleichung u(x) = − 1 2 x 2 + 2x mit der Galerkin-Methode 0 0.5 1 1.5 2 2.5 3 3.5 4 x 0 0.5 1 1.5 2 u(x) analytische Lösung numerische Lösung Abbildung 6.3: Gegen¨ uberstellung des analytischen mit numerischen Funktionsverlaufs 6.2.1 Wahl der Basis- und Wichtungsfunktion Gem¨aß der Galerkin-Methode ist W (x) = u(x) zu setzen. Geeignete Gewichts- und Basisfunktionen m¨ ussen gesucht und gepr¨ uft werden. Die gesuchte Funktion muss stets die geforderten Randbedingungen W(x) = u(x) = 0 bei x ∂ Ω erf¨ ullen. Gesucht, gefunden und gepr¨ uft wurde u(x) = N ∑ n =1 α n ( x n − 1 4 x n +1 ) W (x) = x m − 1 4 x m +1 . 6.2.2 Formulierung der schwachen Form mit Basis- und Wichtungsfunktion Die gefundene Basis- und Wichtungsfunktion wird in die oben stehende Gleichung zur Berechnung des inneren Produkts eingesetzt und umgeformt. Daraus folgt <?page no="67"?> 6.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion 53 ∫ Ω ( x m − 1 4 x m +1 ) [ d 2 dx 2 N ∑ n =1 α n ( x n − 1 4 x n +1 ) + 1 ] = 0 α n N ∑ n =1 ∫ Ω ( x m − 1 4 x m +1 ) d 2 dx 2 ( x n − 1 4 x n +1 ) dx ︸ ︷︷ ︸ [ l mn ] = − ∫ Ω ( x m − 1 4 x m +1 ) dx ︸ ︷︷ ︸ [ g m ] die schwache Form der Differenzialgleichung. 6.2.3 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Die ¨ Uberf¨ uhrung der schwachen Form der Gleichung in die Matrizengleichung erfolgt mittels partieller Integration. Eine anschließende Bildung der Stammfunktion und Einsetzen der R¨ander x(0) = x(4) = 0 ergeben die beiden Matrizen • Matrix [l mn ] : ( − n x m + n −1 ( − n 3 x 2 + 8 n 3 x − 16 n 3 + n x 2 − 8 n x + 16 n) 16 (m 3 + 3 m 2 n + 3 m n 2 − m + n 3 − n) − m 2 n x m + n −1 (16 n − 8 n x + n x 2 + x 2 − 16) 16 (m 3 + 3 m 2 n + 3 m n 2 − m + n 3 − n) + m n x m + n −1 ( − 2 n 2 x 2 + 16 n 2 x − 32 n 2 − n x 2 + 16 n + x 2 + 16) 16 (m 3 + 3 m 2 n + 3 m n 2 − m + n 3 − n) ) ∣ ∣ ∣ ∣ 4 0 = − 2 4 m + n −1 m n m 3 + 3 m 2 n + 3 m n 2 − m + n 3 − n . F¨ ur m =[1 2 3] und n = [1 2 3] folgt die Matrix l mn [l mn ] = ⎛ ⎜ ⎝ −4 3 −8 3 −32 5 −8 3 −128 15 −128 5 −32 5 −128 5 −3072 35 ⎞ ⎟ ⎠ . Dabei ist det(l mn ) = -65536/ 2625 = 0. <?page no="68"?> 54 L¨osung der Gleichung u(x) = − 1 2 x 2 + 2x mit der Galerkin-Methode • Matrix [g m ] : − ∫ Ω ( x m − 1 4 x m +1 ) dx = − ( ∫ Ω x m dx − ∫ Ω 1 4 x m +1 dx ) = − ( 1 m + 1 x m +1 ∣ ∣ ∣ ∣ 4 0 − 1 4(m + 2) x m +2 ∣ ∣ ∣ ∣ 4 0 ) = − 1 m + 1 4 m +1 + 1 4(m + 2) 4 m +2 . F¨ ur m = [1,2,3] folgt [g m ] = ⎛ ⎜ ⎝ −8 3 −16 3 −64 5 ⎞ ⎟ ⎠ . • Matrix [α n ] : Die Matrix der zu l¨osenden Variabel α ist f¨ ur n = [1,2,3] [α n ] = ⎛ ⎜ ⎝ α 1 α 2 α 3 ⎞ ⎟ ⎠ . 6.2.4 L¨osung des linearen Gleichungssystems Das so entstandene lineare Gleichungssystem [α] [l mn ] = [g m ] wird nach [α] mit ⎛ ⎜ ⎝ α 1 α 2 α 3 ⎞ ⎟ ⎠ = ⎛ ⎜ ⎝ −4 3 −8 3 −32 5 −8 3 −128 15 −128 5 −32 5 −128 5 −3072 35 ⎞ ⎟ ⎠ −1 ︸ ︷︷ ︸ [ l mn ] · ⎛ ⎜ ⎝ −8 3 −16 3 −64 5 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ [ g m ] = ⎛ ⎜ ⎝ 2 0 0 ⎞ ⎟ ⎠ <?page no="69"?> 6.2 L¨osung mit nichtlinearer Basis- und Wichtungsfunktion 55 gel¨ost. Damit folgt die L¨osung der Differenzialgleichung von Gl. (6.2) in der Darstellung von Gl. (6.1) u(x) = α 1 ( x − 1 4 x 1+1 ) + α 2 ( x − 1 4 x 2+1 ) + α 3 ( x − 1 4 x 3+1 ) = 2 ( x − 1 4 x 2 ) = − 1 2 x 2 + 2 x. <?page no="70"?> Kapitel 7 L¨osung physik. Bsp. DGL 1’ter Ordnung mit Galerkin-Methode Folgend wird die ausgew¨ahlte Beispielsgleichungen, die Ampere’sche Gleichung aus Abb. 2.1, mit der Galerkin-Methode in eine Matrizengleichung ¨ uberf¨ uhrt und gel¨ost. Der Rechenweg erfolgt mittels zweier unterschiedlicher Basis- und Wichtungsfunktionen. 7.1 Durchflutungsgesetz gel¨ost mit nichtlinearer Basis- und Wichtungsfunktion Gegeben ist die Integralform des Durchflutungsgesetzes nach Abb. 7.1. Die analytische Herleitung der magnetischen Feldst¨arke f¨ ur den Innenraum, f¨ ur die Oberfl¨ache und dem Außenraum des Leiters wird dargestellt. In Abb. 7.2 ist die L¨osung der Gleichung nach der magnetischen Feldst¨arke f¨ ur den Innenraum und Außenraum des Leiters ersichtlich. Die maximale magnetische Feldst¨arke stellt sich an der Leiteroberfl¨ache ein. Gegeben ist die Differentialform des Durchflutungsgesetzes, die Rotation der magnetischen Feldst¨arke mit rot H = J rot z H = ⎛ ⎜ ⎝ 1 r ∂rH Φ ∂r − 1 r ∂H r ∂Φ ︸ ︷︷ ︸ =0 ⎞ ⎟ ⎠ 1 r ∂rH Φ ∂r = J, (7.1) 56 <?page no="71"?> 7.1 Durchflutungsgesetz gel¨ost mit nichtlinearer Basis- und Wichtungsfunktion 57 Abbildung 7.1: Herleitung der magnetischen Feldst¨arke notiert im Zylinderkoordinatensystem. Die Gleichung soll nach H Φ gel¨ost werden. Das Feld besitzt nur eine Komponente in Umfangsrichtung Φ, daher verschwindet der weitere Term der Gleichung. Die L¨osung des magnetischen Feldes erfolgt f¨ ur das Leiterinnere. Das Vorgehen zur L¨osung entspricht dem Vorgehen nach Kap. 3. 7.1.1 Schwache Formulierung der Differenzialgleichung Durch Umstellen von Gl. (7.1) folgt ∂rH Φ ∂r = r J. Die Gleichung ist nach H Φ zu l¨osen. Unter Verwendung der geeigneten Basisfunktion H Φ (r) = N ∑ n =1 α n r n −1 (7.2) <?page no="72"?> 58 L¨osung physik. Bsp. DGL 1’ter Ordnung mit Galerkin-Methode Abbildung 7.2: Magnetischer Feldverlauf im und außerhalb des Leiterstabs folgt N ∑ n =1 α n d dr ( r · r n −1 ) = r J. Unter Einbezug der Wichtungsfunktion w m folgt die Darstellung mit Hilfe des inneren Produkts N ∑ n =1 α n 〈 w m , d dr ( r · r n −1 ) 〉 = 〈 w m , r J 〉 und aus dem Galerkin-Ansatz w m = r m −1 folgt die schwache Formulierung des Durchflutungsgesetzes N ∑ n =1 α n 〈 r m −1 , d dr ( r · r n −1 ) 〉 = 〈 r m −1 , r J 〉 in der Summenschreibweise mit Anzahl N noch zu bestimmende Summanden. <?page no="73"?> 7.1 Durchflutungsgesetz gel¨ost mit nichtlinearer Basis- und Wichtungsfunktion 59 7.1.2 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Die schwache Form des Durchflutungsgesetzes wird mittels Integral ¨ uber dem Raum Ω N ∑ n =1 α n ∫ Ω r m −1 · d dr ( r · r n −1 ) dr ︸ ︷︷ ︸ T erm 1 = ∫ Ω r m −1 · r · J dr ︸ ︷︷ ︸ T erm 2 . formuliert. Dabei ist Ω = r ∈ [0, R], mit Leiterradius R. Die Entwicklung des ersten Terms erfolgt mittels partieller Integration zu ∫ R 0 r m −1 · d dr ( r · r n −1 ) dr = [ r m + n −1 ] R 0 − ∫ R 0 (m − 1) r m −2 r m dr = R m + n −1 − m − 1 m + n − 1 R m + n −1 = n m + n − 1 R m + n −1 . Die Entwicklung des zweiten Terms erfolgt durch Integration zu ∫ R 0 r m −1 · r · J dr = J m + 1 R m +1 . Die schwache Formulierung wird damit erneut als Summe N ∑ n =1 α n n m + n − 1 R m + n −1 ︸ ︷︷ ︸ [ l mn ] = J m + 1 R m +1 ︸ ︷︷ ︸ [ g m ] dargestellt. 7.1.3 L¨osung des linearen Gleichungssystems Die Summenschreibweise der schwachen Form wird in die Matrixschreibweise f¨ ur N = 3 [α] · [l mn ] = [g m ] <?page no="74"?> 60 L¨osung physik. Bsp. DGL 1’ter Ordnung mit Galerkin-Methode [α] · ⎛ ⎜ ⎝ R R 2 R 3 R 2 2 2 R 2 3 3 R 4 4 R 3 3 R 4 2 3 R 5 5 ⎞ ⎟ ⎠ = [g m ] [α] = [l mn ] −1 · [g m ] = ⎛ ⎜ ⎝ 9 R −36 R 2 30 R 3 −18 R 2 96 R 3 −90 R 4 10 R 3 −60 R 4 60 R 5 ⎞ ⎟ ⎠ · ⎛ ⎜ ⎝ J R 2 2 J R 3 3 J R 4 4 ⎞ ⎟ ⎠ ⎛ ⎜ ⎝ α 1 α 2 α 3 ⎞ ⎟ ⎠ = ⎛ ⎜ ⎝ 0 J 2 0 ⎞ ⎟ ⎠ ¨ uberf¨ uhrt. Durch Einsetzen in die Basisfunktion Gl. (7.2) folgt H Φ (r) = 0 · r 0 + J 2 · r 1 + 0 · r 2 = J 2 · r = 0, 509 · 10 6 A 2 m 2 · 2, 5 · 10 −3 m = 636, 25 A/ m das Ergebnis der Feldst¨arke, was dem analytischen Ergebnis aus Abb. 7.1 an der Leiteroberfl¨ache entspricht. 7.2 Gegen¨ uberstellung FEMmit Galerkin-Ergebnis Die Ergebnisgegen¨ uberstellung findet anhand eines gew¨ahlten Leiters statt. F¨ ur den Ergebnisvergleich zwischen numerischer Berechnung nach Galerkin und der FEM-Software COMSOL Multiphysics wurde in Tab. 7.1 die Nr. 2 als Referenz herangezogen. In Tab. 7.2 sind beide Ergebnisse gegen¨ ubergestellt. Eine Abweichung zwischen den beiden Ergebnissen stellt sich im Rahmen der numerischen Genauigkeit und der ber¨ ucksichtigten Nachkommastellen ein. Nachfolgend sind die FEM-Ergebnisse mit COMSOL Multiphysics ersichtlich. In Abb. 7.3 folgt das dazugeh¨orige magnetische Feld im und außerhalb des Leiters. Das Maximum der magnetischen Feldst¨arke befindet sich immer auf der Leiteroberfl¨ache. Die magnetische Feldst¨arke des Leiters mit dem Durchmesser von 5 mm ist an der Leiteroberfl¨ache mit 625 A/ m abzulesen. <?page no="75"?> 7.2 Gegen¨ uberstellung FEMmit Galerkin-Ergebnis 61 Tabelle 7.1: Simulationsdaten Leiter- Leiter- Strom I Strom- Felddurchmesser d fl¨ache A dichte J st¨arke H o Nr. [mm] [(mm) 2 ] [A] [A/ (mm) 2 ] [A/ m] 1 2 3,14 10 3,183 1592,0 2 5 19,63 10 0,509 636,4 3 8 50,27 10 0,199 398,0 4 11 95,03 10 0,105 289,6 5 14 153,93 10 0,065 227,8 6 17 226,98 10 0,044 188,0 7 20 314,16 10 0,032 160,4 Tabelle 7.2: Ergebnisvergleich der magnetischen Feldst¨arke H o an der Leiteroberfl¨ache Nr. COMSOL Multiphysics Galerkin traditionell 2 636,4 A/ m 636,25 A/ m <?page no="76"?> 62 L¨osung physik. Bsp. DGL 1’ter Ordnung mit Galerkin-Methode 0 5 10 15 20 25 Radius r [mm] 0 200 400 600 800 1000 1200 1400 1600 Magnetische Feldstärke H [A/ m] d Stab = 2 mm d Stab = 5 mm d Stab = 8 mm d Stab = 11 mm d Stab = 14 mm d Stab = 17 mm d Stab = 20 mm X: 10 Y: 160.4 X: 2.5 Y: 636.4 X: 3.99 Y: 398.9 Abbildung 7.3: FEM-Simulationsergebnis der magnetischen Feldst¨arke im und außerhalb des Leiterstabs mit Leiterstabdurchmesser als Scharparameter <?page no="77"?> Kapitel 8 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode Im Fortgang werden ausgew¨ahlte Beispielsgleichungen aus Abb. 2.2 gel¨ost. Dies sind im Einzelnen 1. Elektrostatische Feldgleichung 2. W¨armediffusionsgleichung 3. Felddiffusionsgleichung 8.1 Elektrostatische Feldberechnung Das elektrostatische Feld eines Plattenkondensators nach Abb. 8.2 a) soll berechnet werden. Am Kondensator wird eine Spannung von U c = 100 V angelegt. Das Vorgehen erfolgt gem¨aß der Auflistung aus Kap. 4. Im Anschluss wird das elektrostatische Feld durch Gradientenbildung des Potenzials errechnet. 8.1.1 Schwache Formulierung der Differenzialgleichung Zur Berechnung des elektrostatischen Feldes muss zuerst die Differenzialgleichung in Abb. 2.2 c) Gl. (2.1) nach dem Potenzial gel¨ost werden. Die Differenzialgleichung der Elektrostatik in ihrer starken Form ist der Abb. 2.2 c) Gl. (1) zu entnehmen. Durch Umstellen folgt 63 <?page no="78"?> 64 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode ∇ (ε ∇ ϕ) + ρ = 0 = R. Durch Wichtung mit Funktion W des Residuums und Integration ¨ uber das zu betrachtende Gebiet Ω (Plattenabstand) folgt ∫ Ω W R dx = 0 ∫ Ω W [ ∇ (ε ∇ ϕ) + ρ] dx = 0 ∫ Ω W ( ∂ϕ 2 ∂x 2 + ρ ε ) dx = 0 ∫ Ω W ( ∂ϕ 2 ∂x 2 ) dx + ∫ Ω W ( ρ ε ) dx = 0. Durch partielle Integration des ersten Terms der linken Gleichungsh¨alfte folgt die schwache Formulierung der Differenzialgleichung W dϕ dx − ∫ Ω ( dϕ dx dW dx ) dx + ∫ Ω W ( ρ ε ) dx = 0 W ∇ ϕ − ∫ Ω ( ∇ ϕ ∇ W ) dx + ∫ Ω W ( ρ ε ) dx = 0 (8.1) 8.1.2 Diskretisierung des zu l¨osenden Gebiets Ω Die Diskretisierung des Gebietes Ω erfolgt gem. Abb. 5.2. Die Elementl¨ange ist h = 2 mm. 8.1.3 Wahl der Basis- und Wichtungsfunktion Als Basis- und Wichtungsfunktionen werden gem. Kap. 5.1.3 Dreiecksfunktionen φ(x) gew¨ahlt. 8.1.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) Die Vorgehensweise erfolgt in Anlehnung an Kap. 5.1.4. Das Gebiet Ω wurde zwischenzeitlich in n = 5 Teilgebiete (Elemente) Ω n eingeteilt, der Wichtungs- und Basisfunktion <?page no="79"?> 8.1 Elektrostatische Feldberechnung 65 die Dreiecksfunktion zugewiesen. Beim ersten Term der Gl. (8.1) sind gem. Kap. 1.2.3 f¨ ur das bestimmte Integral als Randbedingung jeweils der linke und rechte Rand einzusetzen. An den Berandungen nehmen die Dreiecksfunktionen und damit verbunden auch der Term den Wert Null an. Die Berechnung reduziert sich auf die inneren Knoten. Durch Summenbildung ¨ uber alle Teilgebiete folgt n ∑ i =1 [ ∫ Ω ( ∇ ϕ ∇ φ(x)) dx − ρ ε ∫ Ω φ(x) dx ] = 0. Unter Einbezug der Ansatzfunktion ϕ h (x) = 2 ∑ i =1 ϕ i φ i (x) = ϕ 1 φ 1 (x) + ϕ 2 φ 2 (x), durch Umstellen und Einsetzen folgt ∫ Ω [ ϕ 1 φ 1 (x) dx φ(x) dx + ϕ 2 φ 2 (x) dx φ(x) dx ] dx = ρ ε ∫ Ω φ(x) dx ︸ ︷︷ ︸ Quellterm . Zur Bestimmung von ϕ 1 wird der Funktion φ(x) die Funktion φ 1 (x) und zur Bestimmung von ϕ 2 die Funktion φ 2 (x) zugeordnet. F¨ ur ϕ 1 am Knoten x 1 folgt ∫ Ω [ ϕ 1 φ 1 (x) dx φ 1 (x) dx + ϕ 2 φ 2 (x) dx φ 1 (x) dx ] dx = ρ ε ∫ Ω φ 1 (x) dx und f¨ ur ϕ 2 am Knoten x 2 folgt ∫ Ω [ ϕ 1 φ 1 (x) dx φ 2 (x) dx + ϕ 2 φ 2 (x) dx φ 2 (x) dx ] dx = ρ ε ∫ Ω φ 2 (x) dx. In Matrizenschreibweise zusammengefasst lautet die Gleichung wie folgt: ∫ Ω ( dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx ) dx ( ϕ 1 ϕ 2 ) = ρ ε ∫ Ω φ(x) dx ( 1 1 ) . Beim rechten Term der Gleichung ist eine Unterscheidung der Funktion φ mittels Indizes nicht erforderlich, da diese keinen Einfluss auf das Integrationsergebnis hat. <?page no="80"?> 66 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode 8.1.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Die Vorgehensweise erfolgt in Anlehnung an Kap. 5.1.5. Mit der schwachen Form der Differenzialgleichung folgten die Knotengleichungen eines Elements • Knotenmatrix des ersten inneren Knotens x i mit φ(x) = φ 1 (x). Hiermit folgt das Intervall [x 0 , x 2 ]: ϕ 1 ∫ Ω ( dφ 1 (x) dx dφ 1 (x) dx ) dx + ϕ 2 ∫ Ω ( dφ 2 (x) dx dφ 1 (x) dx ) dx = ρ ε ∫ Ω φ(x) dx und durch Anpassung des Integrationsintervalls ϕ 1 ∫ x 1 x 0 φ 1 dx φ 1 dx dx + ϕ 2 ∫ x 1 x 0 φ 2 dx φ 1 dx dx+ϕ 1 ∫ x 2 x 1 φ 1 dx φ 1 dx dx + ϕ 2 ∫ x 2 x 1 φ 2 dx φ 1 dx dx = ρ ε ∫ Ω φdx. Unter Einbezug der Gegebenheiten von Kap. 5.1.3 folgt gem¨aß Vorgehen nach Kap. 5.1.5 ϕ 1 ∫ x 1 x 0 1 h 1 h dx + ϕ 2 ∫ x 1 x 0 0 1 h dx+ϕ 1 ∫ x 2 x 1 − 1 h − 1 h dx + ϕ 2 ∫ x 2 x 1 1 h − 1 h dx = ρ ε ∫ Ω φ dx. Im zweiten Term ist die Ableitung f¨ ur φ 2 im Intervall [x 0 , x 1 ] nicht definiert und nimmt daher den Wert Null an. Eine weitere Zusammenfassung f¨ uhrt zu ϕ 1 ∫ x 1 x 0 1 h 2 dx + 0 + ϕ 1 ∫ x 2 x 1 1 h 2 dx + ϕ 2 ∫ x 2 x 1 − 1 h 2 dx = ρ ε ∫ Ω φ dx. Nach gliedweiser Integration und Umstellung folgt die Knotenmatrix f¨ ur den Knoten x 1 ϕ 1 1 h + 0 + ϕ 1 1 h − ϕ 2 1 h = ρ ε h 1 h ( 2 − 1 ) ( ϕ 1 ϕ 2 ) = ρ ε h. <?page no="81"?> 8.1 Elektrostatische Feldberechnung 67 • Knotenmatrix des zweiten inneren Knotens x i +1 mit φ(x) = φ 2 (x). Hiermit folgt das Intervall [x 1 , x 3 ]: ϕ 1 ∫ Ω ( dφ 1 (x) dx dφ 2 (x) dx ) dx + ϕ 2 ∫ Ω ( dφ 2 (x) dx dφ 2 (x) dx ) dx = ρ ε ∫ Ω φ(x) dx und durch Anpassung des Integrationsintervalls ϕ 1 ∫ x 2 x 1 φ 1 dx φ 2 dx dx + ϕ 2 ∫ x 2 x 1 φ 2 dx φ 2 dx dx+ϕ 1 ∫ x 3 x 2 φ 1 dx φ 2 dx dx + ϕ 2 ∫ x 3 x 2 φ 2 dx φ 2 dx dx = ρ ε ∫ Ω φdx. Wie zuvor, folgt unter Einbezug der Gegebenheiten von Kap. 5.1.3 und nach dem Vorgehen in Kap. 5.1.5 ϕ 1 ∫ x 2 x 1 − 1 h 1 h dx + ϕ 2 ∫ x 2 x 1 1 h 1 h dx+ϕ 1 ∫ x 3 x 2 0 − 1 h dx + ϕ 2 ∫ x 3 x 2 − 1 h − 1 h dx = ρ ε ∫ Ω φ dx. Im dritten Term ist die erste Ableitung f¨ ur φ 2 im Intervall [x 0 , x 1 ] nicht definiert und nimmt daher den Wert Null an. Eine weitere Zusammenfassung f¨ uhrt zu ϕ 1 ∫ x 2 x 1 − 1 h 2 dx + ∫ x 2 x 1 1 h 2 dx + 0 + ϕ 2 ∫ x 3 x 2 1 h 2 dx = ρ ε ∫ Ω φ dx. Nach gliedweiser Integration und Umstellung folgt die Knotenmatrix f¨ ur den Knoten x 2 ϕ 1 − 1 h + ϕ 2 1 h + 0 + ϕ 2 1 h = ρ ε h 1 h ( − 1 2 ) ( ϕ 1 ϕ 2 ) = ρ ε h. • Elementmatrix: Die Knotengleichungen der beiden inneren Knoten x i und x i +1 werden als Elementgleichung 1 h ( 2 − 1 − 1 2 ) · ( ϕ 1 ϕ 2 ) = ρ ε h ( 1 1 ) zusammengefasst. <?page no="82"?> 68 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode • Koeffizientenmatrix: Alle Knotengleichungen werden in der Koeffizientenmatrix zu einem linearen Gleichungssystem (Koeffizientenmatrizengleichung) ⎛ ⎜ ⎜ ⎜ ⎝ 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 ⎞ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎝ ϕ 1 ϕ 2 ϕ 3 ϕ 4 ⎞ ⎟ ⎟ ⎟ ⎠ = ρ h 2 ε ⎛ ⎜ ⎜ ⎜ ⎝ 1 1 1 1 ⎞ ⎟ ⎟ ⎟ ⎠ zusammengef¨ uhrt. • Gem¨aß Aufgabenstellung liegt am Kondensator eine Spannung an, welche durch die Dirichlet-Randbedingungen ber¨ ucksichtigt wird (elliptische Differenzialgleichung). Der Quellterm ist deshalb mit ρ = 0 gleich null zu setzen, um eine ¨ Uberbestimmtheit zu vermeiden. Unter Einbezug der Dirichlet-Randbedingungen (ϕ(x 0 ) = 0 V, ϕ(x 5 ) = 100 V) folgt das lineare Gleichungssystem mit ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ϕ 0 ϕ 1 ϕ 2 ϕ 3 ϕ 4 ϕ 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 0 0 0 0 0 100 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . 8.1.6 L¨osung des linearen Gleichungssystems Die L¨osung des linearen Gleichungssystems erfolgt gem. Kap. 5.1.6. Der Ergebnisvektor des Potenzials f¨ ur alle Knoten ist ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ϕ 0 ϕ 1 ϕ 2 ϕ 3 ϕ 4 ϕ 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 0 20 40 60 80 100 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ und den inneren Knoten <?page no="83"?> 8.1 Elektrostatische Feldberechnung 69 Tabelle 8.1: Numerische Ergebnisse des Potenzialverlaufs Knotennummer x i , i = 0 1 2 3 4 5 Position [mm] 0 2 4 6 8 10 ϕ h (x i ) [V ] 0 20 40 60 80 100 ϕ h = 4 ∑ i =1 ϕ i (x) φ i (x) = 20 φ 1 + 40 φ 2 + 60 φ 3 + 80 φ 4 . In Tab. 8.1 sind die Ergebnisse zusammengefasst. In Abb. 8.1 ist die grafische Darstellung der L¨osung einschließlich der vier Basisfunktionen und Dirichlet-Randbedingungen ersichtlich. Abbildung 8.1: Ergebnisdarstellung des Potenzialverlaufs mittels Galerkin-Methode In Abb. 8.2 a) ist der Plattenkondensator mit dem elektrostatischen Feld und den Potenzialen der inneren Knoten ϕ 1 bis ϕ 4 ersichtlich. Die Dirichlet-Randbedingungen wurden den ¨außeren Knoten ϕ 0 = 0 V und ϕ 5 = 100 V auferlegt. In Abb. 8.2 b) ist <?page no="84"?> 70 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode Abbildung 8.2: Plattenkondensator und Potenzialverlauf das Ergebnis des Potenzialverlaufs ¨ uber dem Plattenabstand (x-Achse) dargestellt. Aus dem Potenzialverlauf der Abb. 8.2 b) wird das elektrostatische Feld mit E = − gradϕ = − dϕ dx n = − 20 V 2 mm = − 10 kV / m berechnet. Zur Probe erfolgt die Berechnung der Spannung U c ¨ uber dem Kondensator durch Integration entlang des Plattenabstandes mit ∫ x 5 x 0 E d l = 10 kV / m · 10 mm = 100 V. 8.2 Ortsabh¨angige Temperaturberechnung Der W¨armedurchgang durch K¨orper wird mittels der W¨armediffusionsgleichung beschrieben. Diese ist durch zwei Zeit- und zwei Ortsableitungen gekennzeichnet. Zur anschaulichen Deutung der eindimensionalen W¨armediffusion dient Abb. 8.3. Ein K¨orper mit der W¨armeleitf¨ahigkeit Λ aus Kupfer wird an der Stirnfl¨ache einseitig auf 100 ◦ C erw¨armt. Der W¨armestrom breitet sich dabei nur in x-Richtung aus. Die Temperatur nimmt entlang zunehmender x-Achse bei einem betrachteten Zeitpunkt t als Funktion des Ortes ab. Die Temperatur erh¨oht sich entlang zunehmender x-Achse als Funktion <?page no="85"?> 8.2 Ortsabh¨angige Temperaturberechnung 71 der Zeit. Gesucht wird, zu einem Zeitpunkt t, die ¨ortliche Temperaturverteilung im K¨orper, die mittels den gr¨ un dargestellten Pfeilen im Balken symbolisiert wird. Abbildung 8.3: Beispiel eines eindimensionalen W¨armediffusionsvorgangs 8.2.1 Schwache Formulierung der Differenzialgleichung Die eindimensionale W¨armediffusionsgleichung nach Abb. 2.2, b2) wird in die Poisson’sche Differenzialgleichung der Form d 2 υ(x) dx 2 = ρ c λ dυ dt ︸ ︷︷ ︸ K (8.2) d 2 υ(x) dx 2 = K, x ∈ Ω mit den Randbedingungen υ(x 0 ) und υ(x 5 ) ¨ uberf¨ uhrt. Die L¨osung der Differenzialgleichung erfolgt f¨ ur ein angenommenes dυ/ dt in der N¨ahe eines angenommenen Zeitschritts t, gefolgt von der Bildung des inneren Produkts durch Integration des gewichteten Residuums ¨ uber das Gebiet Ω <?page no="86"?> 72 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode ∫ Ω R W dx = 0 ∫ Ω ( d 2 υ(x) dx 2 − K ) W dx = 0 ∫ Ω d 2 υ(x) dx 2 W dΩ − K ∫ Ω W dx = 0. Nach erfolgter partieller Integration des ersten Terms folgt die schwache Form der Diffusionsgleichung W dυ(x) dx − ∫ Ω dυ(x) dx dW dx dx − K ∫ Ω W dx = 0 W ∇ υ(x) − ∫ Ω ∇ υ(x) ∇ W dx − K ∫ Ω W dx = 0. (8.3) 8.2.2 Diskretisierung des zu l¨osenden Gebiets Ω Die Diskretisierung des Gebietes Ω erfolgt gem. Abb. 5.2. Die Elementl¨ange betr¨agt h = 10 mm. 8.2.3 Wahl der Basis- und Wichtungsfunktion Als Basis- und Wichtungsfunktionen werden gem. Kap. 5.1.3 Dreiecksfunktionen definiert. 8.2.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) Die Vorgehensweise erfolgt in Anlehnung an Kap. 5.1.4. Das Gebiet Ω wurde zwischenzeitlich in n = 5 Teilgebiete (Elemente) Ω n eingeteilt. Die Wichtungsfunktion und die Basisfunktion wurden mit der Dreiecksfunktion φ(x) gleichgesetzt. Aufgrund der Randbedingungen wird der erste Term von Gl. (8.3) gleich Null gesetzt, da die Dreiecksfunktionen am Rande gleich null sind. Damit reduziert sich die Berechnung auf die inneren Knoten. Unter Einbezug der Ansatzfunktion υ h (x) = 2 ∑ i υ i φ i (x) = υ 1 φ 1 (x) + υ 2 φ 2 (x) <?page no="87"?> 8.2 Ortsabh¨angige Temperaturberechnung 73 und Einsetzen in Gl. (8.3) mit anschließendem Umformen folgt ∫ Ω [ υ 1 dφ 1 dx dφ dx + υ 2 dφ 2 dx dφ dx ] dx = K ∫ Ω φ(x) dx ︸ ︷︷ ︸ Quellterm . Zur Bestimmung der Temperaturen υ 1 und υ 2 an den Knoten x 1 und x 2 kann erneut die Matrizenschreibweise ∫ Ω ( dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx ) dx ( υ 1 υ 2 ) = K ∫ Ω φ(x) dx ( 1 1 ) angewendet werden. 8.2.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Die Vorgehensweise erfolgt in Anlehnung an Kap. 5.1.5. • Knotenmatrix des ersten inneren Knotens x i mit φ(x) = φ 1 (x). Hiermit folgt das Intervall [x 0 , x 2 ]: 1 h ( 2 − 1 ) · ( υ 1 υ 2 ) = K h • Knotenmatrix des zweiten inneren Knotens x i +1 mit φ(x) = φ 2 (x). Hiermit folgt das Intervall [x 1 , x 3 ]: 1 h ( − 1 2 ) · ( υ 1 υ 2 ) = K h • Elementmatrix: Die beiden inneren Knoten x i und x i +1 werden als Elementgleichung zusammengefasst 1 h ( 2 − 1 − 1 2 ) · ( υ 1 υ 2 ) = K h ( 1 1 ) . <?page no="88"?> 74 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode Tabelle 8.2: Werkstoffangaben/ Koeffizienten/ Randbedingungen Angaben f¨ ur Werkstoff Kupfer Dichte ρ, [kg/ m 3 ] 8933 Spez. W¨armekapazit¨at c, [J/ (kgK)] 383 W¨armeleitf¨ahigkeit λ, [W/ (mK)] 384 ρ c/ λ, [s/ m 2 ] 8937 dυ/ dt, [K/ s]; (48,41-51,61) ◦ C/ 0,5 s, bei t= 4 s -6,4 K, [K/ m 2 ] -57 196,8 h, [m] 0,01 K h 2 , [K] -5,72 Dirichlet-Randbedingungen υ(x 0 ), [ ◦ C] 100 υ(x 5 ), [ ◦ C] 20 • Koeffizientenmatrix: Die beiden Knotengleichungen werden in der Koeffizientenmatrix zu einem linearen Gleichungssystem zusammengef¨ uhrt ⎛ ⎜ ⎜ ⎜ ⎝ 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ S · ⎛ ⎜ ⎜ ⎜ ⎝ υ 1 υ 2 υ 3 υ 4 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ υ h = K h 2 ⎛ ⎜ ⎜ ⎜ ⎝ 1 1 1 1 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ f . • Unter Einbezug der Dirichlet-Randbedingungen nach Tab. 8.2 folgt erneut das Gleichungssystem f¨ ur die inneren und ¨außeren Knoten ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ϕ 0 ϕ 1 ϕ 2 ϕ 3 ϕ 4 ϕ 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 100 − 5, 72 − 5, 72 − 5, 72 − 5, 72 20 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . <?page no="89"?> 8.2 Ortsabh¨angige Temperaturberechnung 75 8.2.6 L¨osung des linearen Gleichungssystems Die L¨osung des linearen Gleichungssystems erfolgt gem. Kap. 5.1.6 mit S υ = f S −1 (S υ h ) = S −1 f ( S −1 S ) ︸ ︷︷ ︸ E υ h = S −1 f υ h = S −1 f = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 100, 0 72, 6 50, 8 34, 8 24, 6 20, 0 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . 4 Zeit t [s] 2 0 20 40 (Linker Rand) Weg x [m] (rechter Rand) 60 [°C] 80 0 100 0 0.005 0.01 120 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0.05 Abbildung 8.4: MATLAB-Ergebnis eines eindimensionalen W¨armediffusionsvorgangs F¨ ur die inneren Knoten folgt <?page no="90"?> 76 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode υ h = 4 ∑ i =1 υ i (x) φ i (x) = 72, 6 φ 1 + 50, 8 φ 2 + 34, 8 φ 3 + 24, 6 φ 4 . In Tab. 8.2 sind die erforderlichen Angaben zur Berechnung angegeben. Die Werkstoffdaten wurden den Tabellen aus [18], die Angaben zur Berechnung des Koeffizienten K aus Abb. 8.5 entnommen. In Abb. 8.4 ist das MATLAB-Ergebnis einer Temperaturverteilung ¨ uber Ort und Zeit dargestellt, erstellt mit der PDE-Toolbox. Der dazugeh¨orige MATLAB-Code kann dem Anhang A.1 entnommen werden. Der ¨ Ubergang zur Darstellung der Temperatur ¨ uber dem Ort wurde in Abb. 8.5 vollzogen. Diese Darstellung erm¨oglicht eine Gegen¨ uberstellung mit den oben erzielten Ergebnissen. Die Abweichungen sind auf Rundungen und Ablesegenauigkeit (MATLAB-Data-Cursor steht nicht genau auf den Knoten) zur¨ uckzuf¨ uhren. Zudem wurde dυ/ dt bei einem mittleren Weg (19,8 mm) abgelesen. In Tab. 8.3 sind die Ergebnisse einander gegen¨ ubergestellt. 0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0.05 (Linker Rand) Weg x [m] (rechter Rand) 0 20 40 60 80 100 120 [°C] 4,0 s 3,5 s 3,0 s 2,5 s 2,0 s 1,5 s 1,0 s 0,5 s 0,0 s X: 0.01007 Y: 73.93 X: 0.0198 Y: 51.61 X: 0.0198 Y: 48.41 X: 0.02987 Y: 33.81 X: 0.03993 Y: 22.76 X: 0 Y: 100 Abbildung 8.5: MATLAB-Ergebnis eines eindimensionalen W¨armediffusionsvorgangs mit der Zeit t als Scharparameter <?page no="91"?> 8.2 Ortsabh¨angige Temperaturberechnung 77 Abbildung 8.6: Grafische Ergebnisdarstellung der ¨ortlichen Temperaturverteilung In Abb. 8.6 sind schematisch die Basisfunktionen einschl. deren Ansatzfunktionen der inneren Knoten ¨ uber dem Ort (Ω) dargestellt. Die Dirichlet-Randbedingungen werden durch die ¨außeren Knoten (υ 0 , υ 5 ) verk¨orpert. Eine Gegen¨ uberstellung der mit MATLAB- und mit Galerkin-Methode erzielten Ergebnisse ist der Abb. 8.7 zu entnehmen. Die dazu erforderlichen Wertepaare wurden aus Tab. 8.3 ¨ ubernommen. Tabelle 8.3: Gegen¨ uberstellung der Ergebnisse des Temperatur- und Ortsverlaufs m. H. Abb. 8.5 Knotennummer x i , i = 0 1 2 3 4 5 Position [mm] 0 10 20 30 40 50 υ h [ ◦ C] bei t = 4 s; Galerkin 100 72,6 50,8 34,8 24,6 20 υ [ ◦ C] bei t = 4 s; MATLAB 100 73,9 51,6 33,8 22,8 20 υ h [ ◦ C] bei t = ∞ s; Galerkin 100 84 68 52 36 20 <?page no="92"?> 78 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode 8.2.7 Diffusionsvorgang vollendet Der Diffusionsvorgang wird als vollendet betrachtet, wenn die zeitliche Temperatur¨anderung des rechten Terms der Gl. (8.2) f¨ ur t → ∞ 0 10 20 30 40 50 x [mm] 20 40 60 80 100 [°C] Galerkin, t = 4 s MATLAB PDE, t = 4 s Galerkin, t = s Abbildung 8.7: Ergebnisgegen¨ uberstellung ρ c λ dυ dt → 0 gegen den Wert Null strebt und damit ein station¨arer Zustand eintritt. Die Poisson’sche Differenzialgleichung wird damit in Laplace’sche Differenzialgleichung ¨ uberf¨ uhrt. Das lineare Gleichungssystem wird zu ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 − 1 2 − 1 0 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ υ 0 υ 1 υ 2 υ 3 υ 4 υ 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 100 0 0 0 0 20 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . Die L¨osung ergibt die ¨ortliche Temperaturverteilung des Beharrungszustandes <?page no="93"?> 8.3 Ortsabh¨angige Magnetfeldberechnung 79 ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ υ 0 υ 1 υ 2 υ 3 υ 4 υ 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 100 84 68 52 36 20 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . Der Beharrungszustand ist in Abb. 8.7 und in Tab. 8.3 dokumentiert. 8.3 Ortsabh¨angige Magnetfeldberechnung In Abb. 8.8 a) ist der Schnitt eines Topfanker-Magnetkreises, bestehend aus dem Anker, Joch, R¨ uckstellfeder, Federvorspannh¨ ulse, Wicklungstr¨ager und Wicklung, ersichtlich. W¨ahrend dem Einschaltvorgang durchdringen die geschlossenen B-Feldlinien die Polfl¨achen von innen nach außen. Abbildung 8.8: Beispiel einer eindimensionalen Magnetfelddiffusion Die weiterf¨ uhrende Betrachtung wird auf den Ausschnitt im Innenpol der Abb. 8.8 b) beschr¨ankt. In aller Regel finden bei Elektromagneten ferromagnetische Werkstoffe mit nichtlinearer B(H)-Kennlinie Anwendung. Im Fortgang wird dagegen ein Werkstoffmit linearer B(H)-Kennlinie angenommen. Als K¨orper-Werkstoffwurde der linearen paramagnetischen Werkstoff Kupfer gew¨ahlt, um eine konstante Permeabilit¨at zu erhalten. <?page no="94"?> 80 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode Am linken Rand des Polausschnitts wird eine hohe Flussdichte eingepr¨agt. Die Intensit¨at der Flussdichte nimmt in radialer Ausdehnung ab und erreicht am Innenradius des Innenpols (rechter Rand) ein Minimum. Die Flussdichte breitet sich damit normal zur Flussrichtung aus (transversale Ausbreitung). Dieser Diffusionsvorgang unterliegt zudem einer zeitlichen ¨ Anderung. 8.3.1 Schwache Formulierung der Differenzialgleichung Die eindimensionale Felddiffusionsgleichung ist in 2.2 b) Gl. (1) entnommen. Da im Fortgang nur eine Dimension des Vektors B betrachtet wird, ist der Vektor gleich dessen Betrag B. Die eindimensionale Felddiffusionsgleichung wird in die Form d 2 B(x) dx 2 = μ 0 κ dB dt ︸ ︷︷ ︸ K d 2 B dx 2 = K, x ∈ Ω ¨ uberf¨ uhrt und ihr die Randbedingungen B(x 0 ), B(x 5 ) zugeordnet. Die L¨osung erfolgt f¨ ur ein angenommenes B/ dt. Es folgt die Integration ¨ uber das mit W gewichtete Residuum ∫ Ω R W dx = 0 ∫ Ω ( d 2 B(x) dx 2 − K ) W dx = 0 ∫ Ω d 2 B(x) dx 2 W dx − K ∫ Ω W dx = 0. Nach erfolgter partieller Integration des ersten Terms folgt die schwache Form der Felddiffusionsgleichung W dB(x) dx − ∫ Ω dB(x) dx W dx dx − K ∫ Ω W dx = 0 W ∇ B(x) − ∫ Ω ∇ B(x) ∇ W dx − K ∫ Ω W dx = 0. (8.4) 8.3.2 Diskretisierung des zu l¨osenden Gebiets Ω Die Diskretisierung des Gebietes Ω erfolgt gem. Abb. 5.2. Die Elementl¨ange betr¨agt h = 2 mm. <?page no="95"?> 8.3 Ortsabh¨angige Magnetfeldberechnung 81 8.3.3 Wahl der Basis- und Wichtungsfunktion Die Basisfunktionen werden gem. Kap. 5.1.3 definiert. 8.3.4 Formulierung der schwachen Form mit Dreiecksfunktionen φ(x) Die Vorgehensweise erfolgt in Anlehnung an Kap. 5.1.4. Das Gebiet Ω wurde zwischenzeitlich in n = 5 Teilgebiete (Elemente) Ω n eingeteilt und die Wichtungsfunktion W zusammen mit der Basisfunktion durch die Dreiecksfunktion φ(x) ersetzt. Aufgrund der Randbedingungen wird der erste Term von Gl. (8.4) gleich Null gesetzt, da die Basisfunktionen am Rande gleich null sind. Damit reduziert sich die Berechnung auf die inneren Knoten. Unter Einbezug der Ansatzfunktion B h (x) = 2 ∑ i =1 B i φ i (x) = B 1 φ 1 (x) + B 2 φ 2 (x) folgt durch Einsetzten und Umstellen ∫ Ω [ B 1 dφ 1 dx dφ dx + B 2 dφ 2 dx dφ dx ] dx = K ∫ Ω φ i (x) dx ︸ ︷︷ ︸ Quellterm . Zur Bestimmung der Flussdichte B 1 und B 2 an den Knoten x 1 und x 2 kann erneut die Matrizenschreibweise ∫ Ω ( dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 1 ( x ) dx dφ 1 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx dφ 2 ( x ) dx ) dx ( B 1 B 2 ) = K ∫ Ω φ(x) dx ( 1 1 ) angewendet werden. 8.3.5 ¨ Uberf¨ uhrung des Gleichungssystems in eine Matrizengleichung Die Vorgehensweise erfolgt in Anlehnung an Kap. 5.1.5. <?page no="96"?> 82 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode • Knotenmatrix des ersten inneren Knotens x i mit φ(x) = φ 1 (x). Hiermit folgt das Intervall [x 0 , x 2 ]: 1 h ( 2 − 1 ) · ( B 1 B 2 ) = K h • Knotenmatrix des zweiten inneren Knotens x i +1 mit φ(x) = φ 2 (x). Hiermit folgt das Intervall [x 1 , x 3 ]: 1 h ( − 1 2 ) · ( B 1 B 2 ) = K h • Elementmatrix: Die beiden inneren Knoten x i und x i +1 werden als Elementgleichung zusammengefasst 1 h ( 2 − 1 − 1 2 ) · ( B 1 B 2 ) = K h. • Koeffizientenmatrix: Die beiden Knotengleichungen werden in der Koeffizientenmatrix zu einem linearen Gleichungssystem zusammengef¨ uhrt 1 K h 2 ⎛ ⎜ ⎜ ⎜ ⎝ 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 − 1 0 0 − 1 2 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ S · ⎛ ⎜ ⎜ ⎜ ⎝ B 1 B 2 B 3 B 4 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ B h = ⎛ ⎜ ⎜ ⎜ ⎝ 1 1 1 1 ⎞ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ f . • Unter Einbezug der Dirichlet-Randbedingungen ( B(x 0 ) = 1 T, B(x 5 ) = 0,2 T) folgt erneut das Gleichungssystem ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0 0 17, 7 − 35, 4 17, 7 0 0 0 0 17, 7 − 35, 4 17, 7 0 0 0 0 17, 7 − 35, 4 17, 7 0 0 0 0 17, 7 − 35, 4 17, 7 0 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ B 0 B 1 B 2 B 3 B 4 B 5 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 1 1 1 1 0, 2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . <?page no="97"?> 8.3 Ortsabh¨angige Magnetfeldberechnung 83 Tabelle 8.4: Werkstoffangaben/ Koeffizienten/ Randbedingungen Angaben f¨ ur Werkstoff Kupfer Permeabilit¨at μ 0 , [V s/ (Am)] 4 π 10 −7 = 1,2 10 −6 Spez. elektr. Leitf¨ahigkeit κ, [A/ (V m)] 56,2 10 6 dB/ dt, [T / s]; (0,49-0,53) T/ 0,2 ms -200 K, [V s/ m 4 ] -14 124,6 h, [m] 0,002 K h 2 , [V s/ m 2 ] -0,0565 1/ (K h 2 ), [m 2 / (V s)] -17,7 Dirichlet-Randbedingungen B(x 0 ), [T ] 1 B(x 5 ), [T ] 0,2 8.3.6 L¨osung des linearen Gleichungssystems Die L¨osung des linearen Gleichungssystems erfolgt gem. Kap. 5.1.6 mit S B h = f S −1 (S B h ) = S −1 f ( S −1 S ) ︸ ︷︷ ︸ E B h = S −1 f B h = S −1 f = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1, 00 0, 73 0, 51 0, 35 0, 25 0, 20 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . F¨ ur die inneren Knoten folgt B h = 4 ∑ i =1 B i (x) φ i (x) = 0, 73 φ 1 + 0, 51 φ 2 + 0, 35 φ 3 + 0, 25 φ 4 . <?page no="98"?> 84 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode In Tab. 8.4 sind die erforderlichen Angaben zur Berechnung angegeben. Die Werkstoffdaten wurden den Tabellen aus [18], die Angaben zur Berechnung des Koeffizienten K aus Abb. 8.10 entnommen. 2 Zeit t [s] 10 -3 1 0 0.2 0.4 (Linker Rand) Weg x [m] (rechter Rand) 0.6 B [T] 0.8 0 1 0 0.001 0.002 1.2 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01 Abbildung 8.9: MATLAB-Ergebnis eines eindimensionalen Magnetfelddiffusionsvorgangs In Abb. 8.9 ist das MATLAB-Ergebnis einer Flussdichteverteilung ¨ uber Ort und Zeit dargestellt, erstellt mit der PDE-Toolbox. Der dazugeh¨orige MATLAB-Code ist dem Anhang A.2 zu entnehmen. Das Ergebnis wurde mit COMSOL Multiphysics best¨atigt. Der ¨ Ubergang in die Darstellung der Flussdichte ¨ uber dem Ort mit der Zeit t als Scharparameter wurde in Abb. 8.10 vollzogen. Diese Darstellung erm¨oglicht eine Gegen¨ uberstellung mit den oben erzielten Ergebnissen. Die Abweichungen sind auf Rundungen und Ableseungenauigkeit (MATLAB-Data-Cursor steht nicht genau auf den Knoten) zur¨ uckzuf¨ uhren. Zudem wurde dB/ dt bei einem mittleren Weg (4,08 mm) abgelesen. In Tab. 8.5 sind die Ergebnisse einander gegen¨ ubergestellt. In Abb. 8.11 sind schematisch die Basisfunktionen einschl. deren Ansatzfunktionen der inneren Knoten ¨ uber dem Ort (Ω) dargestellt. Die Dirichlet-Randbedingungen werden durch die ¨außeren Knoten (B 0 , <?page no="99"?> 8.3 Ortsabh¨angige Magnetfeldberechnung 85 B 5 ) verk¨orpert. Eine Gegen¨ uberstellung der MATLAB- und mit Galerkin-Methode erzielten Ergebnisse ist der Abb. 8.12 zu entnehmen. Die dazu erforderlichen Wertepaare wurden aus Tab. 8.5 ¨ ubernommen. Im Anhang A.3 wurden die mit der MATLAB- PDE-Toolbox erzielten Ergebnisse mit den Ergebnissen von COMSOL-Multiphysics gegen¨ ubergestellt. 0 0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01 (Linker Rand) Weg x [m] (rechter Rand) 0 0.2 0.4 0.6 0.8 1 1.2 B [T] 1,4 ms 1,2 ms 1,0 ms 0,8 ms 0,6 ms 0,4 ms 0,2 ms 0,0 ms X: 0.002041 Y: 0.7496 X: 0.004082 Y: 0.5281 X: 0.004082 Y: 0.49 X: 0.006122 Y: 0.3583 X: 0.007959 Y: 0.2619 X: 0 Y: 1 Abbildung 8.10: MATLAB-Ergebnis eines eindimensionalen Felddiffusionsvorgangs mit der Zeit t als Scharparameter <?page no="100"?> 86 L¨osung physik. Bsp. DGL 2’ter Ordnung mit Galerkin-Methode Tabelle 8.5: Gegen¨ uberstellung der Ergebnisse des Flussdichte- und Ortsverlaufs m. H. Abb. 8.10 Knotennummer x i , i = 0 1 2 3 4 5 Position [mm] 0 2 4 6 8 10 B(x i ) [T ] bei t = 1,4 ms; Galerkin 1 0,73 0,51 0,35 0,25 0,2 B(x i ) [T ] bei t = 1,4 ms; MATLAB 1 0,75 0,53 0,36 0,26 0,2 Abbildung 8.11: Grafische Ergebnisdarstellung der ¨ortlichen Flussdichteverteilung <?page no="101"?> 8.3 Ortsabh¨angige Magnetfeldberechnung 87 0 0.002 0.004 0.006 0.008 0.01 x [mm] 0.2 0.4 0.6 0.8 1 B [T] Galerkin MATLAB PDE Abbildung 8.12: Ergebnisgegen¨ uberstellung mit Wertepaare aus Tab. 8.5 <?page no="102"?> Kapitel 9 Einf¨ uhrung in die Finite-Differenzen-Methode Diese Methode sieht das Ersetzen des Differenzialquotienten durch einen Differenzenquotienten vor. Daraus folgt: Diff erenzialquotient = Diff erenzenquotient + Diskretisationsf ehler. Mit dieser Diskretisation wird eine Gleichung in eine algebraische Gleichung ¨ uberf¨ uhrt (algebraisiert). In Abb. 2.2 b) ist mit Gl. (1) die Felddiffusionsgleichung in Deltaoperator- Schreibweise ersichtlich. Die Differentialoperator-Schreibweise lautet: ∂ 2 B ∂x 2 = μ 0 κ ∂ B ∂t . (9.1) Die L¨osung erfolgt mittels der Finite-Differenzen-Methode (FDM) mit impliziter und expliziter Methode. Als empfehlenswerte Literatur sei hier [22], Kap. 3: ”Finite Difference Methods“ mit Aufgaben und L¨osungen genannt. 9.1 Numerische Notation der linearen Felddiffusionsgleichung Gl. (9.1) wird als eindimensionale Felddiffusionsgleichung mittels numerischer Methoden gel¨ost. Im weiteren Verlauf werden die Ableitungen mit Hilfe von Differenzenquotienten ausgedr¨ uckt. Hierzu wird Gl. (9.2) als Vorw¨artsdifferenzenquotient und Gl. (9.3) als Zentraldifferenzenquotient geschrieben ([22], S. 126 f.) 88 <?page no="103"?> 9.2 L¨osung mit impliziter Methode nach Crank-Nicolson 89 ∂B ∂t = B j +1 − B j Δt (9.2) ∂ 2 B ∂x 2 = B i +1 − 2B i + B i −1 (Δx) 2 . (9.3) Werden die Gleichungen (9.2) und (9.3) in Gl. (9.1) eingesetzt, so folgt B i +1 ,j − 2B i,j + B i −1 ,j (Δx) 2 = μκ B i,j +1 − B i,j Δt . (9.4) Die L¨osung von Gl. (9.4) erfolgt mittels impliziter und expliziter Methode. 9.2 L¨osung mit impliziter Methode nach Crank- Nicolson In diesem Anwendungsfall muss die Zeit und die Strecke in kleinere Einheiten unterteilt (diskretisiert) werden. Dies macht den Einsatz der Finiten-Differenzen-Methode vorteilhaft. Die Differenzenquotienten der Diffusionsgleichung (Differenzialgleichung 2’ter Ordnung) werden durch Differenzenquotienten ersetzt. Im weiteren Verlauf wird die implizite Methode nach Crank-Nicolson angewendet. Es folgt die ¨ Uberf¨ uhrung der Gleichung in ein lineares (n,n)-Gleichungssystem mit anschließender L¨osung, gefolgt von einem Anwendungsbeispiel. 9.2.1 ¨ Uberf¨ uhrung der Diffusionsgleichung in eine Matrizengleichung Hierzu wird der linke Term der Gl. (9.4) durch den Mittelwert des zentralen Differenzenquotienten der j’ten und j + 1’ten Zeitreihe ersetzt 1 2 ( B i +1 ,j − 2B i,j + B i −1 ,j (Δx) 2 + B i +1 ,j +1 − 2B i,j +1 + B i −1 ,j +1 (Δx) 2 ) = μ κ B i,j +1 − B i,j Δt . Durch Umstellen folgt (B i +1 ,j − 2B i,j + B i −1 ,j + B i +1 ,j +1 − 2B i,j +1 + B i −1 ,j +1 ) Δt 2(Δx) 2 μκ = B i,j +1 − B i,j . <?page no="104"?> 90 Einf¨ uhrung in die Finite-Differenzen-Methode Mit der Substitution k = 1 2μ κ Δt (Δx) 2 (9.5) wird die Lesbarkeit erleichtert kB i +1 ,j − 2kB i,j + kB i −1 ,j + kB i +1 ,j +1 − 2kB i,j +1 + kB i −1 ,j +1 = B i,j +1 − B i,j . Durch Umsortieren und Trennung der einzelnen Terme nach den j’ten und j+1’ten Zeitschritt folgt kB i −1 ,j − 2kB i,j + B i,j + kB i +1 ,j = − kB i −1 ,j +1 + B i,j +1 + 2kB i,j +1 − kB i +1 ,j +1 . Eine anschließende Zusammenfassung erm¨oglicht kB i −1 ,j + (1 − 2k) B i,j + kB i +1 ,j = − kB i −1 ,j +1 + (1 + 2k) B i,j +1 − kB i +1 ,j +1 (9.6) das Umstellen der Gleichung nach dem j’ten und j+1’ten Zeitschritt. Die linke Seite der Gl. (9.6) ist bekannt, da diese den gegenw¨artigen j’ten Zeitschritt und bekannte geometrische Schritte enth¨alt. Dem gegen¨ uber steht die rechte Seite der Gl. (9.6) mit zwar bekannten geometrischen Schritten, aber dem unbekannten Zustand zur Zeit j + 1 (vgl. Abb. 9.1). Mit der Substitution der linken Seite von Gl. (9.6) b 1 = k B i −1 ,j + (1 − 2k) B i,j + k B i +1 ,j = [k (1 − 2k) k] ⎛ ⎜ ⎝ B i −1 ,j B i,j B i +1 ,j ⎞ ⎟ ⎠ (9.7) und Substitution der rechten Seite von Gl. (9.6) folgt jeweils in Matrixschreibweise b 1 = [ − k (1 + 2k) − k] ︸ ︷︷ ︸ A ⎛ ⎜ ⎝ B i −1 ,j +1 B i,j +1 B i +1 ,j +1 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ x . (9.8) <?page no="105"?> 9.2 L¨osung mit impliziter Methode nach Crank-Nicolson 91 Abbildung 9.1: Schritteinteilung der impliziten Methode Die Gleichungen haben damit die Form b = A x, (9.9) wobei A und x Matrizen sind. 9.2.2 L¨osung der Matrizengleichung Eine Voraussetzung f¨ ur die L¨osbarkeit eines inhomogenen linearen (n,n)-Systems b = A x ist die quadratische Matrix, was bei Gl. (9.8) noch nicht der Fall ist. Um diese L¨osungsmethode dennoch anwenden zu k¨onnen, werden die geometrischen Schritte von Gl. (9.7) mit b 2 = [k (1 − 2k) k] ⎛ ⎜ ⎝ B i,j B i +1 ,j B i +2 ,j ⎞ ⎟ ⎠ (9.10) b 3 = [k (1 − 2k) k] ⎛ ⎜ ⎝ B i +1 ,j B i +2 ,j B i +3 ,j ⎞ ⎟ ⎠ (9.11) und die der Gl. (9.8) mit <?page no="106"?> 92 Einf¨ uhrung in die Finite-Differenzen-Methode b 4 = [ − k (1 + 2k) − k] ⎛ ⎜ ⎝ B i,j +1 B i +1 ,j +1 B i +2 ,j +1 ⎞ ⎟ ⎠ (9.12) b 5 = [ − k (1 + 2k) − k] ⎛ ⎜ ⎝ B i +1 ,j +1 B i +2 ,j +1 B i +3 ,j +1 ⎞ ⎟ ⎠ (9.13) Abbildung 9.2: Erweiterte Schritteinteilung der impliziten Methode erweitert. Gl. (9.10) und Gl. (9.11) werden in Gl. (9.7) eingef¨ ugt, b 1 = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ k (1 − 2k) k 0 0 0 k (1 − 2k) k 0 0 0 k (1 − 2k) k 0 0 0 k (1 − 2k) 0 0 0 0 k ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ B i −1 ,j B i,j B i +1 ,j B i +2 ,j B i +3 ,j ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ k B i −1 ,j + (1 − 2k) B i,j + k B i +1 ,j k B i,j + (1 − 2k)B i +1 ,j + k B i +2 ,j k B i +1 ,j + (1 − 2k)B i +2 ,j + k B i +3 ,j k B i +2 ,j + (1 − 2k)B i +3 ,j k B i +3 ,j ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ b <?page no="107"?> 9.2 L¨osung mit impliziter Methode nach Crank-Nicolson 93 alle Elemente enthalten bekannte Gr¨oßen. Gl. (9.12) und Gl. (9.13) werden in Gl. (9.8) eingef¨ ugt. Damit folgt b 1 = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ − k (1 + 2k) − k 0 0 0 − k (1 + 2k) − k 0 0 0 − k (1 + 2k) − k 0 0 0 − k (1 + 2k) 0 0 0 0 − k ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ A · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ B i −1 ,j +1 B i,j +1 B i +1 ,j +1 B i +2 ,j +1 B i +3 ,j +1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ x . Die Abb. 9.2 zeigt die Erweiterung der geometrischen Schritte. Die Matrix A ist nun quadratisch. Damit gen¨ ugt die Gleichung den Anforderungen zur L¨osung linearer Gleichungssysteme mit Hilfe der Cramer’schen Regel. F¨ ur den darauf folgenden Zeitschritt wird j zu j+1 und j+1 zu j+2. Der L¨osungsvorgang bei Matrizengleichungen beginnt mit der Pr¨ ufung auf Nichtsingularit¨at der Matrix A, was mit det(A) = 0 = − k 5 gegeben ist. Im Fortgang erfolgt die Bildung der inversen Matrix A −1 A −1 = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ − 1 k − 2 k +1 k 2 − 3 k 2 +4 k +1 k 3 − 4 k 3 +10 k 2 +6 k +1 k 4 − 5 k 4 +20 k 3 +21 k 2 +8 k +1 k 5 0 − 1 k − 2 k +1 k 2 − 3 k 2 +4 k +1 k 3 − 4 k 3 +10 k 2 +6 k +1 k 4 0 0 − 1 k − 2 k +1 k 2 − 3 k 2 +4 k +1 k 3 0 0 0 − 1 k − 2 k +1 k 2 0 0 0 0 − 1 k ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ mit bekannten Methoden. Damit verbleibt die Multiplikation mit der Matrix b zur L¨osung des Flussdichte-Spaltenvektors x des j+1-ten Zeitschrittes gem¨aß A −1 b = A −1 A ︸ ︷︷ ︸ E x A −1 b = x fort. E ist die Einheitsmatrix. <?page no="108"?> 94 Einf¨ uhrung in die Finite-Differenzen-Methode 9.2.3 Anwendungsbeispiel Als Anwendungsbeispiel dient das Nachrechnen eines Magnetfeld-Diffusionsvorgangs nach Abb. 9.3. Hierzu wurden f¨ unf geometrische Schritte x 0 bis x 4 festgelegt. Die zu vergleichenden Schritte sind mit dem MATLAB Data Cursor in der Abbildung gekennzeichnet. Der 0’te Zeitschritt (t = 0 ms) ist dadurch gekennzeichnet, dass am linken Rand, bei x = 0 m die Randbedingung B(x 0 ) gesetzt wurde. Die Gegen¨ uberstellung erfolgt mit dem zweiten Zeitschritt t = 0,2 ms. Die zur Nachrechnung erforderlichen Daten sind den Tabellen 9.1 und 9.2 zu entnehmen. 0 0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01 (Linker Rand) Weg x [m] (rechter Rand) 0 0.2 0.4 0.6 0.8 1 1.2 B [T] 1,4 ms 1,2 ms 1,0 ms 0,8 ms 0,6 ms 0,4 ms 0,2 ms 0,0 ms X: 0.002449 Y: 0.3035 X: 0.005102 Y: 0.03224 X: 0.007551 Y: 0.001555 X: 0.01 Y: 5.795e-05 X: 0 Y: 1 X: 0.002449 Y: 0.4668 X: 0.005102 Y: 0.1296 Abbildung 9.3: MATLAB-Ergebnis eines eindimensionalen Felddiffusionsvorgangs mit der Zeit t als Scharparameter f¨ ur Faktor k = 0,2272 in Tab. 9.2 Unter Einbezug der Dirichlet-Randbedingungen wird die Matrixgleichung mit <?page no="109"?> 9.2 L¨osung mit impliziter Methode nach Crank-Nicolson 95 Tabelle 9.1: Werkstoffangaben Angaben f¨ ur Werkstoff Kupfer: Permeabilit¨at μ 0 , [V s/ (Am)] 4 π 10 −7 = 1,2 10 −6 Spez. elektr. Leitf¨ahigkeit κ, [A/ (V m)] 56,2 10 6 1/ (2 μ κ), [m 2 / s] 0,0071 ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0 0 − 0, 2272 1, 4544 − 0, 2272 0 0 0 − 0, 2272 1, 4544 − 0, 2272 0 0 0 − 0, 2272 1, 4544 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ B i −1 ,j +1 B i,j +1 B i +1 ,j +1 B i +2 ,j +1 B i +3 ,j +1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0 0 0 0, 001 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ berechnet. Die Orts- und Zeitangaben der Flussdichte B wurden dem MATLAB-Simulationsergebnis Abb. 9.3 entnommen und in Tab. 9.2 eingef¨ ugt. 0 0.002 0.004 0.006 0.008 0.01 x [mm] 0 0.2 0.4 0.6 0.8 1 B [T] Crank-Nicolson B(x 4 ) = 57,9 10 -6 T Crank-Nicolson B(x 4 ) = 0,001 T MATLAB PDE Abbildung 9.4: Vergleich zwischen MATLAB-PDE- und Crank-Nicolson-Ergebnis Die Ergebniswertepaare aus Tab. 9.2 wurden in Abb. 9.4 grafisch gegen¨ ubergestellt. Die Ann¨aherung des Crank-Nicolson-Ergebnisses an das MATLAB-Ergebnis wird durch <?page no="110"?> 96 Einf¨ uhrung in die Finite-Differenzen-Methode Tabelle 9.2: Gegen¨ uberstellung der Ergebnisse des Flussdichte- und Ortsverlaufs m. H. Abb. 9.3 Dirichlet-Randbedingungen: Linker Rand B(x 0 ), [T ] 1 Rechter Rand B(x 4 ), [T ] 0,001 und 57,9 10 −6 Zeitliche, r¨aumliche Diskretisierung, k-Faktor: Δt, [s] 0,0002 Δx, [m] 0,0025 Δt/ Δx 2 , [m/ s 2 ] 32 k, [1] 0,2272 x-Position B(t=0 s) B(t=0,2 ms) B(t=0,2 ms) B(t=0,2 ms) [m] [T ] MATLAB Crank Nicolson Crank Nicolson 0 1,0 1,0 1,0 1,0 0,0025 0 0,3035 0,25 0,0145 0,005 0 0,0322 0,04 0,0023 0,0075 0 0,0016 0,0064 0,37 10 −3 0,01 0 57,9 10 − 6 0,001 57,9 10 −6 Variation der Randbedingung B(x 4 ) ver¨andert. Wird in der Crank-Nicolson-Methode die Randbedingung B(x 4 ) = 57,9 10 −6 T (Wert aus Abb. 9.3) durch die Randbedingung B(x 4 ) = 1 mT ersetzt, so wird eine deutliche Ann¨aherung an das MATLAB-PDE- Ergebnis erzielt. 9.3 L¨osung mit expliziter Methode Im Fortgang wird die explizite Methode zur L¨osung der Felddiffusionsgleichung Gl. (9.1), vorgestellt. <?page no="111"?> 9.3 L¨osung mit expliziter Methode 97 9.3.1 ¨ Uberf¨ uhrung der Diffusionsgleichung in eine Matrizengleichung Hierzu wird die Gl. (9.4) B i +1 ,j − 2B i,j + B i −1 ,j (Δx) 2 = μκ B i,j +1 − B i,j Δt herangezogen. Mit der Substitution des Faktors k gem. Gl. (9.5) folgt k B i +1 ,j − 2k B i,j + k B i −1 ,j = B i,j +1 − B i,j . Mit dem Umsortieren der einzelnen Terme nach dem Zeitschritt j und j+1 wird k B i −1 ,j + (1 − 2k) B i,j + k B i +1 ,j = B i,j +1 erreicht. In Matrixschreibweise folgt Abbildung 9.5: Schritteinteilung der expliziten Methode [k (1 − 2k) k] ︸ ︷︷ ︸ A ⎛ ⎜ ⎝ B i −1 ,j B i,j B i +1 ,j ⎞ ⎟ ⎠ ︸ ︷︷ ︸ b = B i,j +1 ︸ ︷︷ ︸ x . (9.14) <?page no="112"?> 98 Einf¨ uhrung in die Finite-Differenzen-Methode Der linke Term von Gl. (9.14) ist bekannt, da dieser den gegenw¨artigen Zeitschritt und bekannte geometrische Schritte beinhaltet. Abb. 9.5 ist die Schritteinteilung zu entnehmen. 9.3.2 L¨osung der Matrizengleichung Nach erfolgter Berechnung der Funktionswerte des gegenw¨artigen Zeitschrittes j kann daraus der Funktionswert des zuk¨ unftigen Zeitschrittes j + 1 explizit berechnet werden. Die lineare Gl. (9.14) wird auf ein m/ n-System ⎛ ⎜ ⎝ k (1 − 2k) k 0 0 0 k (1 − 2k) k 0 0 0 k (1 − 2k) k ⎞ ⎟ ⎠ ︸ ︷︷ ︸ A · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ B i −1 ,j B i,j B i +1 ,j B i +2 ,j B i +3 ,j ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ︸ ︷︷ ︸ b = ⎛ ⎜ ⎝ B i,j +1 B i +1 ,j +1 B i +2 ,j +1 ⎞ ⎟ ⎠ ︸ ︷︷ ︸ x erweitert. Abb. 9.6 zeigt die Erweiterung der geometrischen Schritte. Die Anzahl der Knoten im j+1’ten Zeitschritt entspricht der Anzahl der Knoten im j-ten Zeitschritt vermindert um 2. Abbildung 9.6: Erweiterte Schritteinteilung der expliziten Methode <?page no="113"?> 9.3 L¨osung mit expliziter Methode 99 9.3.3 Anwendungsbeispiel In Abb. 9.7 ist das MATLAB-Ergebnis eines Felddiffusionsvorgangs ersichtlich. Nachgerechnet wird der Diffusionsvorgang von der Zeit 2 ms auf die Zeit 4 ms. Die zur Nachrechnung erforderlichen Daten sind in Tab. 9.3 zusammengefasst. Mit dem dort errechneten k-Faktor konnte die beste Ann¨aherung an das MATLAB-Ergebnis erzielt werden. Die zu l¨osende Matrixgleichung ist 0 0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01 (Linker Rand) Weg x [m] (rechter Rand) 0 0.2 0.4 0.6 0.8 1 1.2 B [T] 1,4 ms 1,2 ms 1,0 ms 0,8 ms 0,6 ms 0,4 ms 0,2 ms 0,0 ms X: 0.00102 Y: 0.7617 X: 0.00102 Y: 0.668 X: 0.002041 Y: 0.3911 X: 0.002041 Y: 0.5442 X: 0.003061 Y: 0.1984 X: 0.003061 Y: 0.3631 X: 0.004082 Y: 0.08652 X: 0 Y: 1 Abbildung 9.7: MATLAB-Ergebnis eines eindimensionalen Felddiffusionsvorgangs mit der Zeit t als Scharparameter f¨ ur Faktor k = 1,42 in Tab. 9.3 ⎛ ⎜ ⎝ 1, 42 − 1, 84 1, 42 0 0 0 1, 42 − 1, 84 1, 42 0 0 0 1, 42 − 1, 84 1, 42 ⎞ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0, 668 0, 391 0, 198 0, 087 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎝ 0, 74 0, 51 0, 31 ⎞ ⎟ ⎠ . <?page no="114"?> 100 Einf¨ uhrung in die Finite-Differenzen-Methode Tabelle 9.3: Beispiel 1: Gegen¨ uberstellung der Ergebnisse des Flussdichte- und Ortsverlaufs m. H. Abb. 9.7 Dirichlet-Randbedingungen: Linker Rand B(x 0 ), [T ] 1 Zeitliche, r¨aumliche Diskretisierung, k-Faktor: Δt, [s] 0,0002 Δx, [m] 0,001 Δt/ Δx 2 , [m/ s 2 ] 200 k, [1] 1,42 x-Position B(t=0,2 ms) [T ] B(t=0,4 ms) [T ] B(t=0,4 ms) [T ] [m] MATLAB MATLAB Explizit 0 1,0 - - 0,001 0,67 0,76 0,74 0,002 0,39 0,54 0,51 0,003 0,20 0,36 0,31 0,004 0,09 - - <?page no="115"?> 9.3 L¨osung mit expliziter Methode 101 Im Fortgang wird eine weitere Berechnung durchgef¨ uhrt und dabei der Faktor k so gew¨ahlt, dass (1-2k) = 0 wird. Die Matrixgleichung wird zu ⎛ ⎜ ⎝ 0, 5 0 0, 5 0 0 0 0, 5 0 0, 5 0 0 0 0, 5 0 0, 5 ⎞ ⎟ ⎠ · ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 0, 499 0, 157 0, 0347 0, 0069 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = ⎛ ⎜ ⎝ 0, 579 0, 267 0, 082 ⎞ ⎟ ⎠ . Die dazu erforderlichen Daten zur Nachrechnung sind der Abb. 9.8 zu entnehmen und in Tab. 9.4 zusammengefasst. Das dritte Rechenbeispiel ist in Tab. 9.5 ersichtlich. Diese abschließende Rechnung erlaubt eine Einsch¨atzung des k-Faktor-Einflusses (Wahl des k-Faktors) auf die Berechnungsergebnisse. 0 0.001 0.002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01 (Linker Rand) Weg x [m] (rechter Rand) 0 0.2 0.4 0.6 0.8 1 1.2 B [T] 1,4 ms 1,2 ms 1,0 ms 0,8 ms 0,6 ms 0,4 ms 0,2 ms 0,0 ms X: 0.001608 Y: 0.4992 X: 0.001608 Y: 0.6328 X: 0.005025 Y: 0.1354 X: 0.005025 Y: 0.03474 X: 0.003367 Y: 0.3171 X: 0.003367 Y: 0.1572 X: 0.006432 Y: 0.05606 X: 0.006432 Y: 0.006885 X: 0 Y: 1 Abbildung 9.8: MATLAB-Ergebnis eines eindimensionalen Felddiffusionsvorgangs mit der Zeit t als Scharparameter f¨ ur Faktor k = 0,5 in Tab. 9.4 <?page no="116"?> 102 Einf¨ uhrung in die Finite-Differenzen-Methode Tabelle 9.4: Beispiel 2: Gegen¨ uberstellung der Ergebnisse des Flussdichte- und Ortsverlaufs m. H. Abb. 9.8 Dirichlet-Randbedingungen: Linker Rand B(x 0 ), [T ] 1 Zeitliche, r¨aumliche Diskretisierung, k-Faktor: Δt, [s] 0,0002 Δx, [m] 0,001685 Δt/ Δx 2 , [m/ s 2 ] 70,44 k, [1] 0,5 x-Position B(t=0,2 ms) [T ] B(t=0,4 ms) [T ] B(t=0,4 ms) [T ] [m] MATLAB MATLAB Explizit 0 1,0 - - 0,00168 0,499 0,633 0,579 0,00337 0,157 0,317 0,267 0,005 0,035 0,135 0,082 0,0064 0,0069 - - <?page no="117"?> 9.3 L¨osung mit expliziter Methode 103 Tabelle 9.5: Beispiel 3: Gegen¨ uberstellung der Ergebnisse des Flussdichte- und Ortsverlaufs m. H. Abb. 9.3 Dirichlet-Randbedingungen: Linker Rand B(x 0 ), [T ] 1 Zeitliche, r¨aumliche Diskretisierung, k-Faktor: Δt, [s] 0,0002 Δx, [m] 0,0025 Δt/ Δx 2 , [m/ s 2 ] 32 k, [1] 0,2272 x-Position B(t=0,2 ms) [T ] B(t=0,4 ms) [T ] B(t=0,4 ms) [T ] [m] MATLAB MATLAB Explizit 0 1,0 - - 0,0025 0,3022 0,466 0,398 0,005 0,0032 0,1296 0,086 0,0075 0,0016 0,025 0,0079 0,01 57 10 −6 - - <?page no="118"?> Kapitel 10 Anwendungen der FEM zur Produktentwicklung Die Anwendung der FEM beschr¨ankt sich f¨ ur den Anwender i. d. R. auf die Bedienung geeigneter FEM-Programme. Allen gemeinsam ist die dreischrittige Vorgehensweise in der Programmanwendung, welche in die drei Phasen • Preprocessing • Processing • Postprocessing gegliedert werden kann. Als Beispiele zur Anwendung der FEM wird 1. das Nachrechnen (Analyse) eines bereits vorhandenen Proportionalmagnetaktors 2. das Vorausberechnen (Entwurf) eines planaren Asynchron-Scheibenl¨aufermotors vorgestellt. 10.1 Analyse eines Proportionalmagnetaktors Als Beispiel zur Nachrechnung eines vorhandenen Produkts dient der Proportionalmagnetaktor mit geometrischer Kennlinienbeeinflussung nach Abb. 10.1 von der Fa. Robert Bosch GmbH, welcher an der Hubschieber-Reiheneinspritzpumpe zur Regelung der Kraftstofff¨orderung eingesetzt wird. Die Benennung der Magnetelemente erfolgt in Tab. 10.1. 104 <?page no="119"?> 10.1 Analyse eines Proportionalmagnetaktors 105 Abbildung 10.1: Teilschnitt eines Proportionalmagnetaktors Tabelle 10.1: Bezeichnungen des Proportionalmagnetaktors Nr. Bezeichnung Nr. Bezeichnung 1 magnetischer R¨ uckschluss (Stator) 4 Anker mit Wirkelement 2 nichtmagnetische H¨ ulse 5 orthozyklische Wicklung 3 Gleitlager 6 Umspritzung 10.1.1 Preprocessing Die Phase des Preprocessing beinhaltet das Zeichnen der gew¨ unschten Kontur mit anschließender Vernetzung. Das Zeichnen erfolgt entweder mittels des in der FEM- Software integrierten Grafikeditors oder durch Zeichnen mittels eines externen Grafikeditors und anschließender Importierung der Grafik in die FEM-Software. Desgleichen erfolgt die Vernetzung entweder direkt mit der FEM-Software oder durch Vernetzung der extern erstellten Grafik mit anschließendem Import in die FEM-Software. In Abb. 10.2 a) ist ein 2D-L¨angsschnitt des Proportionalmagnetaktors nach Abb. 10.1 mit dem umgebenden Luftraum abgebildet. Die Vernetzung (Meshing) aller Bauelemente einschließlich des Luftraums ist in Abb. 10.2 b) erfolgt. Die Rotationsachse befindet sich jeweils im Bild links der Geometrie. Den einzelnen Bauelementen werden stoffliche Eigenschaften sowie Randbedingungen zugewiesen. Siehe hierzu auch Kap. 1.2.6. <?page no="120"?> 106 Anwendungen der FEM zur Produktentwicklung Abbildung 10.2: Preprocessing am Beispiel eines Proportionalmagnetaktors 10.1.2 Processing Die Phase des Processings umfasst die eigentliche L¨osung der bei der FEM verwendeten Gleichungen. Unterschieden werden zwei grundlegende L¨osungsmethoden: • Direkte Methode: Die L¨osung wird in einem (einzigen) riesigen Rechenschritt erarbeitet. Als Vertreter sei die Galerkin-Methode genannt. Die Solver arbeiten nach der vom Polnischen Mathematiker, Astronom und Geodesist Tadeusz Banachiewicz 1938 ver¨offentlichten ”LU decomposition-method“. Dabei steht LU f¨ ur ”lower“ und ”upper“ und ”decomposition“ f¨ ur ”Zerlegung“. Die LU-Zerlegung bezieht sich auf die quadratische Matrix A in eine Matrix L, in welcher die Elemente links unten bis zur Diagonalen besetzt sind und der Matrix U, in welcher die Elemente rechts oberhalb der Diagonalen besetzt sind. Damit sei A = L · U. Als n¨ utzliche Literaturstelle sei hierzu auf [24], S. 405 ff. verwiesen. • Iterative Methode: Die iterativen Methoden n¨ahern sich der L¨osung schrittweise an. Zu nennen sind hier ”conjugate gradient method“, ”generalized minimum residual method“ sowie ”biconjugate gradient stabilized method“. Diese Methoden erlauben w¨ahrend einer konvergierenden Berechnung die Beobachtung des kleiner werdenden Fehlers (Konvergenz) und der Anzahl der bereits berechneten Schritte. <?page no="121"?> 10.1 Analyse eines Proportionalmagnetaktors 107 Gut konditionierte Rechenprobleme weisen eine monotone Konvergenz auf. Zeigt sich jedoch nur ein verlangsamtes oder ein oszillierendes Konvergenzverhalten, deutet dies auf ein weniger gut konditioniertes Rechenproblem hin. Der Software-Anwender kann bei kommerzieller Software i. d. R. keinen Einfluss auf den L¨osungsalgorithmus nehmen. Selbst die erforderlichen Randbedingungen werden seitens der Software bereits vorgeschlagen oder angepasst. 10.1.3 Postprocessing Die Phase des Postprocessings beinhaltet die eigentliche Analysephase, in welcher das Ergebnis der gel¨osten Variabel mittels farblicher Codierung dargestellt wird und interpretiert werden muss. In Abb. 10.3 sind die magnetische Flussdichte B farblich und die Isolinien des magnetischen Flusses Φ grau hinterlegt. Abbildung 10.3: Postprocessing am Beispiel eines Proportionalmagnetaktors Die Farben rot und gelb verk¨orpern hohe, die Farben blau und gr¨ un niedrige Flussdichten. <?page no="122"?> 108 Anwendungen der FEM zur Produktentwicklung 10.2 Synthese eines planaren Asynchron-Scheibenl¨aufermotors Im Rahmen eines ¨offentlich gef¨orderten Projekts der AiF Projekt GmbH (ZIM-Kooperationsprojekt) wurde ein Dreiphasen-Asynchron-Scheibenl¨aufermotor mit Hilfe der FEM entwickelt und im Anschluss gem¨aß den errechneten Daten gebaut. Das magnetische Wanderfeld der Statorwicklung induziert im L¨aufer einen Strom, welcher seinerseits ein Magnetfeld ausbildet und mit dem umlaufenden Wanderfeld interagiert. Als Werkstofff¨ ur die magnetische Flussf¨ uhrung wird ein Pulververbundwerkstoffeingesetzt. Dieser zeichnet sich durch eine sehr niedrige spezifische elektrische Leitf¨ahigkeit aus. Das Anlaufmoment betr¨agt 0,3 Nm. Im Fokus der Entwicklung stand das zu reduzierende Motorgewicht. Der Motor besticht durch seine flache Bauweise. Als FEM-Tool wurde JAMG der Fa. Powersys erfolgreich angewendet. 10.2.1 Preprocessing Die Phase des Preprocessings ist in Abb. 10.4 dargestellt. In Abb. 10.4 a) werden die Motorelemente mittels JMAG-Editors gezeichnet. In Abb. 10.4 b) ist die automatische Vernetzung der Kontur erfolgt. Abbildung 10.4: Preprocessing am Beispiel des planaren Asynchronmotors 10.2.2 Processing Die Vorgehensweise des Processings erfolgt gem. Kap. 10.1.2. <?page no="123"?> 10.2 Synthese eines planaren Asynchron-Scheibenl¨aufermotors 109 10.2.3 Postprocessing In Abb. 10.5 ist die magnetische Flussdichte B in farblicher Codierung dargestellt. Hier bedeutet rot eine hohe Flussdichte und violett eine niedrige Flussdichte. Dargestellt ist der Stator und der Kurzschlussl¨aufer (Kurzschlussst¨abe mit Kurzschlussringen). Abbildung 10.5: Postprocessing am Beispiel des planaren Asynchronmotors 10.2.4 Musterbau des planaren Asynchronmotors In Abb. 10.6 ist der Prototyp des planaren Asynchron-Scheibenl¨aufermotors ersichtlich. Der Motor besteht aus zwei Statoren mit einem Doppelk¨afigrotor. In der Abbildung wurde ein Stator entfernt, um den Blick auf den Rotor zu erm¨oglichen. Die erforderlichen Bezeichnungen sind in Tab. 10.2 hinterlegt. Bez¨ uglich der Fertigung sei angemerkt, dass die beiden Statoren wie auch der Rotor aus Scheiben, bestehend aus Pulververbundwerkstoff, hergestellt wurden. Die erforderlichen Konturen f¨ ur die Leiterbahnen und Kurzschlussst¨abe wurden gefr¨ast. <?page no="124"?> 110 Anwendungen der FEM zur Produktentwicklung Tabelle 10.2: Bezeichnungen des Asynchron-Scheibenl¨aufermotors Nr. Bezeichnung Nr. Bezeichnung 1 Pulververbundwerkstoff (Rotor) 4 Kurzschlussring innen (Rotor) 2 Kurzschlussring außen (Rotor) 5 Distanzh¨ ulse 3 Kurzschlussstab (Rotor) 6 Rotor-Abtriebswelle Abbildung 10.6: Musteraufbau des planaren Asynchron-Scheibenl¨aufermotors <?page no="125"?> Kapitel 11 Anwendung der FEM zur Produktoptimierung In Abb. 11.1 ist die Einbettung einer FEM-Software in eine Optimierungssoftware ersichtlich. Dieses Vorgehen erlaubt die Optimierung mit vielen unabh¨angigen Variablen (multikriterielle Optimierung). Bei der Mehrzieloptimierung (multikriterielle Optimierung) werden an eine Probleml¨osung mehrere Anforderungen gestellt, die es bestm¨oglich zu erf¨ ullen gilt. Diese Anforderungen sind typischerweise gegens¨atzlich, sodass ein Optimum bez¨ uglich aller Funktionen nicht mit einer L¨osung zu erreichen ist. Ein Beispiel ist die Pareto-Optimierung, benannt nach dem italienischen ¨ Okonom Vilfredo Pareto (1848 - 1923). Ein Pareto-Optimum ist ein Zustand, in dem es nicht m¨oglich ist, eine L¨osung besser zu stellen, ohne zugleich eine andere L¨osung schlechter zu gestalten. Die Menge der Pareto-optimalen Punkte heißt Pareto-Front. Ein Berechnungsergebnis einer Mehrzieloptimierung ist der Abb. 11.2 zu entnehmen. Ersichtlich ist eine Motormoment-Optimierung in Abh¨angigkeit der geometrischen Parameter d 1 und d 2 . Alle Punkte und Kreise kennzeichnen dabei einer L¨osung. Die m¨oglichen L¨osungen verteilen sich entlang der Pareto-Front. Mehrzieloptimierverfahren suchen bei der Pareto-Optimierung nicht nur nach einer m¨oglichst guten L¨osung, sondern nach einer Menge von Kompromissl¨osungen, aus denen der Anwender eine zur Realisierung ausw¨ahlt. Weitere Beispiele zur Optimierung eines rotatorischen Antriebs (permanentmagnetisch erregten Synchronmaschine) sind der Dissertationsschrift [17] zu entnehmen. 111 <?page no="126"?> 112 Anwendung der FEM zur Produktoptimierung Abbildung 11.1: Vorgehen bei der Optimierung mit FEM <?page no="127"?> 113 Abbildung 11.2: Simulationsergebnis einer Mehrzieloptimierung <?page no="128"?> Anhang A A.1 MATLAB-Code - W¨armediffusionsskript % -------------------------------------------------------------- % Reinhold-Wuerth Hochschule Campus Kuenzelsau % Autor: Prof. Dr.-Ing. J. Ulm % MATLAB-Programm zur Loesung der thermischen Diffusionsgleichung % Datum: Sommerurlaub 2016 % -------------------------------------------------------------function pde_Diffusion close all; clear all; clc; t = [0: 0.5: 4]; % Linearvektor %t = linspace(0,1.0E-3,10); % Linearvektor x = linspace(0,0.05,20); % Linearvektor<<< %x = linspace(0,0.05,150); % Linearvektor m = 0; % Ordnungszahl Symmetrieproblem der Ebene % Loest Randwertprobleme elliptischer und parabolischer PDE’s u = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t); % Matrix u: Zeilen beinhalten Weginfo, beginnend bei null % Spalten beinhalten Zeitinfo, beginnend bei null 114 <?page no="129"?> A.1 MATLAB-Code - W¨armediffusionsskript 115 % Rufe Funktion auf um Werkstoffname zu erhalten Mat_Name = ’Kupfer’; % Figure-Titel % Position Figure links = 20; % Bezugskoordinate linker Bildrand (vertikal) unten = 300; % Bezugskoordinate unterer Bildrand (horizontal) breite = 600; % Bildbreite augehend von ’Bezugskoordinate linker Bildrand’ hoehe = 350; % Bildhoehe ausgehend von ’Bezugskoordinate unterer Bildrand’ figure(’position’,[links,unten,breite,hoehe]); surf(x,t,u,’Linewidth’,1.5); set(gca,’Fontsize’,12); view(8,12); %Titel = strcat(’Werkstoff: ’, Mat_Name, ’ \kappa / \kappa_0 = ’,num2str(Faktor)); Titel = strcat(’Werkstoff: ’, Mat_Name); title(Titel); xlabel(’(Linker Rand) Weg x [m] (rechter Rand)’); ylabel(’Zeit t [s]’); zlabel(’[ ◦ C]’); %print -depsc2 -tiff Diff1.eps links = 700; % Bezugskoordinate linker Bildrand (vertikal) unten = 300; % Bezugskoordinate unterer Bildrand (horizontal) breite = 600; % Bildbreite augehend von ’Bezugskoordinate linker Bildrand’ hoehe = 350; % Bildhoehe ausgehend von ’Bezugskoordinate unterer Bildrand’ figure(’position’,[links,unten,breite,hoehe]); %plot(rot90(u,3),’Linewidth’,1.2); plot(x,rot90(u,3),’Linewidth’,1.2); <?page no="130"?> 116 Titel = strcat(’Werkstoff: ’, Mat_Name); set(gca,’Fontsize’,12); xlabel(’(Linker Rand) Weg x [m] (rechter Rand)’); ylabel(’[ ◦ C]’); title(Titel); legend(’4,0 s’,’3,5 s’,’3,0 s’,’2,5 s’,’2,0 s’,’1,5 s’, ’1,0 s’,’0,5 s’,’0,0 s’) grid on; % Beispiele fuer Bildausgaben %print -depsc2 -tiff Diff2.eps %print -depsc2 -tiff Name.eps %print -dpng Name.png %print -dbmp Name.bmp % -------------------------------------------------------------- % Partielle Differenzialgleichung % Angaben aus Kuchling; Taschenbuch der Physik function [c,f,s] = pdex1pde(x,t,u,DuDx) Rho = 8933; % kg/ m^3; cth = 383; % J/ (kg K); Lamda = 384; % W/ (m K); c = Rho * cth / Lamda; f = DuDx; s = 0; % -------------------------------------------------------------- % Anfangsbedingung (initial conditions) function u0 = pdex1ic(x) u0 = 0; % -------------------------------------------------------------- % Randbedingungen <?page no="131"?> A.1 MATLAB-Code - W¨armediffusionsskript 117 function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t) % l = linker Rand, Anfang des Probekoerpers bei x = 0 % r = rechter Rand, Ende des Probekoerpers bei x = x_max % Randbedingung Beispiel 1- Kupfer - Referenz COMSOL % mit Randbedingung u0=0 % Temperatur = 100; % pl = ul-Temperatur; % ql = 0; % pr = ur-0; % qr = 0; %print -depsc2 -tiff Diffusion_Beispiel_1.eps % Randbedingung Beispiel 2 - Diff-Laenge gegnueber Bsp. 1 halbiert % mit Randbedingung u0=0 Temperatur = 100.0; pl = ul-Temperatur; ql = 0; pr = ur; qr = 20; %print -depsc2 -tiff Diffusion_Beispiel_2.eps % Randbedingung Beispiel 3 - Diff-Laenge gegnueber Bsp. 1 halbiert % mit Randbedingung u0=0 % Temperatur = 100; % pl = ul-Temperatur; % ql = 0; % pr = 0; % qr = ur-1 ; %print -depsc2 -tiff Diffusion_Beispiel_3.eps % -------------------------------------------------------------- <?page no="132"?> 118 A.2 MATLAB-Code - Magnetfelddiffusionsskript % -------------------------------------------------------------- % Reinhold-Wuerth Hochschule Campus Kuenzelsau % Autor: Prof. Dr.-Ing. J. Ulm % MATLAB-Programm zur Loesung der Magnetfeld-Diffusionsgleichung % Datum: Sommerurlaub 2016 % -------------------------------------------------------------function pde_Diffusion close all; clear all; clc; t = [0E-4: 2E-4: 1.4E-3]; % Linearvektor %t = linspace(0,1.0E-3,10); % Linearvektor x = linspace(0,0.010,50); % Linearvektor m = 0; % Ordnungszahl Symmetrieproblem der Ebene % Loest Randwertprobleme elliptischer und parabolischer PDE’s u = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t); % Matrix u: Zeilen beinhalten Weginfo, beginnend bei null % Spalten beinhalten Zeitinfo, beginnend bei null % rufe Funktion auf um Werkstoffname zu erhalten [my_kappa,Mat_Name] = Permeabilitaet_kappa(0.1); % Position Figure links = 20; % Bezugskoordinate linker Bildrand (vertikal) unten = 300; % Bezugskoordinate unterer Bildrand (horizontal) breite = 600; % Bildbreite augehend von ’Bezugskoordinate linker Bildrand’ hoehe = 350; % Bildhoehe ausgehend von ’Bezugskoordinate unterer Bildrand’ <?page no="133"?> A.2 MATLAB-Code - Magnetfelddiffusionsskript 119 figure(’position’,[links,unten,breite,hoehe]); surf(x,t,u,’Linewidth’,1.5); set(gca,’Fontsize’,12); view(8,12); %Titel = strcat(’Werkstoff: ’, Mat_Name, ’ \kappa / \kappa_0 = ’,num2str(Faktor)); Titel = strcat(’Werkstoff: ’, Mat_Name); %title(Titel); xlabel(’(Linker Rand) Weg x [m] (rechter Rand)’); ylabel(’Zeit t [s]’); zlabel(’B [T]’); %print -depsc2 -tiff Diff_B_Galerkin_01.eps links = 700; % Bezugskoordinate linker Bildrand (vertikal) unten = 300; % Bezugskoordinate unterer Bildrand (horizontal) breite = 600; % Bildbreite augehend von ’Bezugskoordinate linker Bildrand’ hoehe = 350; % Bildhoehe ausgehend von ’Bezugskoordinate unterer Bildrand’ figure(’position’,[links,unten,breite,hoehe]); %plot(rot90(u,3),’Linewidth’,1.2); plot(x,rot90(u,3),’Linewidth’,1.2); Titel = strcat(’Werkstoff: ’, Mat_Name); set(gca,’Fontsize’,12); xlabel(’(Linker Rand) Weg x [m] (rechter Rand)’); ylabel(’B [T]’); legend(’1,4 ms’, ’1,2 ms’, ’1,0 ms’, ’0,8 ms’, ’0,6 ms’, ’0,4 ms’, ’0,2 ms’, ’0,0 ms’); %title(Titel); grid on; % Beispiele fuer Bildausgaben %print -depsc2 -tiff Diff_B_Galerkin_02.eps %print -depsc2 -tiff Diff_B_Galerkin_02.eps <?page no="134"?> 120 %print -depsc2 -tiff Name.eps %print -dpng Name.png %print -dbmp Name.bmp % -------------------------------------------------------------- % Partielle Differenzialgleichung function [c,f,s] = pdex1pde(x,t,u,DuDx) [my_kappa,Mat_Name] = Permeabilitaet_kappa(u); c = my_kappa; f = DuDx; s = 0; % -------------------------------------------------------------- % Anfangsbedingung (initial conditions) function u0 = pdex1ic(x) u0 = 0; % -------------------------------------------------------------- % Randbedingungen function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t) % l = linker Rand, Anfang des Probekoerpers bei x = 0 % r = rechter Rand, Ende des Probekoerpers bei x = x_max % Randbedingung Beispiel 1- Kupfer - Referenz COMSOL % mit Randbedingung u0=0 % Tesla = 1.0; % pl = ul-Tesla; % ql = 0; % pr = ur-0.2; % qr = 0; %print -depsc2 -tiff Diffusion_Beispiel_1.eps % Randbedingung Beispiel 2 - Diff-Laenge gegnueber Bsp. 1 halbiert <?page no="135"?> A.2 MATLAB-Code - Magnetfelddiffusionsskript 121 % mit Randbedingung u0=0 Tesla = 1.0; pl = ul-Tesla; ql = 0; pr = ur; qr = 1; %print -depsc2 -tiff Diffusion_Beispiel_2.eps % Randbedingung Beispiel 3 - Diff-Laenge gegnueber Bsp. 1 halbiert % mit Randbedingung u0=0 % Tesla = 1.0; % pl = ul-Tesla; % ql = 0; % pr = 0; % qr = ur-1 ; %print -depsc2 -tiff Diffusion_Beispiel_3.eps % -------------------------------------------------------------- % -------------------------------------------------------------function [my_kappa,Mat_Name] = Permeabilitaet_kappa(B) % -------------------------------------------------------------- % zur Verfuegung stehende Werkstoffe Werkstoff_M = 2; % waehle Werkstoff aus (1 ... 5) perm = 3; % 1 = differenzeille Permeab.; % 2 = Permeabilitaet % 3 = Permeabilitaet des Vakuums if Werkstoff_M == 1 % Vacoflux 50 Mat_Name= [’Vacoflux50’]; H_mat = [0 13 20 22 26 30 34 39 45 52 60 69 79 91 104 138 242 447 782 2393 5514 15000 50000 ]; B_mat = [0 0.04 0.1 0.16 0.26 0.43 0.65 0.86 1.09 1.28 1.42 1.51 1.59 1.65 1.71 1.81 1.95 2.07 2.15 2.25 2.28 2.30 2.40]; <?page no="136"?> 122 kappa = 2.5.*10.^6; elseif Werkstoff_M == 2 % Kupfer Mat_Name= [’Kupfer’]; H_mat = [0 10 25000 50000]; B_mat = [0 1.2566E-5 0.0314159 0.062831853]; kappa = 56.2.*10.^6; elseif Werkstoff_M == 3 Mat_Name= [’9SMn28K’]; H_mat = [ 0 396 793 1031 1190 1587 3174 5952 9920 19841 30158 50000]; B_mat = [ 0 0.6 0.97 1.15 1.24 1.39 1.53 1.63 1.7 1.75 1.77 1.8 ]; kappa = 4.5.*10.^6; % Leitfaehigkeit geschaetzt elseif Werkstoff_M == 4 Mat_Name= [’11SMn 30’]; H_mat = [0 1250 2500 5000 10000 15000 20000 25000 30000 35000 40000 45000 50000 ]; B_mat = [0 0.95 1.29 1.65 1.81 1.89 1.95 2 2.02 2.04 2.05 2.07 2.09]; kappa = 4.55.*10.^6; elseif Werkstoff_M == 5 Mat_Name= [’100Cr6’]; H_mat = [0 100 300 500 1250 2500 5000 10000 15000 20000 25000 30000 35000 40000 45000 50000]; B_mat = [0 0.01 0.10 0.21 0.63 1.17 1.37 1.47 1.54 1.58 1.60 1.63 1.65 1.66 1.67 1.69]; kappa = 4.65.*10.^6; else error(’Waehle Werkstoff aus! ! ! ’) end my_0 = 4*pi*10^-7; if perm == 1 % Differenzielle Permeabilitaet dB = diff(B_mat); <?page no="137"?> A.2 MATLAB-Code - Magnetfelddiffusionsskript 123 dH = diff(H_mat); my_mat = zeros(size(H_mat)); my_mat(1: end-1) = dB./ dH; my_kappa = interp1(B_mat,my_mat,B,’pchip’) .* kappa; elseif perm == 2 % Permeabilitaet my_mat = B_mat./ (H_mat+eps); my_mat(1) = my_mat(2); % unterdruecke NAN-Ausgabe my_kappa = interp1(B_mat,my_mat,B,’pchip’) .* kappa; elseif perm == 3 % Permeabilitaet des Vakuums my_kappa = kappa.*my_0; else error(’Waehle Permeabilitaet aus! ! ! ’) end % -------------------------------------------------------------- <?page no="138"?> 124 A.3 MATLAB-Ergebnisse vs. COMSOL Multiphysics-Ergebnisse Gegen¨ ubergestellt werden die beiden Simulationsergebnisse eines eindimensionalen Felddiffusionsvorgangs in einem unendlich langen Kupferk¨orper mit einer spezifischen elektrischen Leitf¨ahigkeit κ = 56 10 6 1/ (Ωm) und konstanter Permeabilit¨at μ 0 . Als Randbedingung wurde einseitig die Flussdichte von 1 T vorgegeben. In der Abb. A.1 ist das mit der MATLAB PDE-Toolbox und in Abb. A.2 das mit COMSOL Multiphysics erzielte Ergebnis ersichtlich. In beiden Ergebnisdarstellungen wurde die Zeit t als Scharparameter gew¨ahlt. Die Ergebnisse beider Tools stimmen sehr gut ¨ uberein. Der MATLAB-Code kann dem Anhang A.2 entnommen werden. 0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 x [m] 0 0.2 0.4 0.6 0.8 1.0 Flussdichte B [T] 1,0 ms 0,9 ms 0,8 ms 0,7 ms 0,6 ms 0,5 ms 0,4 ms 0,3 ms 0,2 ms Abbildung A.1: MATLAB-PDE-Toolbox-Simulationsergebnis <?page no="139"?> A.3 MATLAB-Ergebnisse vs. COMSOL Multiphysics-Ergebnisse 125 0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 x [m] 0 0.2 0.4 0.6 0.8 1.0 Flussdichte B [T] 1,0 ms 0,9 ms 0,8 ms 0,7 ms 0,6 ms 0,5 ms 0,4 ms 0,3 ms 0,2 ms Abbildung A.2: COMSOL Multiphysics-Simulationsergebnis <?page no="140"?> Literaturverzeichnis [1] Bastos, J. ; Ida, N. : Electromagnetics and Calculation of Fields; 2 Auflage. Springer, 1997 [2] Bastos, J. ; Sadowski, N. : Electromagnetic Modeling by Finite Element Methods. Marcel Dekker, Inc., 2003 [3] Bronstein, I. N. ; Semendjajew, K. A. ; Musilo, G. ; M¨ uhlig, H. : Taschenbuch der Mathematik. 5. Auflage. Verlag Harri Deutsch, 2000 [4] Burke-Hubbard, B. : Wavelets. Die Mathematik der kleinen Wellen. Birkh¨auser Verlag, 1997 [5] DIN-1301 : Einheiten - Einheiten¨ahnliche Namen und Zeichen; Beiblatt zu Teil 1. Beuth Verlag GmbH, 1982 [6] DIN-1301 : Einheiten - Teil 1: Einheitennamen, Einheitenzeichen. Beuth Verlag GmbH, 2002 [7] DIN-1302 : Allgemeine mathematische Zeichen und Begriffe. Beuth Verlag GmbH, 1999 [8] DIN-1303 : Vektoren, Matrizen, Tensoren - Zeichen und Begriffe. Beuth Verlag GmbH, 1987 [9] DIN-1304 : Formelzeichen - Teil 1: Allgemeine Formelzeichen. Beuth Verlag GmbH, 1994 [10] DIN-1338 : Formelschreibweise und Formelsatz. Beuth Verlag GmbH, 2011 [11] DIN-4895 : Orthogonale Koordinatensysteme - Teil 2: Differentialoperatoren der Vektoranalysis. Beuth Verlag GmbH, 1977 126 <?page no="141"?> LITERATURVERZEICHNIS 127 [12] DIN-4895 : Orthogonale Koordinatensysteme - Teil 1: Allgemeine Begriffe. Beuth Verlag GmbH, 1997 [13] Fetzer, J. ; Haas, M. ; Kurz, S : Numerische Berechnung elektromagnetischer Felder - Band 627. expert verlag GmbH, 2002 [14] Fletcher, C.A.J. : Computational Galerkin Methods. Springer, 1984 [15] Harrington, R. F. : Field Computation by Moment Methods. New York; Macmillan, 1968 [16] Jung, M. ; Langer, U. : Methode der finiten Elemente f¨ur Ingenieure; 2. Auflage. Springer/ Vieweg, 2013 [17] Krotsch, J. : Mehrkriterielle Optimierung permanentmagneterregter Synchronmotoren in Außenl¨auferbauweise unter besonderer Ber¨ucksichtigung der Radialkr¨afte. Dissertation Friedrich-Alexander-Universit¨at, Erlangen-N¨ urnberg, 2016 [18] Kuchling, H. : Taschenbuch der Physik. Verlag Harri Deutsch, 1991 [19] Liedl, R. ; Kuhnert, K. : Analysis in einer Variablen. B. I. Wissenschaftsverlag M¨ unchen, 1992 [20] Marsden, J. E. ; Tromba, A. J. : Vector Calculus. W. H. Freeman and Company, 2000 [21] Munz, C.D. ; Westermann, T. : Numerische Behandlung gew¨ohnlicher und partieller Differenzialgleichungen; 2. Auflage. Springer, 2009 [22] Sadiku, M. N. O. : Numerical Techniques in Elektromagnetics, 2nd ed. CRC Press LLC, 2001 [23] Sawitzki, A. : Zur Berechnung der elektromagnetischen Felder von ausgedehnten komplexen Systemen durch Erweiterung der Momentenmethode um eine effiziente Rasterung. Dissertation Technische Universit¨at Ilmenau, 2017 [24] Schwarzenberg-Czerny, A. : On Matrix Factorization and Efficient Least Squares Solution. Astronomy and Astrophysics Supplement, v.110, p.405, 1995 [25] Weller, F. : Numerische Mathematik f¨ur Ingenieure und Naturwissenschaftler. Verlag Vieweg, 1996 <?page no="142"?> Dipl.-In Sim Grundl 2., erw. A (Kontakt ISBN 97 Zum Buc Das Fachb wie man e das System Einsatz g Simulink e Inhalt: Der Autor aus Wisse verschiede der Autor GPSS für behandelt der Leser behandelte Die Intere Das Buch Produktion alle, die In die Simula möchten. Rezensio »Die ausf Systemsim Der Auto Bernd Ack Bremsrege Federungs seit vielen ng. Bernd mulat agen un Aufl. 2011, & Studium 8-3-8169-2 ch: buch führt in ein System b m simuliert w geeigneter P ein Beispiel k stellt einleite enschaft un enen Anwen systemorien die zeitdisk und praktisc ein vertieftes en und weite essenten: richtet sich n, an Techni nteresse an ation als Too onen: führliche Int mulation« or: ker simuliert elsysteme u ssysteme, vo Jahren als L d Acker tions d praktis , 160 S., 38 m, 690) 2999-4 n die Welt d eschreibt un werden kann Programmier kennen. end Simulat d Technik v dungsgebiet ntierte Sprac kreten Syste ch erprobt. D s Verständni ere einfache an Mitarbeit iker und Kau der Simulat ol bei ihrer tä erpretation te zunächst nd entwicke or allem der Lehrbeauftra Be Tel: 0715 E-Mail: ex stech sche Anw 8,80 €, 64, der Simulatio nd daraus ein n. Er erhält G rsprachen u ion im Rahm veranschaul ten. Zur Sim chen ein: AC eme. Die Sp Durch eine au s der System Matlab Simu ter aus Entw ufleute, insbe tionstechnik äglichen Arbe der Simulat Teilkompon elte geeignet aktiven Fahr gter weiter. estellhotl 59 / 92 65-0 xpert@expe hnik wendung ,50 CHF on ein. Der L n Modell able Grundkenntni und lernt m men von Sys icht er die mulation zeitk CSL und Ma prachen we usführliche In msimulation. ulink-Beispie wicklung und esondere an haben bzw eit einsetzen tionsergebnis nenten von F te Regelalgo rzeugfederun line: 0 • Fax: -20 ertverlag.de gen Leser erfähr eitet, mit dem isse über de mit MATLAB stemstudien große Bede kontinuierlich atlab Simulin rden jeweils nterpretation Zur eigenen le als Simulin d n . n sse vermitte Fahrzeugen orithmen. Sp ng. Seine Er e t, m n B dar. Anhand eutung der er und zeitd nk für die ze s anhand an der Simulat n Anwendung nk Programm elt ein vertie wie z.B. Au päter wurde rfahrung auf d zahlreicher Simulations diskreter Sys eitkontinuierl nschaulicher tionsergebnis g erhält er d mfile. eftes Verstä Instan utomatik-Get er Experte diesen Feld r Beispiele stechnik in steme führt lichen und r Beispiele sse erfährt ie im Buch ändnis der ndhaltung triebe und für PKWern gibt er <?page no="143"?> Dipl.-Ing. Dr. sc. techn. ETH/ SIA Yasar Deger pg Die Methode der Finiten Elemente Grundlagen und Einsatz in der Praxis 8. Aufl. 2017, 160 S., 131 Abb., 29,80 €, 38,80 CHF (Kontakt & Studium, 551) ISBN 978-3-8169-3399-1 Zum Buch: Das Buch vermittelt Neueinsteigern, Anwendern und Entscheidungsträgern einen Überblick über Grundlagen, Möglichkeiten und Grenzen der FE-(Finite-Elemente-)Methode. Es erklärt die Arbeitsweise der zugehörigen Programme auf leicht verständliche Art und beschreibt die Voraussetzungen und Vorgehensschritte für den erfolgreichen und effizienten Einsatz in der Ingenieurpraxis. Dabei wird besonders darauf Wert gelegt, das physikalische Problem als ein mechanisch/ mathematisches Modell so einfach wie möglich, aber so genau wie nötig zu simulieren, die aussagekräftigen Ergebnisse aus der Berechnung zu selektieren und diese kritisch auszuwerten bzw. kompetent zu interpretieren. Die zahlreichen, einfach gehaltenen, anwendungsspezifischen Beispiele aus einer breiten Palette von Problemen mit Praxisbezug regen den Leser zum selbständigen Üben an. Inhalt: Theoretische Grundlagen - Finite Elemente als Ersatz für elastische Körper - Rotationssymmetrie - Modellierung des Materialverhaltens - Stabilitätsuntersuchungen / Nichtlinearitäten - Dynamische FE-Berechnungen - Thermische FE-Berechnungen - Regeln für den Umgang mit der FE-Methode - Übungen Die Interessenten: Entwicklungsingenieure, Konstrukteure, Berechnungsingenieure aus Maschinenbau, Fahrzeugbau, Luft- und Raumfahrt, Anlagenbau, Apparatebau, Bauwesen, Forschung und Entwicklung - Entscheidungsträger mit entsprechendem Zuständigkeitsbereich - Lehrer und Studenten an Fachhochschulen Rezensionen: »Im Buch werden die Grundlagen und der Einsatz der FEM in der Praxis detailliert vermittelt und anhand von 131 Bildern verdeutlicht.« NAFEM, Online-Magazin »Endlich mal ein praxisnahes Lehrbuch zur Thematik der Finiten Elemente! Zwar kommt die Theorie nicht zu kurz, doch unterfüttert der Autor das Ganze mit konkreten Anwendungen aus der technischen Mechanik. Dies macht die schwierige Materie extrem les- und lernbar. Eines der Bücher, nach deren Lektüre dauerhaftes Wissen zurückbleibt.« Design&Elektronik Der Autor: Dr. Y. Deger war über 30 Jahre als Dozent für Technische Mechanik und Finite Elemente Methode an der HSR, Hochschule für Technik, Rapperswil in der Schweiz tätig. Während dieser Zeit war er zugleich als Berechnungsingenieur und Experte in der Industrie engagiert und blickt auf breitgefächerte Erfahrung in der Anwendung der FEM in mehreren Bereichen zurück. Seit 1992 erteilt er zudem Weiterbildungskurse für die Ingenieure und Konstrukteure in der Praxis. Nach den ersten Seminarveranstaltungen an der Technischen Akademie Esslingen (TAE) folgten jene für NAFEMS sowie letztlich für das VDI- Wissensforum, welche sich großen positiven Echos seitens der Teilnehmer erfreuen. Blätterbare Leseprobe und einfache Bestellung unter: www.expertverlag.de/ 3399 Bestellhotline: Tel: 07159 / 92 65-0 • Fax: -20 E-Mail: expert@expertverlag.de <?page no="144"?> Dipl.-Ing. Nick Eckert P: \AK jpg Innovationskraft steigern mit LOBIM Eine praxisnahe Methodenkopplung von TRIZ und Bionik - Entwicklung und Konstruktion; Maschinen und Maschinenelemente 2017, 234 S., 143 farb. Abb., 39,00 €, 51,00 CHF (Reihe Technik) ISBN 978-3-8169-3325-0 Zum Buch: Das vorliegende Buch bietet dem Leser eine Vielzahl von Lösungsangeboten anhand innovativer Grundprinzipien aus Natur und Technik. Die klassischen 40 TRIZ Grundprinzipien werden durch bionische Prinzipien zu einer sehr einfachen, neuen Methode erweitert und anhand von Beispielen erklärt. Das Buch soll außerdem zu Innovation motivieren und Mut machen. Dazu wird erklärt, dass Innovation ein immanenter Bestandteil der natürlichen und menschlichen Evolution ist. Durch neueste Erkenntnisse der Kreativitätsforschung und anhand von Erfinderbiografien wird versucht, die innovative Persönlichkeit und den kreativen Prozess zur erklären. Es wird auf die Rahmenbedingungen für ein innovatives Klima ebenso eingegangen, wie auf die Notwendigkeit neuer Innovationsmethoden. Inhalt: Chemische, biologische und physikalische Evolution und Innovation - Kurze Geschichte der menschlichen Innovation - Merkmale innovativer Persönlichkeiten mit beispielhaften Biografien - Erkenntnisse aus der Kreativitätsforschung zum kreativen Prozess - Rahmenbedingungen für Innovation in Unternehmen - Kurzer Abriss bekannter Innovationsmethoden - LOBIM als neuartiger, praxisnaher Methodenverbund aus TRIZ und Bionik - Beispiele zur schnellen Inspiration bei der Lösungssuche - Anwendungsbeispiele von LOBIM aus der Automobilindustrie Die Interessenten: Das Buch ist ein Gewinn für Entwicklungsingenieure, private und professionelle Erfinder, Manager und Mitarbeiter von F&E, ebenso für technisch, naturwissenschaftlich und auch wirtschaftlich interessierte Studenten und Lehrer Rezensionen: »Das vorgelegte Buch ist hervorragend! Sorgfältig strukturiert, mit vielen konkreten Beispielen und Bildern visualisiert, regt es die Phantasie zum Erfinden an. Es sollte sowohl in der Lehre als auch in Entwicklungs- und Konstruktionsabteilungen gute Impulse generieren« GFPMagazin - Gesellschaft für Produktionsmanagement e.V. Der Autor: Nick Eckert ist Maschinenbauingenieur und seit 2002 Entwicklungsingenieur bei einem globalen Zulieferer für automobile Sicherheit. Er war als Vertriebsingenieur und Produktingenieur bei verschiedenen Industrieunternehmen sowie als Konstruktionsingenieur in der Automobilindustrie tätig. Zahlreiche Patente auf dem Gebiet der Fahrzeugsicherheit. Seit 2004 Tätigkeit als TRIZ-Moderator. Er beschäftigt sich nebenberuflich mit Bionik und ist Mitglied des Berliner Arbeitskreises TRIZ Blätterbare Leseprobe und einfache Bestellung unter: www.expertverlag.de/ 3325 Be Tel: 071 E-Mail: ex estellhot 159 / 92 65xpert@exp tline: 0 • Fax: -20 pertverlag.d 0 de