Gelöste Aufgaben/FEB4

Aufgabenstellung

Wie in FEB3 geht es hier um die Ausnutzung von Symmetrien mit einer skalaren Gleichgewichtsformulierung.

Lageplan

Das skizzierte System besteht aus Euler-Bernoulli-Balken (jeweils Länge a) der Biegesteifigkeit EI. Die Konstruktion wird in A durch ein gelenkiges Festlager gehalten und in B durch die Kraft F belastet.

Gesucht ist die Ersatzfeder-Steifigkeit des Systems, das Sie aus Finiten Elementen mit Euler-Bernoulli-Balken zusammensetzen.

Gesucht ist die Auslenkung des Punktes B. Das Modell soll aus wenigen Finiten Elementen bestehen, überlegen Sie, welche Modelle sinnvoll sind. Die neutralen Fasern der Balken seien in Längsrichtung undehnbar – die Querschnitts-Schwerpunkte verschieben sich also nicht in Balken-Längsrichtung.

Gegeben: a, E I, F

Lösung mit Maxima

Header

Für das FE-Modell wollen wir nicht alle Sektionen der Feder modellieren - wir schneiden die Feder entlang ihrer Längsachse durch und nehmen eins der so entstandenen Elemente exemplarisch heraus.

Federelement.

Von diesen Elementen haben wir fünf in der Feder - die Längung dieses Elements der Feder entspricht also 1/5 der Feder-Längung in B.

Im Schnitt der horizontalen Sektionen mit der Längsachse der Feder ist das Schnittmoment jeweils Null - wir dürfen uns also an den Rändern unseres Sub-Systems ein Gelenk "denken" und mit dem hier gezeichneten System arbeiten.

Aber Achtung: so wie gezeichnet ist es statisch unbestimmt: wir können im Gleichgewicht beliebige Verschiebungen des Systems in Feder-Längsrichtung angeben.

✔ Geht's noch einfacher?:

* Wie müsste ein System aussehen, das noch einfacher gestrickt ist? Welche neuen Randbedingungen ergeben sich?

/*******************************************************/
/* MAXIMA script                                       */
/* version: wxMaxima 16.04.2                           */
/* author: Andreas Baumgart                            */
/* last updated: 2020-01-10                            */
/* ref: TM-C                                           */
/* description: FEM-solution for spring section        */
/*******************************************************/

tmp

Declarations

Wir haben also ein Ersatz-System mit drei Finiten Elementen und vier Knoten:

Knoten und Stäbe.

Die Knoten 0 und 3 entstehen in unserem Ersatz-System, weil hier das Biegemoment in den horizontalen Sektionen (I und III) Null ist. Diese beiden Sektionen sind nur a/2 lang, die vertikale Sektion ist a lang.

Damit das System in Längsrichtung festgelegt ist, führen wir in Knoten "0" ein gelenkiges Festlager ein - die Verschiebung in Knoten "3" ist dann die Auslenkung von Knoten "3" relativ zu Knoten "0".

Dementsprechnd haben die drei Finiten Elementen die Längen

${\begin{array}{l}\ell _{1}=a/2\\\ell _{2}=a\\\ell _{3}=a/2\\\end{array}}$ .

Die Element-Steifigkeitsmatrix kopieren wir aus Abschnitt FEM-Formulierung für den Euler-Bernoulli-Balken zu:

$\displaystyle K_{i}(\ell _{i})={\frac {EI}{\ell _{i}^{3}}}\cdot {\begin{pmatrix}12&6\cdot \ell &-12&6\cdot \ell _{i}\\6\cdot \ell _{i}&4\cdot {{\ell _{i}}^{2}}&-6\cdot \ell _{i}&2\cdot {{\ell _{i}}^{2}}\\-12&-6\cdot \ell _{i}&12&-6\cdot \ell _{i}\\6\cdot \ell _{i}&2\cdot {{\ell _{i}}^{2}}&-6\cdot \ell _{i}&4\cdot {{\ell _{i}}^{2}}\end{pmatrix}}$ .

/*******************************************************/
/* coordinates and boundary conditions                 */
 l   : [a/2,a,a/2];
norm : ŵ = F*a^3/(3*EI) * 5;

 X : append(makelist(Q[i]= [W[i-1], Φ[i-1], W[i], Φ[i]],i,1,3),[Q[0]= [U[1], Ψ[1], Ψ[2], U[3]]]);
δX : subst([Q=δQ,W=δW,Φ=δΦ],X);

/* boundary conditions */
kinematics: [[W[0]=  0  ,        δW[0]=  0           ],
             [W[1]=  0  , W[2]=0,δW[1]=  0   ,δW[2]=0],
             [W[2]=-W[1],        δW[2]=-δW[1]        ]];

tmp

Für jedes Finite Element führen wir zuerst das standardisierte lokale x/z-Koordinatensystem ein. Und wir definieren neue, anschauliche - globale - Koordinaten des Systems, um sie von den Standard-Koordinaten W_i , Φ_i unterscheiden zu können:

Lokale Orts-Koordinaten xi, zi	Globale Knoten-Koordinaten Ui, Ψi

