U 2 – Die statische Impulsbilanz – Mathematical Engineering

Das gegebene System besteht aus einem Stab der Länge $l = 3L$ und konstanter Dehnsteifigkeit $EA$ , der im mittleren Drittel durch die konstante Streckenlast ${p_0}$ belastet ist. Am linken Rand ist der Stab fest eingespannt und am rechten Rand wird die Verschiebungsrandbedingung $u_1^*$ aufgebracht.

stab

Aufgaben:

Diskretisieren Sie das System mit drei finiten Stabelementen. Stellen Sie die Element-Steifigkeitsmatrizen und die Elementlastvektoren unter Verwendung linearer Ansatzfunktionen auf.
Führen Sie unter Verwendung Ihrer Ergebnisse aus Aufgabenteil a) die Ensemblierung des Systems durch: Stellen Sie für das zu lösende Gleichungssystem die Systemsteifigkeitsmatrix, den Systemlastvektor und den Systemverschiebungsvektor auf.
Was verstehen Sie allgemein unter statischer Kondensation? Führen Sie die statische Kondensation für das System durch und bestimmen Sie die Lösung der entstehenden Gleichungen.
Sie verwenden zur Approximation den Polynomgrad $p = 1$ . Wie viele Gauß-Punkte würden Sie zur exakten Integration eines Lastvektors bei einer quadratischen Streckenlast bzw. im Falle eines transienten Problems bei der Massenmatrix benötigen? Begründen Sie Ihr Ergebnis kurz.

Gegeben: $EA$ , $L$ , ${p_0}$ , $u_1^* = 0,5\frac{{{p_0}{L^2}}}{{EA}}$

Lösung 2

a) Diskretisierung, Element-Steifigkeitsmatrizen und Elementlastvektoren

Diskretisierung mit drei finiten Stabelementen:

finite-elemente-diskretisierung-des-stabes

Zur Herleitung der schwachen Form nutzen wir ein differentielles Element des Stabes:

differentielles-stabelement

Kräftegleichgewicht:

$-{N_1}+{N_1}+{N_{1,1}}d{x_1}+{p_1}-\rho A{\ddot u_1} = 0$

mit dem Konstitutivgesetz

${N_{1,1}} = {\sigma _{11,1}}A = EA{\varepsilon _{11,1}}$

Differentialgleichung (Kinetik):

$0 = {\ddot u_1}\left( {{x_1}} \right)\rho A-{p_1}\left( {{x_1}} \right)-{N_{1,1}}d{x_1}$

Kinematik:

${\varepsilon _{11,1}}\left( {{x_1}} \right) = {u_{1,11}}\left( {{x_1}} \right)$

Testfunktion:

$\delta {u_1},\quad \delta {u_{1,1}} = \delta {\varepsilon _{11}}$

Neumann-Randbedingung:

$\left. {\begin{array}{*{20}{c}} {-N_1^{e_1^*}-N_1^{{e_1}} = 0} \\ {-N_1^{e_2^*}+N_1^{{e_2}} = 0} \end{array}} \right\}\quad \forall {x_1} \in {\Gamma _L}$

