Gelöste Aufgaben/FEAG

Aufgabenstellung

Analog zu FEAF untersuchen wir hier die Schwingungen eines Kontinuums beim Loslassen aus der entspannten Rugelage. Hier nicht mit einem Dehnstab, sondern einem Euler-Bernoulli-Balken.

Lageplan

Gesucht ist die Schwingung eines Euler-Bernoulli-Balkens beim Loslassen aus der Ruhelage. Wir gehen nach dem Standardrezept der Finite Elemente Methode vor, arbeiten also mit dem Sources/Lexikon/Prinzip der virtuellen Arbeit.

Lösung mit Matlab

Interessant ist hier, dass - im Gegensatz zu Stablängsschwingungen - die Eigenfrequenz nicht ein gerades Vielfaches der untersten Eigenfrequenz ist. Falls Sie ein Saiteninstrument spielen, verstehen Sie sofort, warum das wichtig ist.

Maxima können wir hier nicht gut gebrauchen: die Gleichungen werden zu umfangreich. Wir arbeiten also mehr mit numerischen Verfahren, da ist Matlab geeigneter. Allerdings können wir Matlab-Inhalte nicht gut auf dieser Seite unterbringen - deshalb gibt es dafür die Seite FEAG-Matlab, die der gleichen Struktur folgt.

Header

Im Programm arbeiten wir mit einer dimensionslosen Formulierung - wir brauchen dafür eine Bezugszeit t_Bez und eine Bezugslänge l_Bez.

Dazu nehmen wir eine "Anleihe" bei der analytischen Lösung des Schwingungsproblems (vgl. Aufgabe SKEB):

Analytische Lösung: homohener Lösungsanteil

Analytische Lösung: partikularer Lösungsanteil

Die homogene Bewegungsgleichung des Euler-Bernoulli-Balkens $\varrho \,A\cdot {\ddot {w}}+E\,I\cdot {w}^{IV}=0$ hat Lösungen vom Typ $\displaystyle w(x,t)=\sum _{i=1}^{\infty }W_{i}\cdot e^{\displaystyle \kappa _{i}\cdot x}\cdot e^{\displaystyle \omega _{0,i}\cdot t}$ Wir bekommen zu jedem ω0 vier κ ${\begin{array}{l}\displaystyle \kappa _{1}=+{\frac {j\,\alpha }{\ell }},\\\displaystyle \kappa _{2}=-{\frac {\alpha }{\ell }},\\\displaystyle \kappa _{3}=-{\frac {j\,\alpha }{\ell }},\\\displaystyle \kappa _{4}=+{\frac {\alpha }{\ell }}\end{array}}\;\;{\text{ mit }}\;\;\displaystyle \alpha ^{4}={\frac {{{\ell }^{4}}\,A\,\rho }{E\,I}}\cdot {{\omega }_{0}^{2}}$ Das α ist eine praktische Abkürzung, hinter der wir ω0 verstecken. Anders als beim Dehnstab (SKER) finden wir hier keine analytische Beziehung, sondern nur die numerischen Beziehungen, für die unsere Randbedingungen erfüllt sind: ${\begin{array}{c}\alpha _{1}=1.875\ldots ,\\\alpha _{2}=4.694\ldots ,\\\alpha _{3}=7.854\ldots ,\\\vdots \end{array}}$ Die langsamste Eigenmode gehört zu α1 mit der Periodendauer $\displaystyle T^{*}=0.5688\ldots \pi \ell ^{2}{\sqrt {\frac {\varrho \,A}{E\,I}}}$. Also wählen wir $t_{Bez}:=T^{*}$.

Die zugeordnete inhomogene Bewegungsgleichung des Euler-Bernoulli-Balkens $E\,I\cdot {w}^{IV}=\varrho \,A\cdot g$ hat die partikulare Lösung $\displaystyle E\,I\,{w}_{p}\left(x\right)=A\,g\,\varrho \,\ell ^{4}\cdot \left({\frac {{\xi }^{4}}{24}}-{\frac {{\xi }^{3}}{6}}+{\frac {{\xi }^{2}}{4}}\right)$ Die statische Auslenkung am unteren Ende ist demnach $\displaystyle {{w}_{s}}={\frac {{{\ell }^{4}}\,A\,g\,\rho }{8\,E\,I}}$. Wir wählen \ell_{Bez} := w_s.

System-Parameter des FEM-Modells