Wir haben also je Element die eingezeichneten globalen Knoten-Koordinaten

$U_{i-1},\Psi _{i-1},U_{i},\Psi _{i}$

sowie die lokalen FE-Koordinaten  

$W_{i-1},\Phi _{i-1},W_{i},\Phi _{i}$.

die wir miteinander "verheiraten" müssen. Das machen wir beim Komponieren des Gesamt-Gleichungssystems.

Die virtuelle Formänderungsenergie eines Finiten Elements ist nun

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

mit

$\delta \Pi =\displaystyle \sum _{1}^{3}\delta \Pi _{i}$.

Die virtuelle Arbeit der Kraft ist dabei

$\delta W^{a}=F\cdot \delta U_{3}$

Aus den Bildern für die Knoten-Koordinaten oben lesen wir nun die (geometrischen) Rand- und Übergangsbedingungen ab und skizzieren darunter die Einarbeitung dieser Bedingungen in die virtuellen Formänderungsenergien.

Element I	Element II	Element III
${\begin{array}{lcl}W_{0}&=&0\\\end{array}}$	${\begin{array}{lcl}W_{1}&=&0\\W_{2}&=&0\end{array}}$	${\begin{array}{lcl}\%\end{array}}$

Dabei sind folgende Koordinaten sind in den drei Elementen gleich oder unterscheiden sich nur durch ein Vorzeichen:

${\begin{array}{ccccc}I&&II&&III\\\Phi _{1}&=&\Phi _{1}&&\\&&\Phi _{2}&=&\Phi _{2}\\W_{1}&&=&&-W_{2}&{\text{(Sektion II wird nicht länger!)}}\end{array}}$

Was noch fehlt ist die Übersetzung der lokalen in die globalen Koordinaten

local		global
${\begin{array}{c}\Phi _{0}\\W_{1}\\\Phi _{1}\\\Phi _{2}\\W_{3}\\\Phi _{3}\end{array}}$	↔	${\begin{array}{c}\Psi _{0}\\U_{1}\\\Psi _{1}\\\Psi _{2}\\U_{3}\\\Psi _{3}\end{array}}$

1+1

tmp

Assembly of System Matrices

Wir bauen nun das Gleichungssystem

zusammen, indem wir die Anteile der Element-Steifigkeitsmatrizen und die äußeren Lasten in die jeweiligen System Matrizen hineinaddieren.

Globale Koordinaten

Dafür müssen wir die Koordinaten zwischen lokalen und globalen Koordinaten transformieren, nämlich

${\begin{array}{ccr}\Phi _{0}&=&\Psi _{0}\\W_{1}&=&U_{1}\\\Phi _{1}&=&\Psi _{1}\\\Phi _{2}&=&-\Psi _{2}\\W_{3}&=&-U_{3}\\\Phi _{3}&=&-\Psi _{3}&\end{array}}$.

Die virtuelle Formänderungsenergie ist dann

• für Sektion I:
  $\delta \Pi _{1}=\left(\delta \Psi _{0},\delta U_{1},\delta \Psi _{1}\right)\cdot \displaystyle {\frac {EI}{a^{3}}}{\begin{pmatrix}8\cdot {{a}^{2}}&-24\cdot a&4\cdot {{a}^{2}}\\-24\cdot a&96&-24\cdot a\\4\cdot {{a}^{2}}&-24\cdot a&8\cdot {{a}^{2}}\end{pmatrix}}\cdot \left({\begin{array}{l}\Psi _{0}\\U_{1}\\\Psi _{1}\end{array}}\right)$
• für Sektion II:
  $\delta \Pi _{2}=\left(\delta \Psi _{1},\delta \Psi _{2}\right)\cdot \displaystyle {\frac {EI}{a^{3}}}{\begin{pmatrix}4\cdot {{a}^{2}}&-2\cdot {{a}^{2}}\\-2\cdot {{a}^{2}}&4\cdot {{a}^{2}}\end{pmatrix}}\cdot \left({\begin{array}{l}\Psi _{1}\\\Psi _{2}\end{array}}\right)$
• für Sektion III: $\delta \Pi _{3}=\left(\delta U_{1},\delta \Psi _{2},\delta U_{3},\delta \Psi _{3}\right)\cdot \displaystyle {\frac {EI}{a^{3}}}{\begin{pmatrix}96&24\cdot a&-96&24\cdot a\\24\cdot a&8\cdot {{a}^{2}}&-24\cdot a&4\cdot {{a}^{2}}\\-96&-24\cdot a&96&-24\cdot a\\24\cdot a&4\cdot {{a}^{2}}&-24\cdot a&8\cdot {{a}^{2}}\end{pmatrix}}\cdot \left({\begin{array}{l}U_{1}\\\Psi _{2}\\U_{3}\\\Psi _{3}\end{array}}\right)$

Die Gesamt-Steifigkeitsmatrix aus Sektion 1 (rot) , Sektion 2 (grün) und Sektion 3 (blau) ist damit