$0 = \int\limits_{{\Omega ^e}} {\delta {u_1}\left( {\rho A{{\ddot u}_1}-{p_1}} \right)d{x_1}} -\int\limits_{{\Omega ^e}} {\delta {u_1}\underbrace {{N_{1,1}}}_{EA{u_{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}} {\delta {u_1}EA{u_{1,11}}d{x_1}} = \int\limits_{{\Gamma ^e}} {\delta {u_1}EA{u_{1,1}}d\Gamma } -\int\limits_{{\Omega ^e}} {\delta {u_{1,1}}EA{u_{1,1}}d{x_1}}$

Durch Einsetzen ergibt sich:

$\int\limits_{{\Omega ^e}} {\delta {u_1}\rho A{{\ddot u}_1}d{x_1}} -\int\limits_{{\Omega ^e}} {\delta {u_1}{p_1}d{x_1}} +\int\limits_{{\Omega ^e}} {\underbrace {\delta {u_{1,1}}}_{\delta {\varepsilon _{11}}}EA\underbrace {{u_{1,1}}}_{{\varepsilon _{11}}}d{x_1}} = \int\limits_{{\Gamma ^e}} {\delta {u_1}\underbrace {EA{u_{1,1}}}_{{N_1}}d\Gamma }$

$\quad \Rightarrow \quad \int\limits_{{\Omega ^e}} {\delta {u_1}\rho A\ddot ud{x_1}} +\int\limits_{{\Omega ^e}} {\delta {\varepsilon _{11}}EA{\varepsilon _{11}}d{x_1}} = \delta u_1^{{e_2}}N_1^{{e_2}}-\delta u_1^{{e_1}}N_1^{{e_1}}+\int\limits_{{\Omega ^e}} {\delta {u_1}{p_1}d{x_1}}$

Wir setzen nun die Neumann-Randbedingung ein:

$\underbrace {\int\limits_{{\Omega ^e}} {\delta {u_1}\rho A\ddot ud{x_1}} }_{\delta W_{dyn}^e}+\underbrace {\int\limits_{{\Omega ^e}} {\delta {\varepsilon _{11}}EA{\varepsilon _{11}}d{x_1}} }_{\delta W_{int}^e \Rightarrow {{\underline k }^e}} = \underbrace {\delta u_1^{e_2^*}N_1^{e_2^*}+\delta u_1^{e_1^*}N_1^{e_1^*}+\int\limits_{{\Omega ^e}} {\delta {u_1}{p_1}d{x_1}} }_{\Omega _{ext}^e}$

Jetzt wollen wir die hergeleitete schwache Form approximieren:

$\int\limits_{{\Omega ^e}} {\delta {{\tilde u}_1}\rho A\ddot {\tilde u}d{x_1}} +\int\limits_{{\Omega ^e}} {\delta {{\tilde \varepsilon }_{11}}EA{{\tilde \varepsilon }_{11}}d{x_1}} = \delta \tilde u_1^{e_2^*}N_1^{e_2^*}+\delta \tilde u_1^{e_1^*}N_1^{e_1^*}+\int\limits_{{\Omega ^e}} {\delta {{\tilde u}_1}{p_1}d{x_1}}$

Jetzt approximieren wir die das Verschiebungsfeld, die Verzerrungen, sowie deren Variationen:

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

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

${\varepsilon _{11}}\left( {{\xi _1}} \right) \approx {{\tilde \varepsilon }_{11}}\left( {{\xi _1}} \right) = \frac{{\partial {{\tilde u}_1}\left( {{\xi _1}} \right)}}{{\partial {x_1}}} = \frac{{\partial \sum\limits_{i = 1}^{p+1} {{N^i}\left( {{\xi _1}} \right)u_1^{{e_i}}} }}{{\partial {x_1}}} = \sum\limits_{i = 1}^{p+1} {\frac{{\partial {N^i}\left( {{\xi _1}} \right)}}{{\partial {x_1}}}u_1^{{e_i}}}$

$\quad \Rightarrow \quad {{\tilde \varepsilon }_{11}}\left( {{\xi _1}} \right) = \sum\limits_{i = 1}^{p+1} {\partial N_{,1}^i\left( {{\xi _1}} \right)u_1^{{e_i}}}$

Mit der inneren virtuellen Arbeit bestimmen wir die Elementsteifigkeitsmatrix:

$\delta W_{int}^{{e}} = \sum\limits_{i = 1}^2 {\sum\limits_{j = 1}^2 {\int\limits_{-1}^1 {\left[ {\delta u_1^{{e_i}}N_{,1}^i} \right]EA\left[ {u_1^{{e_j}}N_{,1}^j} \right]\frac{L}{2}d{\xi _1}} } }$

$\quad = \sum\limits_{i = 1}^2 {\sum\limits_{j = 1}^2 {\delta u_1^{{e_i}}u_1^{{e_j}}\underbrace {EA\frac{L}{2}\int\limits_{-1}^1 {N_{,1}^iN_{,1}^jd{\xi _1}} }_{{{\underline k }^{{e_{ij}}}}}} }$

Wir verwenden nun wieder die linearen Ansatzfunktionen, die bereits aus Übung 1 bekannt sind:

${N^1} = \frac{1}{2}\left( {1-{\xi _1}} \right),\quad N_{,1}^1 = \frac{{\partial {N^i}}}{{\partial {\xi _1}}}\frac{{\partial {\xi _1}}}{{\partial {x_1}}} = -\frac{1}{L}$

${N^2} = \frac{1}{2}\left( {1+{\xi _1}} \right),\quad N_{,1}^2 = \frac{{\partial {N^i}}}{{\partial {\xi _1}}}\frac{{\partial {\xi _1}}}{{\partial {x_1}}} = +\frac{1}{L}$

Für die Elementsteifigkeitsmatrix ergibt sich damit:

${\underline k ^e} = \frac{{EA}}{L}\left[ {\begin{array}{*{20}{c}} 1&{-1} \\ {-1}&1 \end{array}} \right] = {\underline k ^{{e_1}}} = {\underline k ^{{e_2}}} = {\underline k ^{{e_3}}}$

Da alle Elemente die gleiche Geometrie und Dehnsteifigkeit haben, sind die Elementsteifigkeitsmatrizen ebenfalls gleich.

Wir stellen nun den Elementlastvektor auf:

$\int\limits_{{\Omega ^e}} {\delta {{\tilde u}_1}{p_1}d{x_1}} = \int\limits_{-1}^1 {\left[ {\sum\limits_{i = 1}^{p+1} {{N^i}\left( {{\xi _1}} \right)u_1^{{e_i}}} } \right]{p_1}\frac{L}{2}d{\xi _1}} = \sum\limits_{i = 1}^{p+1} {u_1^{{e_i}}\underbrace {\int\limits_{-1}^1 {{N^i}\left( {{\xi _1}} \right){p_1}\frac{L}{2}d{\xi _1}} }_{{r^{{e_i}}}}}$

$r_p^{{e_i}} = \int\limits_{-1}^1 {{N^i}{p_1}\frac{L}{2}d{\xi _1}}$

$r_p^{{e_1}} = {p_1}\frac{L}{2}\int\limits_{-1}^1 {\frac{1}{2}\left( {1-{\xi _1}} \right)d{\xi _1}} = {p_1}\frac{L}{4}\left[ {{\xi _1}-\frac{{\xi _1^2}}{2}} \right]_{-1}^1 = {p_1}\frac{L}{2}$

$r_p^{{e_2}} = {p_1}\frac{L}{2}\int\limits_{-1}^1 {\frac{1}{2}\left( {1+{\xi _1}} \right)d{\xi _1}} = {p_1}\frac{L}{4}\left[ {{\xi _1}+\frac{{\xi _1^2}}{2}} \right]_{-1}^1 = {p_1}\frac{L}{2}$

$\quad \Rightarrow \quad \vec r_p^e = \left\{ {\begin{array}{*{20}{c}} {r_p^{{e_1}}} \\ {r_p^{{e_2}}} \end{array}} \right\} = \frac{L}{2}\left\{ {\begin{array}{*{20}{c}} {{p_1}} \\ {{p_1}} \end{array}} \right\}$

Für unsere konkreten Elemente gilt also:

$\vec r_p^2 = \frac{{{p_0}L}}{2}\left\{ {\begin{array}{*{20}{c}} 1 \\ 1 \end{array}} \right\}\quad ;\quad \vec r_p^1 = \vec r_p^3 = \left\{ {\begin{array}{*{20}{c}} 0 \\ 0 \end{array}} \right\}$

b) Ensemblierung und Aufstellen des Gleichungssystems

Systemsteifigkeitsmatrix:

$\underline K = \frac{{EA}}{L}\left[ {\begin{array}{*{20}{c}} 1&{-1}&0&0 \\ {-1}&{1+1}&{-1}&0 \\ 0&{-1}&{1+1}&{-1} \\ 0&0&{-1}&1 \end{array}} \right]$

${{\vec r}_p} = \frac{{{p_0}L}}{2}\left\{ {\begin{array}{*{20}{c}} 0 \\ 1 \\ 1 \\ 0 \end{array}} \right\},\quad {{\vec r}_u} = \left\{ {\begin{array}{*{20}{c}} {r_u^1} \\ 0 \\ 0 \\ {r_u^4} \end{array}} \right\}\quad \Rightarrow \quad \vec r = {{\vec r}_p}+{{\vec r}_u} = \left\{ {\begin{array}{*{20}{c}} {r_u^1} \\ {\frac{{{p_0}L}}{2}} \\ {\frac{{{p_0}L}}{2}} \\ {r_u^4} \end{array}} \right\}$

Das zu lösende Gleichungssystem mit eingesetzten Randbedingungen lautet also:

$\underline K {\vec u^e} = \vec r\quad \Rightarrow \quad \frac{{EA}}{L}\left[ {\begin{array}{*{20}{c}} 1&{-1}&0&0 \\ {-1}&2&{-1}&0 \\ 0&{-1}&2&{-1} \\ 0&0&{-1}&1 \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} 0 \\ {u_1^2} \\ {u_1^3} \\ {u_1^*} \end{array}} \right\} = \left\{ {\begin{array}{*{20}{c}} {r_u^1} \\ {\frac{{{p_0}L}}{2}} \\ {\frac{{{p_0}L}}{2}} \\ {r_u^4} \end{array}} \right\}$