Für die Diskretisierung wählen wir als Anzahl der Finiten Elemente (Number Of Elements):

$\displaystyle I=4\;\;\;({\text{ im Matlab-Code: }}NOE)$

und damit je Element

$\displaystyle \ell _{i}={\frac {\ell _{0}}{I}}$.

Wir wählen die kubische Ansatzfunktionen aus dem Abschnitt Finite Elemente Methode je Element, also

${\begin{array}{crc}\phi _{1}=&\left(\xi -1\right)^{2}\cdot \left(1+2\cdot \xi \right)\\\phi _{2}=&\left(\xi -1\right)^{2}\cdot \xi \cdot \ell _{i}\\\phi _{3}=&\left(\xi ^{2}\cdot \left(3-2\cdot \xi \right)\right)\\\phi _{4}=&\left(\left(\xi -1\right)\cdot \xi ^{2}\cdot \ell _{i}\right)\end{array}}$

Für die numerisch Implementierung stört in dieser Darstellung das Element-spezifische ℓ_i - die Elementlänge. Wir behelfen uns, indem wir die Ansatzfunktionen schreiben als

$\displaystyle \underbrace {\left({\begin{array}{c}\phi _{1}\\\phi _{2}\\\phi _{3}\\\phi _{4}\end{array}}\right)} _{\displaystyle =:{\underline {\phi }}}=\underbrace {\left({\begin{array}{cccc}1&0&0&0\\0&\ell _{i}&0&0\\0&0&1&0\\0&0&0&\ell _{i}\\\end{array}}\right)} _{\displaystyle =:{\underline {\underline {d}}}}\cdot \underbrace {\left({\begin{array}{c}\left(\xi -1\right)^{2}\cdot \left(1+2\cdot \xi \right)\\\left(\xi -1\right)^{2}\cdot \xi \\\left(\xi ^{2}\cdot \left(3-2\cdot \xi \right)\right)\\\left(\left(\xi -1\right)\cdot \xi ^{2}\right)\end{array}}\right)} _{\displaystyle =:{\underline {\tilde {\phi }}}}$,

wobei wir die Ansatz-Polynome in ${\underline {\tilde {\phi }}}$ verpacken.

Das schaut unnötig komplex aus - allerdings stecken wir die Diagonal-Matrizen d (für jedes Element eine) in eine Matlab-Variable, so dass sie dort nicht weiter auffällt.

Die abhängigen Koordinaten des FEM-Modells sind dann zunächst - bis zum Einarbeiten der Randbedingungen - diese:

${\underline {Q}}_{t}=\left({\begin{array}{c}W_{0}(t)\\\Phi _{0}(t)\\W_{1}(t)\\\vdots \\\Phi _{I}(t)\\W_{I}(t)\\\end{array}}\right)$.

Equilibrium Conditions

xffsd

Die Gleichgewichtsbedingung mit dem Prinzip der virtuellen Verrückungen setzen sich additiv aus

$\displaystyle \delta W=\sum _{i=1}^{I}\left(\delta W_{i}^{a}-\delta \Pi _{i}\right)+\delta W_{R}^{a}$,

zusammen - also den virtuellen Arbeiten je Element zuzüglich von virtuellen Arbeiten am Rand. In unserem Beispiel haben wir allerdings keine eingeprägten, äußeren Lasten am Rand, also ist

$\delta W_{R}^{a}=0$.

Je Element erfassen wir die virtuellen Arbeiten durch Element-Matrizen. So ist z.B. für die virtuelle Formänderungsenergie (vgl. Finite Elemente Methode)

$\delta \Pi _{i}=\left(\delta W_{i-1},\delta \Phi _{i-1},\delta W_{i},\delta \Phi _{i}\right)\cdot {\underline {\underline {K}}}_{i}\cdot \left({\begin{array}{c}W_{i-1}(t)\\\Phi _{i-1}(t)\\W_{i}(t)\\\Phi _{i}(t)\end{array}}\right)$

mit der Element-Steifigkeitsmatrix ${\underline {\underline {K}}}_{i}$.

Aus der FEM-Formulierung für den Euler-Bernoulli-Balken könnten wir die Element-Matrizen der virtuellen Arbeit der D'Alembert'schen Trägheitskraft sowie der Virtuelle Formänderungsenergie herauskopieren. Auch die "rechte Seite" mit der virtuelle Arbeit der Gewichtskraft aus

