Blog

U 3 – Die stationäre Wärmeleitungsgleichung

Die stationäre Wärmeleitungsgleichung lautet

\nabla \cdot \vec Q+r = 0.

Hierin beschreibt \vec Q den Wärmefluss über den Körperrand und r eine innere Wärmequelle. Die Primärvariable ist die Temperatur \theta, die mit dem Wärmeflussvektor über das isotrope Fourier’sche Wärmeleitgesetz

\vec Q = -\lambda \nabla \theta

in Verbindung steht. Der Materialparameter \lambda ist die Wärmeleitfähigkeit.

  1. Formulieren Sie die Differentialgleichung und die Randbedingungen für den eindimensionalen Fall des Wärmeleitstabes. Welche Randbedingungen werden im Kontext der FEM unterschieden, welche Restriktionen in Bezug auf die Randbedingungen müssen erfüllt sein? Warum sind keine Anfangsbedingungen zu formulieren?
  2. Leiten Sie die schwache Form der Differentialgleichung her. Welchen Bedingungen muss dabei die Testfunktion genügen? Bemerkung: Für die kinematische Beziehung kann {\gamma _{11}} = {\theta _{1,1}} benutzt werden.
  3. Approximieren Sie die unbekannten Größen in der schwachen Form für ein repräsentatives Element e mit linearen Ansätzen. Welche Bedingungen müssen Ansatzfunktionen allgemein erfüllen? Was verstehen Sie unter dem isoparametrischen Konzept?
  4. Schreiben Sie die approximierte, schwache Form nieder. Berechnen Sie analog zur Impulsbilanz die in der approximierten, schwachen Form der Wärmeleitungsgleichung auftretende Konduktivitätsmatrix (Steifigkeitsmatrix).

Lösung 3

a) Differentialgleichung und Randbedingungen

Differentialgleichung für den eindimensionalen Fall:

\frac{{\partial Q}}{{\partial x}}+r = 0\quad ,\quad Q = -\lambda \nabla \theta

\quad \Rightarrow \quad -\lambda \frac{{{\partial ^2}\theta }}{{\partial {x^2}}}+r = 0

Als Randbedingungen werden im Kontext der FEM die Dirichlet- und die Neumann-Randbedingung unterschieden. Für den eindimensionalen Fall des Wärmeleitstabes sind diese im Folgenden dargestellt.

Randbedingung 1. Art (Dirichlet)

dirichlet-randbedingung

\left. {\begin{array}{*{20}{c}} {\theta _1^{{e_1}} = \theta _1^{e_1^*}} \\ {\theta _1^{{e_2}} = \theta _1^{e_2^*}} \end{array}} \right\}\quad \forall {x_1} \in \Gamma _\theta ^e

Randbedingung 2. Art (Neumann):

neumann-randbedingung

\left. {\begin{array}{*{20}{c}} {Q_1^{{e_1}} = -Q_1^{e_1^*}} \\ {Q_1^{{e_2}} = Q_1^{e_2^*}} \end{array}} \right\}\quad \forall {x_1} \in \Gamma _Q^e

Es darf allerdings nur entweder eine RB 1. oder 2. Art gegeben sein.

Anfangsbedingungen werden nicht benötigt, da das Problem stationär ist.

b) Herleitung der schwachen Form

Für die Herleitung der schwachen Formulierung benötigen wir:

-\lambda {\theta ^{\prime \prime }}+r = 0

\left. {\begin{array}{*{20}{c}} {\theta _1^{{e_1}} = \theta _1^{e_1^*}} \\ {\theta _1^{{e_2}} = \theta _1^{e_2^*}} \end{array}} \right\}\quad \forall {x_1} \in \Gamma _\theta ^e

Als Testfunktion verwenden wir:

\delta {\theta _1}