$\displaystyle {\underline {\underline {K}}}={\frac {E\,I}{a^{3}}}\cdot {\begin{pmatrix}{\color {red}{8a^{2}}}&{\color {red}{-24a}}&{\color {red}{4a^{2}}}&0&0&0\\{\color {red}{-24a}}&{\color {red}{96}}+{\color {blue}{96}}&{\color {red}{-24a}}&{\color {blue}{24a}}&{\color {blue}{-96}}&{\color {blue}{24a}}\\{\color {red}{4a^{2}}}&{\color {red}{-24a}}&{\color {red}{8a^{2}}}+{\color {green}{4a^{2}}}&{\color {green}{-2a^{2}}}&0&0\\0&{\color {blue}{24a}}&{\color {green}{-2a^{2}}}&{\color {green}{4a^{2}}}+{\color {blue}{8a^{2}}}&{\color {blue}{-24a}}&{\color {blue}{4a^{2}}}\\0&{\color {blue}{-96}}&0&{\color {blue}{-24a}}&{\color {blue}{96}}&{\color {blue}{-24a}}\\0&{\color {blue}{24a}}&0&{\color {blue}{4a^{2}}}&{\color {blue}{-24a}}&{\color {blue}{8a^{2}}}\end{pmatrix}}$.

Die Spaltenmatrix der eingeprägten, äußeren Lasten auf das System ist dementsprechend

$\displaystyle {\underline {P}}=\left({\begin{array}{r}0\\0\\0\\0\\F\\0\end{array}}\right)$.

1+1

Solving

Das lineare Gleichungssystem lautet also

$\displaystyle {\frac {EI}{a^{3}}}\cdot {\begin{pmatrix}8\cdot {{a}^{2}}&-24\cdot a&4\cdot {{a}^{2}}&0&0&0\\-24\cdot a&192&-24\cdot a&24\cdot a&-96&24\cdot a\\4\cdot {{a}^{2}}&-24\cdot a&12\cdot {{a}^{2}}&-2\cdot {{a}^{2}}&0&0\\0&24\cdot a&-2\cdot {{a}^{2}}&12\cdot {{a}^{2}}&-24\cdot a&4\cdot {{a}^{2}}\\0&-96&0&-24\cdot a&96&-24\cdot a\\0&24\cdot a&0&4\cdot {{a}^{2}}&-24\cdot a&8\cdot {{a}^{2}}\end{pmatrix}}\cdot {\begin{pmatrix}{{\Psi }_{0}}\\{{U}_{1}}\\{{\Psi }_{1}}\\{{\Psi }_{2}}\\{{U}_{3}}\\{{\Psi }_{3}}\end{pmatrix}}={\begin{pmatrix}0\\0\\0\\0\\F\\0\end{pmatrix}}$

Und dessen Lösung ist

$\left({\begin{array}{c}\Psi _{0}\\U_{1}\\\Psi _{1}\\\Psi _{2}\\U_{3}\\\Psi _{3}\end{array}}\right)=\left({\begin{array}{c}\displaystyle {\frac {3}{8a}}\\\displaystyle {\frac {1}{6}}\\\displaystyle {\frac {1}{4a}}\\\displaystyle {\frac {1}{4a}}\\\displaystyle {\frac {1}{3}}\\\displaystyle {\frac {1}{8a}}\end{array}}\right)\cdot \displaystyle {\frac {a^{3}F}{EI}}$.

/* solve .....                                        */
sol: solve(eom,Y)[1];

ACM : augcoefmatrix(eom, Y);
ole : [K[0] = submatrix(ACM,7), p = -col(ACM,7)];
print(subst(ole, K[0]),"*",transpose(Y),"=",subst(ole,p))$

Post-Processing

Verformung des Elements.

Wir normieren die Auslenkung des Systems mit

$\displaystyle {\hat {w}}=5{\frac {Fa^{3}}{3EI}}$

und plotten das System damit im ausgelenkten Zustand:

Die Auslenkung in Knoten 3 ist gleichzeitig die Längung der insgesamt 5 Sub-Systeme, die Auslenkung des Gesamtsystems in B ist demnach

${\begin{array}{ll}u_{B}&=5\cdot U_{3}\\&=\displaystyle {\frac {5a^{3}}{3EI}}F\end{array}}$

und die gesuchte Ersatzfedersteifigkeit des Systems ist

$\displaystyle k={\frac {3EI}{5a^{3}}}$.

/* post-processing                                     */

displ: makelist(expand(subst(norm,(subst(sol,subst(globals,subst(kinematics[j],subst(X,Q[j])))).transpose(subst([ℓ[i]=l[j]],phi))/ŵ))),j,1,3);
offset: subst(norm,subst(sol,U[1]/ŵ));

plot2d([[parametric, (1-t)/2,  +subst([xi=t],displ[1]), [t,0,1]],
        [parametric, +subst([xi=t],displ[2]), offset+t, [t,0,1]],
        [parametric, (  t)/2, 1-subst([xi=t],displ[3]), [t,0,1]]],
        [x,-1/2,2/2],[y,-0.1,1.4], [same_xy , true],
        [legend,"I","II","III"])$

Links

• ...

Literature

• ...

g