$\displaystyle \delta W_{i}^{a,g}=\int _{\ell _{i}}\varrho \;A\;g\cdot \delta w(x)\;dx_{i}$,

könnten wir als Element-Lastmatrix (Spaltenmatrix) abspeichern.

$\displaystyle {\underline {G}}_{i}={\frac {\varrho \,A\,g\,\ell _{i}}{2}}\cdot \left({\begin{array}{c}+1\\\displaystyle +{\frac {\ell _{i}}{6}}\\+1\\\displaystyle -{\frac {\ell _{i}}{6}}\end{array}}\right)$.

Dass es auch anders geht,sehen Sie in der Definition der Matlab-Klasse FEM_Element_Modell, bei der die Element direkt aus den Trial-Functions - hier den

${\underline {\tilde {\phi }}}$ -

hergeleitet. werden:

: K(row,col) = diff(... polyval(... polyint(... conv(polyder(self.phi(row,:)),polyder(self.phi(col,:)))),... [0,1]));  :

Komponieren der System-Matrizen

Die System-Matrizen des Gesamt–Systems komponieren wir nun durch Hinzuaddieren der Anteile je Element.

Beispiel: die Gesamt-Steifigkeitsmatrix für I=2:

$\displaystyle {\underline {\underline {K}}}={\frac {E\,I}{\ell _{i}^{3}}}\cdot \left({\begin{array}{ccccc}{\color {red}{+12}}&{\color {red}{+6}}&{\color {red}{-12}}&{\color {red}{+6}}&0&0\\{\color {red}{+6\ell _{i}}}&{\color {red}{+4\ell _{i}^{2}}}&{\color {red}{-6\ell _{i}}}&{\color {red}{+2\ell _{i}^{2}}}&0&0\\{\color {red}{-12}}&{\color {red}{-6}}&{\color {red}{+12}}{\color {green}{+12}}&{\color {red}{-6}}{\color {green}{+6}}&{\color {green}{-12}}&{\color {green}{+6}}\\{\color {red}{+6\ell _{i}}}&{\color {red}{2\ell _{i}^{2}}}&{\color {red}{-6\ell _{i}}}{\color {green}{+6\ell _{i}}}&{\color {red}{+2\ell _{i}^{2}}}{\color {green}{+2\ell _{i}^{2}}}&{\color {green}{-6\ell _{i}}}&{\color {green}{+2\ell _{i}^{2}}}\\0&0&{\color {green}{-12}}&{\color {green}{-6}}&{\color {green}{+12}}&{\color {green}{-6}}\\0&0&{\color {green}{+6\ell _{i}}}&{\color {green}{2\ell _{i}^{2}}}&{\color {green}{-6\ell _{i}}}&{\color {green}{+2\ell _{i}^{2}}}\end{array}}\right)$

Die Beiträge der zwei Elemente sind hier in rot bzw. grün eingefärbt.

Dieses "Hinzuaddieren" passiert in Matlab hier:

: e = 0;% counter for all elements of all sections for section = 1: length(element.noe)  : for i = 1:element.noe(section) e = e+1;  : j = 2*(e-1)+1;  : K(j:j+3,j:j+3) = K(j:j+3,j:j+3)+b* d(:,:,e)*element.K*d(:,:,e);  : end  :

In Matrix-Schreibweise lautet die Bewegungsgleichung des Gesamt-Systems jetzt:

${\underline {\underline {M}}}\cdot {\underline {\ddot {Q}}}+{\underline {\underline {K}}}\cdot {\underline {Q}}={\underline {G}}$.

Einarbeiten der Randbedingungen

Die Randbedingungen arbeiten wir ein, indem wir zeilenweise (für δW und δΦ) und spaltenweise (für W und Φ) streichen.

Hier sind es

• die erste Zeile/Spalte für W₀ (blau) und <span style="color:#009999"> This is an example...</span>
• die zweite Zeile / Spalte für Φ₀ (grün).

Solving

xffsd

Adapt to Initial Conditions

xffsd

Post-Processing

xffsd

Header

xffsd

Header

xffsd

Links

• ...

Literature

• ...

Particulare Lösung.

Animation der Lösung

Excel-File mit System-Parametern.

Excel-Parameterdatei: Randbedingungen

Parameter des Modells

Randbedingungen einarbeiten.

Lösung im Zeitbereich.

Mode #1

Mode #3

Mode #5

Modalformen

Animation in Time Domain.