c) Statische Kondensation

Unter statischer Kondensation versteht man die Vereinfachung des Gleichungssystems durch Umsortierung nach bekannten und unbekannten Größen. (siehe auch: hier)

In unserem Fall ergibt sich das vereinfachte System:

$\left[ {\begin{array}{*{20}{c}} {{{\underline K }_{uu}}}&{{{\underline K }_{up}}} \\ {{{\underline K }_{pu}}}&{{{\underline K }_{pp}}} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} {{{\vec u}_u}} \\ {{{\vec u}_p}} \end{array}} \right\} = \left\{ {\begin{array}{*{20}{c}} {{{\vec r}_u}} \\ {{{\vec r}_p}} \end{array}} \right\}$

Durch Spaltentausch wird aus der Matrix:

$\frac{{EA}}{L}\left[ {\begin{array}{*{20}{c}} 1&0&{-1}&0 \\ {-1}&0&2&{-1} \\ 0&{-1}&{-1}&2 \\ 0&1&0&{-1} \end{array}} \right]$

Durch Zeilentausch erhalten wir dann das LGS:

$\frac{{EA}}{L}\left[ {\begin{array}{*{20}{c}} 1&0&\vline & {-1}&0 \\ 0&1&\vline & 0&{-1} \\ \hline {-1}&0&\vline & 2&{-1} \\ 0&{-1}&\vline & {-1}&2 \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} 0 \\ {u_1^*} \\ \hline {u_1^2} \\ {u_1^3} \end{array}} \right\} = \left\{ {\begin{array}{*{20}{c}} {r_u^1} \\ {r_u^4} \\ \hline {\frac{{{p_0}L}}{2}} \\ {\frac{{{p_0}L}}{2}} \end{array}} \right\}$