Diese muss folgende Bedingungen erfüllen:

  • Sie muss die Dirichlet-Randbedingungen genügen
  • Sie muss den Feldgleichungen genügen, d.h. \delta {\gamma _{11}} = \delta {\theta _{1,1}}\quad \forall {x_1} \in {\Omega ^e}
  • Sie ist infinitesimal klein und ansonsten beliebig

Zur Herleitung der schwachen Formulierung benutzen wir die Variante mit der 1. Green’schen Identität:

0 = r-\lambda \Delta {\theta _1}\qquad \left| { \cdot \delta {\theta _1}} \right.\qquad \left| {\int {} } \right.

0 = \int\limits_{{\Omega ^e}} {\delta {\theta _1}rd{x_1}} -\int\limits_{{\Omega ^e}} {\lambda \delta {\theta _1}\underbrace {\Delta {\theta _1}}_{{\theta _{1,11}}}d{x_1}}

Mit der Green’schen Formel

\int\limits_\Omega {div\left( {w\nabla v} \right)dx} = \int\limits_\Omega {\left( {\nabla w \cdot \nabla v+w\Delta v} \right)dx} = \int\limits_{\partial \Omega } {w\frac{{\partial v}}{{\partial n}}ds} = \int\limits_{\partial \Omega } {w\left( {\nabla v \cdot \vec n} \right)dx}

\quad \Rightarrow \quad \int\limits_\Omega {\left( {w\Delta v} \right)dx = } \int\limits_{\partial \Omega } {w\left( {\nabla v \cdot \vec n} \right)ds} -\int\limits_\Omega {\left( {\nabla w \cdot \nabla v} \right)dx}

folgt:

\int\limits_{{\Omega ^e}} {\lambda \delta {\theta _1}{\theta _{1,11}}d{x_1}} = \int\limits_{{\Gamma ^e}} {\lambda \delta {\theta _1}\left( {{\theta _{1,1}} \cdot \vec n} \right)d\Gamma } -\int\limits_{{\Omega ^e}} {\lambda \delta {\theta _{1,1}}{\theta _{1,1}}d{x_1}}

\quad \Rightarrow \quad \int\limits_{{\Omega ^e}} {\lambda \delta {\theta _1}{\theta _{1,11}}d{x_1}} = \int\limits_{{\Gamma ^e}} {\delta {\theta _1}\underbrace {\lambda {\theta _{1,1}}}_{-Q}d\Gamma } -\int\limits_{{\Omega ^e}} {\lambda \underbrace {\delta {\theta _{1,1}}}_{\delta {\gamma _{11}}}\underbrace {{\theta _{1,1}}}_{{\gamma _{11}}}d{x_1}}

Im Eindimensionalen gilt: \vec n\;d\Gamma = d\Gamma

Einsetzen:

0 = \int\limits_{{\Omega ^e}} {\delta {\theta _1}rd{x_1}} -\int\limits_{{\Omega ^e}} {\lambda \delta {\theta _1}{\theta _{1,11}}d{x_1}}

\quad \Rightarrow \quad \int\limits_{{\Omega ^e}} {\delta {\gamma _{11}}\lambda {\gamma _{11}}d{x_1}} +\int\limits_{{\Omega ^e}} {\delta {\theta _1}rd{x_1}} = -\delta \theta _1^{{e_1}}Q_1^{e_1^*}-\delta \theta _1^{{e_2}}Q_1^{e_2^*}

c) Approximation

Wir wollen nun die unbekannten Größen der schwachen Form approximieren.

Allgemeine Bedingungen für die Ansatzfunktionen:

Die Ansatzfunktionen …

  • müssen konform sein (Interaktion und Kompatibilität benachbarter finiter Elemente)
  • müssen mindestens konstante Temperaturableitungen liefern
  • …

Approximationen:

{\theta _1}\left( {{\xi _1}} \right) \approx {\tilde \theta _1}\left( {{\xi _1}} \right) = \sum\limits_{i = 1}^{p+1} {{N^i}\left( {{\xi _1}} \right)\theta _1^{{e_i}}}

