PDGL mittels FEM < partielle < Differentialgl. < Analysis < Hochschule < Mathe < Vorhilfe
|
Status: |
(Frage) überfällig | Datum: | 17:50 Sa 28.03.2009 | Autor: | sabadei |
Aufgabe | Gegeben ist folgende PDGL auf [mm] \Omega:
[/mm]
[mm] \bruch{\partial^2\,p}{\partial\,x^2}+\left(1+3\epsilon\,\cos\left(x\right)\right)\bruch{\partial^2\,p}{\partial\,y^2}=1, \qquad 0\le\epsilon\le\bruch{1}{3}, [/mm]
p(x)=0 [mm] \forall x\in\bruch{\partial\,\Omega}{\partial\,x}
[/mm]
Es ist eine Lösung mittels FEM gesucht. [mm] \Omega [/mm] wird mittels rechteckige 4-Knoten Elemente mit den Formfunktionen
[mm] N=[(1+\xi)(1+\eta) \quad (1-\xi)(1+\eta) \quad (1-\xi)(1-\eta) \quad (1+\xi)(1-\eta)] [/mm]
unterteilt. |
Hallo liebe Forummitglieder,
ich versuche gerade diese Aufgabe zu lösen, komme aber nicht weiter. Ich gehe dabei wie folgt vor:
Übergang zu der Variationsformulierung liefert
[mm] \int_{\Omega}\left(\bruch{\partial\,p}{\partial\,x}\bruch{\partial\,v}{\partial\,x}+\left(1+3\epsilon\,\cos\left(x\right)\right)\bruch{\partial\,p}{\partial\,y}\bruch{\partial\,v}{\partial\,x}\right)\,d\Omega=\int_{\Omega}1\,v\,d\Omega
[/mm]
und mit den Formfunktionen
[mm] \bigcup^{Elemente}\int_{\Omega}\left(N_{,x}^T N_{,x}+\left(1+3\epsilon\,\cos\left(x\right)\right)\N_{,y}^T N_{,y}\right)\,d\Omega P=\bigcup^{Elemente}\int_{\Omega}N^T\,d\Omega.
[/mm]
Wegen den Koordinatentransformationen
[mm] x(\xi)=\bruch{(x_a+x_e)}{2}+\bruch{(x_a-x_e)}{2}\xi
[/mm]
[mm] y(\eta)=\bruch{(y_a+y_e)}{2}+\bruch{(y_a-y_e)}{2}\eta
[/mm]
[mm] (x_a,x_e,y_a,y_e [/mm] - Elementegrenzen)
setze ich dabei
[mm] (a=\bruch{(x_a+x_e)}{2}, b=\bruch{(x_a-x_e)}{2}, c=\bruch{(y_a-y_e)}{2})
[/mm]
[mm] \int_{\Omega}N_{,x}^T N_{,x}\,d\Omega [/mm] = [mm] \int_{-1}^1 \int_{-1}^1 N_{,\xi}^T N_{,\xi}\bruch{bc}{b^2}\,d\xi\,d\eta
[/mm]
[mm] \int_{\Omega}\left(1+3\epsilon\cos\left(x\right)\right)N_{,y}^T N_{,y}\,d\Omega [/mm] = [mm] \int_{-1}^1 \int_{-1}^1\left(1+3\epsilon \cos\left(a+b\xi\right)\right)N_{,\eta}^T N_{,\eta}\bruch{bc}{c^2}\,d\xi\,d\eta
[/mm]
[mm] \int_{\Omega}N^T\,d\Omega=\int_{-1}^1\int_{-1}^1N^T bc\,d\xi\,d\eta.
[/mm]
Dann werden die Integrale ausgewertet. Für den 2. bekomme ich z.B.
[mm] I_1=\int_{-1}^1 \int_{-1}^1\left(1+3\epsilon \cos\left(a+b\xi\right)\right)(1+\xi)(1+\xi)\bruch{bc}{c^2}\,d\xi\,d\eta=\bruch{16}{3}+ \epsilon\bruch{(12 (-2 \cos(a) \sin(b) + 2 b (\cos(a+b) + b \sin(a+b))))}{b^3}
[/mm]
[mm] I_2=\int_{-1}^1 \int_{-1}^1\left(1+3\epsilon \cos\left(a+b\xi\right)\right)(1+\xi)(1-\xi)\bruch{bc}{c^2}\,d\xi\,d\eta=\bruch{8}{3}+\epsilon \bruch{(24 \cos(a) (-b \cos(b) + \sin(b)))}{b^3}
[/mm]
[mm] I_3=\int_{-1}^1 \int_{-1}^1\left(1+3\epsilon \cos\left(a+b\xi\right)\right)(1-\xi)(1-\xi)\bruch{bc}{c^2}\,d\xi\,d\eta=\bruch{16}{3} [/mm] + [mm] \epsilon \bruch{(8 (9 b \cos(a-b) - 9 (b^2 \sin(a-b) + \cos(a)\sin(b))))}{b^3}
[/mm]
und forme sie zu der lokalen Steifigkeitsmatrix zusammen:
[mm] \pmat{ I_1 & I_2& -I_2& -I_1 \\ I_2 & I_3& -I_3& -I_2\\ -I_2 & -I_3& I_3& I_2\\ -I_1 & -I_2& I_2& I_1}
[/mm]
Anschließend wird die globale Steifigkeitsmatrix zusammengesetzt. Die Lösung ist aber falsch. Der Fehler steigt mit dem zunehmenden [mm] \epsilon. [/mm] Für [mm] \epsilon=.6 [/mm] bekommt man etwa 10%, für [mm] \epsilon=.99 [/mm] etwa 40%. Da man für [mm] \epsilon [/mm] gegen 0 exakte Werte bekommt, vermute ich, dass etwas bei der Auswertung des Integrals mit [mm] \cos [/mm] falsch ist, oder irgendwo bei der Herleitung (Übergang zu der Variationsformulierung?) ein Fehler passiert.
(die Lösung habe ich mit der von Matlab (PDE-Toolbox) und Comsol verglichen. Beide Progamme bekommen etwa dasselbe raus. Bei der Erhöhung der Elementenzahl gibt es keine Verbesserungen, meine Lösung konvergiert also gegen einen anderen Wert. Integrale habe ich mit Mathematica und TI ausgewertet (in beiden Fällen gleich) und bereits sehr oft überprüft.)
Hoffe, dass jemand eine Idee hat.
Ich habe diese Frage in keinem Forum auf anderen Internetseiten gestellt.
|
|
|
|
Status: |
(Mitteilung) Reaktion unnötig | Datum: | 19:20 So 05.04.2009 | Autor: | matux |
$MATUXTEXT(ueberfaellige_frage)
|
|
|
|