Zur Lösung schreiben wir:

${\underline K _{uu}}{{\vec u}_u}+{\underline K _{up}}{{\vec u}_p} = {{\vec r}_u}\quad \quad \quad \left( 2 \right)$

${\underline K _{pu}}{{\vec u}_u}+{\underline K _{pp}}{{\vec u}_p} = {{\vec r}_p}\quad \quad \quad \left( 1 \right)$

Aus (1) folgt:

${\vec u_p} = \underline K _{pp}^{-1}\left\{ {{{\vec r}_p}-{{\underline K }_{pu}}{{\vec u}_u}} \right\}$

Für die Inverse gilt:

${\left[ {\begin{array}{*{20}{c}} a&b \\ c&d \end{array}} \right]^{-1}} = \frac{1}{{\det }}\left[ {\begin{array}{*{20}{c}} d&{-b} \\ {-c}&a \end{array}} \right]$

Damit folgt:

$\underline K _{pp}^{-1} = \frac{1}{3}{\left( {\frac{L}{{EA}}} \right)^2}\left[ {\begin{array}{*{20}{c}} 2&1 \\ 1&2 \end{array}} \right]\frac{{EA}}{L} = \frac{L}{{3EA}}\left[ {\begin{array}{*{20}{c}} 2&1 \\ 1&2 \end{array}} \right]$

$\Rightarrow \quad {{\vec u}_p} = \frac{L}{{3EA}}\left[ {\begin{array}{*{20}{c}} 2&1 \\ 1&2 \end{array}} \right]\left\{ {\left\{ {\begin{array}{*{20}{c}} {\frac{{{p_0}L}}{2}} \\ {\frac{{{p_0}L}}{2}} \end{array}} \right\}-\frac{{EA}}{L}\left[ {\begin{array}{*{20}{c}} {-1}&0 \\ 0&{-1} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} 0 \\ {u_1^*} \end{array}} \right\}} \right\}$

$\quad \quad \quad = \frac{L}{{3EA}}\left[ {\begin{array}{*{20}{c}} 2&1 \\ 1&2 \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} {\frac{{{p_0}L}}{2}} \\ {\frac{{{p_0}L}}{2}+\frac{{EA}}{L}u_1^*} \end{array}} \right\}$

$\quad \quad \quad = \frac{L}{{3EA}}\left\{ {\begin{array}{*{20}{c}} {\frac{{3{p_0}L}}{2}+\frac{{EA}}{L}u_1^*} \\ {\frac{{3{p_0}L}}{2}+2\frac{{EA}}{L}u_1^*} \end{array}} \right\}$

In der Aufgabenstellung ist gegeben:

$u_1^* = \frac{{{p_0}{L^2}}}{{2EA}}$

Damit folgt:

${\vec u_p} = \left\{ {\begin{array}{*{20}{c}} {u_1^2} \\ {u_1^3} \end{array}} \right\} = \frac{{{p_0}{L^2}}}{{6EA}}\left\{ {\begin{array}{*{20}{c}} 4 \\ 5 \end{array}} \right\}$

Aus (2) folgt schließlich:

${{\vec r}_u} = {\underline K _{uu}}{{\vec u}_u}+{\underline K _{up}}{{\vec u}_p} = \left[ {\begin{array}{*{20}{c}} 1&0 \\ 0&1 \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} 0 \\ {u_1^*} \end{array}} \right\}+\left[ {\begin{array}{*{20}{c}} {-1}&0 \\ 0&{-1} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} {u_1^2} \\ {u_1^3} \end{array}} \right\}$

$\quad \Rightarrow \quad {{\vec r}_u} = \left[ {\begin{array}{*{20}{c}} 1&0 \\ 0&1 \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} 0 \\ {\frac{{{p_0}{L^2}}}{{2EA}}} \end{array}} \right\}+\frac{{{p_0}{L^2}}}{{6EA}}\left[ {\begin{array}{*{20}{c}} {-1}&0 \\ 0&{-1} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} 4 \\ 5 \end{array}} \right\}$

$\quad \Rightarrow \quad {{\vec r}_u} = -\frac{{{p_0}{L^2}}}{{EA}}\left\{ {\begin{array}{*{20}{c}} {\frac{2}{3}} \\ {\frac{1}{3}} \end{array}} \right\}$

d) Exakte Approximation

Im Falle einer quadratischen Streckenlast würde für den Lastvektor gelten:

${r^{{e_i}}} = \frac{L}{2}\int\limits_{-1}^1 {\underbrace {{N^i}}_{p = 1}\underbrace {{p_1}}_{p = 2}d{\xi _1}} \quad \Rightarrow \quad p = 3$

Hierbei steht p für den jeweiligen Polynomgrad.

Für die Anzahl der benötigten Gaußpunkte gilt:

$p \leq 2 \cdot {N_G}-1\quad \Rightarrow \quad {N_G} \geq \frac{{p+1}}{2} = \frac{{3+1}}{2} = 2$

Wir brächten zur exakten Integration des Lastvektors mit linearen Ansatzfunktion und quadratischer Strecken also 2 Gauß-Punkte.

Im Falle eines transienten Problems benötigen wir zur exakten Integration der Massenmatrix:

$\delta W_{dyn}^e = \int\limits_{{\Omega ^e}} {\delta {u_1}\rho A{{\ddot u}_1}dx}$

$\quad \Rightarrow \quad \delta W_{dyn}^e = \int\limits_{{\Omega ^e}} {\delta {{\tilde u}_1}{N_1}\rho A{{\ddot {\tilde u}}_1}{N_1}\frac{L}{2}d{\xi _1}}$

$\quad \Rightarrow \quad \delta W_{dyn}^e = \sum\limits_{i = 1}^{p+1} {\sum\limits_{j = 1}^{p+1} {\delta u_1^{{e_i}}\ddot u_1^{{e_j}}\underbrace {\int\limits_{-1}^1 {{N^i}\left( {{\xi _1}} \right)\rho A{N^j}\left( {{\xi _1}} \right)\frac{L}{2}d{\xi _1}} }_{Elementmasse\;\;{m^{{e_{ij}}}}}} }$

$\quad \Rightarrow \quad {m^{{e_{ij}}}} = \int\limits_{-1}^1 {\underbrace {{N^i}\left( {{\xi _1}} \right)}_{p = 1}\rho A\underbrace {{N^j}\left( {{\xi _1}} \right)}_{p = 1}\frac{L}{2}d{\xi _1}} \quad \Rightarrow \quad p = 2$

Für die Anzahl der benötigten Gaußpunkte gilt hier:

${N_G} \geq \frac{{p+1}}{2} = \frac{{2+1}}{2} = \frac{3}{2} = 1,5\quad \Rightarrow \quad {N_G} \geq 2$

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

Lösung 2

a) Diskretisierung, Element-Steifigkeitsmatrizen und Elementlastvektoren

b) Ensemblierung und Aufstellen des Gleichungssystems

c) Statische Kondensation

d) Exakte Approximation

You Might Also Like

Aufgabe 1 – Analytische Lösungen für verschiedene Lastfälle

1 – Einführung

U 3 – Die stationäre Wärmeleitungsgleichung

Leave a Reply Cancel reply