\delta {\theta _1}\left( {{\xi _1}} \right) \approx \delta {\tilde \theta _1}\left( {{\xi _1}} \right) = \sum\limits_{i = 1}^{p+1} {{N^i}\left( {{\xi _1}} \right)\delta \theta _1^{{e_i}}}

{\gamma _{11}}\left( {{\xi _1}} \right) \approx {{\tilde \gamma }_{11}}\left( {{\xi _1}} \right) = {{\tilde \theta }_{1,1}}\left( {{\xi _1}} \right) = \frac{{\partial {{\tilde \theta }_1}\left( {{\xi _1}\left( {{x_1}} \right)} \right)}}{{\partial {x_1}}} = \sum\limits_{i = 1}^{p+1} {\frac{{\partial {N^i}\left( {{\xi _1}} \right)}}{{\partial {x_1}}}\theta _1^{{e_i}}}

\quad \Rightarrow \quad {\gamma _{11}}\left( {{\xi _1}} \right) \approx {{\tilde \gamma }_{11}}\left( {{\xi _1}} \right) = \sum\limits_{i = 1}^{p+1} {N_{,1}^i\left( {{\xi _1}} \right)\theta _1^{{e_i}}}

\delta {\gamma _{11}}\left( {{\xi _1}} \right) \approx \delta {{\tilde \gamma }_{11}}\left( {{\xi _1}} \right) = \sum\limits_{i = 1}^{p+1} {N_{,1}^i\left( {{\xi _1}} \right)\delta \theta _1^{{e_i}}}

Für eine genaue Definition von {N^i} siehe Skript und die vorherigen Übungsaufgaben.

Unter isoparametrischem Konzept (isoparametrischer Darstellung) verstehen wir, dass für verschiedene Approximationen die gleichen Formfunktionen verwendet werden.

d) Approximierte schwache Form und Konduktivitätsmatrix

Approximierte schwache Form:

\int\limits_{{\Omega ^e}} {\delta {\gamma _{11}}\lambda {\gamma _{11}}d{x_1}} +\int\limits_{{\Omega ^e}} {\delta {\theta _1}rd{x_1}} = -\delta \theta _1^{{e_1}}Q_1^{e_1^*}-\delta \theta _1^{{e_2}}Q_1^{e_2^*}

\quad \Rightarrow \quad \int\limits_{-1}^1 {\delta {{\tilde \gamma }_{11}}\lambda {{\tilde \gamma }_{11}}\frac{L}{2}d{\xi _1}} +\int\limits_{-1}^1 {\delta {{\tilde \theta }_1}r\frac{L}{2}d{\xi _1}} = -\delta \tilde \theta _1^{{e_1}}Q_1^{e_1^*}-\delta \tilde \theta _1^{{e_2}}Q_1^{e_2^*}

Die Konduktivitätsmatrix bekommen wir wie folgt:

\int\limits_{-1}^1 {\delta {{\tilde \gamma }_{11}}\lambda {{\tilde \gamma }_{11}}\frac{L}{2}d{\xi _1}} = \sum\limits_{i = 1}^2 {\sum\limits_{j = 1}^2 {\delta \theta _1^{{e_i}}\theta _1^{{e_i}}\underbrace {\int\limits_{-1}^1 {\lambda N_{,1}^i\left( {{\xi _1}} \right)N_{,1}^j\left( {{\xi _1}} \right)\frac{L}{2}d{\xi _1}} }_{{k^{{e_{ij}}}}}} }

Für die Konduktivitätsmatrix folgt somit analog zu Übung 2 :

{\underline k ^e} = \frac{\lambda }{L}\left[ {\begin{array}{*{20}{c}} 1&{-1} \\ {-1}&1 \end{array}} \right]

\mathcal{S}\mathcal{W}\& \mathcal{J}\mathcal{K}

Leave a Reply

Close Menu