Gelöste Aufgaben/SKEB

Aufgabenstellung

Die Bewegung des Balkens wird durch das Zusammenspiel von elastischen Verformungen und Trägheitskräften bestimmt. Man nennt das "Schwingungen von Kontinua" - diese untersuchen wir hier. Der zentrale Aufgabe besteht in der Berechnung der homogenen Lösung - und der Anpassung der Lösungsanteile an die Anfangsbedingungen.

Lageplan

Gesucht ist analytische Lösung für Schwingungen des Euler-Bernoulli-Balkens beim Loslassen aus der enspannten Ruhelage.

Lösung mit Maxima

Header

Für die mathematische Behandlung - insbesondere der Auflösung quadratischer Gleichungen - setzen wir in Maxima voraus, dass

$E>0,\;\;\;I>0,\;\;\;\ell _{0}>0,\;\;\;\varrho >0;$.

/*******************************************************/
/* MAXIMA script                                       */
/* version: wxMaxima 16.04.2                           */
/* author: Andreas Baumgart                            */
/* last updated: 2018-01-15                            */
/* ref: oscillations of continuous media               */
/* description: analytic solutioun for EBB             */
/*******************************************************/

/* settings */
workdir: "C:\\Users\\X...X\\plots\\";

/* assumptions */
assume(rho>0, E>0,I>0,A>0,l[0]>0,alpha>0,gamma>0);

/* abbreviations used later .... */
abbrev : [alpha = (sqrt(omega[0])*A^(1/4)*rho^(1/4)*l[0])/(E^(1/4)*I^(1/4)),
          gamma = ((cosh(alpha)+cos(alpha)))/(sinh(alpha)+sin(alpha))];
abbrev : append(abbrev, solve(abbrev,[sqrt(omega[0]),cosh(alpha)])[1]);
abbrev : append(abbrev, [abbrev[3]^2,abbrev[3]^3]);

Equations of Motion

In der Gleichgewichtsbeziehung für den Euler-Bernoulli-Balken setzen wir als eingeprägte, äußere Streckenlast q(x,t) die D'Alembert'sche Trägheitskraft und die Gewichtskraft an, also

$q(x,t)=\varrho \;A\cdot g-\varrho \;A\cdot {\ddot {w}}(x,t)$

Die Lösung der linearen, partielle Bewegungsgleichung

$\varrho \,A\cdot {\ddot {w}}+E\,I\cdot {w}^{IV}=\varrho \,A\cdot g$

setzt sich dann aus zwei Lösungsanteilen, der partikularen Lösung w_p und der homogenen Lösung w_h, zusammen; wir schreiben

$w_{t}(x,t)=w_{p}(x,t)+w_{h}(x,t)$.

/* equation of motion               */
eom: rho*A*diff(w(x,t),t,2) + E*I*diff(w(x,t),x,4) = rho*A*g;

Particular Solution

Die partikulare Lösung w_p erfüllt die "rechte Seite" der Bewegungsgleichung, also ϱ A⋅g:

$E\,I\cdot {w}_{p}^{IV}=\varrho \,A\cdot g$.

Die rechte Seite ist zeit-unveränderlich - so auch die partikulare Lösung.

Wir integrieren die Bewegungsgleichung vier Mal und erhalten

$\displaystyle E\,I\,w\left(x,t\right)={\frac {A\,g\,\rho {{x}^{4}}}{24}}+{{C}_{3}}\,{\frac {{x}^{3}}{6}}+{{C}_{2}}\,{\frac {{x}^{2}}{2}}+{{C}_{1}}\,x+{{C}_{0}}$.

Die vier Integrationskonstanten C_i müssen wir nun an die Randbedingungen

${\begin{array}{rl}w_{p}(0)&=0\\w_{p}'(0)&=0\\E\,I\,w_{p}''(\ell )&=0\\E\,I\,w_{p}'''(\ell )&=0\end{array}}$

anpassen, wir erhalten mit dem linearen Gleichungssystem

${\begin{array}{l}0={{C}_{0}}\\0={{C}_{1}}\\0={\frac {{{\ell }_{0}^{2}}Ag\rho }{2}}+{{\ell }_{0}}\,{{C}_{3}}+{{C}_{2}}\\0={{\ell }_{0}}Ag\rho +{{C}_{3}}\end{array}}$

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)$.

Statische Auslenkung

Die maximale Auslenkung - am rechten Rand - nutzen wir als Bezugslänge

$\displaystyle {{w}_{s}}={\frac {{{\ell }^{4}}\,A\,g\,\rho }{8\,E\,I}}$

Und so sieht w_p aus:

/***************************************************************************/
/* p a r t i c u l a r   s o l u t i o n                                   */
/***************************************************************************/
peom : subst(['diff(w(x,t),t,2)=0],eom);
w[p] : integrate(integrate(integrate(integrate(peom,x),x),x),x);
w[p] : subst([%c1=C[3],%c2=C[2],%c3=C[1],%c4=C[0]],w[p]);

BouCon : [subst([x= 0  ],subst([      w(x,t)     =0],     w[p]     )),
          subst([x= 0  ],subst(['diff(w(x,t),x,1)=0],diff(w[p],x,1))),
          subst([x=l[0]],subst(['diff(w(x,t),x,2)=0],diff(w[p],x,2))),
          subst([x=l[0]],subst(['diff(w(x,t),x,3)=0],diff(w[p],x,3)))];
w[p] : subst(solve(BouCon,makelist(C[i],i,0,3))[1],w[p]);
w[p] : ratsimp(lhs(w[p]/(E*I)) = subst([x=l[0]*xi], rhs(w[p])/(E*I)));

params: [w[s] = subst([xi=1],subst(w[p],w(x,t))), omega[0]=2*%pi/T];

plot2d(ratsimp(subst(params,subst(w[p],w(x,t))/w[s])),[xi,0,1], 
               [legend, "w[p]/w[s]"],
               [style, [lines,3,1]], [xlabel, "ξ →"], [ylabel, "← w"],
               [gnuplot_preamble, "set yrange [1.1:-0.1] reverse"])$

Homogeneous Solution

Die homogene Bewegungsgleichung

$\varrho \,A\cdot {\ddot {w}}+E\,I\cdot {w}^{IV}=0$

lässt Lösungen vom Typ

$\displaystyle {w}_{h}\left(x,t\right)={\hat {W}}\cdot {e}^{\kappa \,x}\cdot {e}^{j\,\omega _{0}\,t}\;\;\;{\text{ mit }}\;\;\;j={\sqrt {-1}}$

zu. Wir bekommen zu jedem ω₀ 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}}$,

also

$w_{h}(x,t)=\left({{\hat {W}}_{1}}\,{{e}^{\displaystyle {\frac {j\,\alpha \,x}{\ell }}}}+{{\hat {W}}_{3}}\,{{e}^{\displaystyle -{\frac {j\,\alpha \,x}{\ell }}}}+{{\hat {W}}_{4}}\,{{e}^{\displaystyle {\frac {\alpha x}{\ell }}}}+{{\hat {W}}_{2}}\,{{e}^{\displaystyle -{\frac {\alpha x}{\ell }}}}\right)\cdot e^{\displaystyle j\;\omega _{0}\,t}$

Die komplexen Exponenten deuten schon darauf hin: im Allgemeinen sind die Koeffizienten W_i ebenfalls komplex. Für dieses Beispiel macht es deshalb viel Sinn, diese allgemeine Lösung der homogenen Differentialgleichung etwas anders hinzuschreiben. Wir nutzen dafür:

• die Hyperbelfunktionen
  

  ${\begin{array}{l}\displaystyle \sinh(z)={\frac {e^{z}-e^{-z}}{2}}\\\displaystyle \cosh(z)={\frac {e^{z}+e^{-z}}{2}}\end{array}}$

• die Euler'sche Formel
  

  $\displaystyle e^{j\,z}=\cos(z)+j\cdot \sin(z)$

Dann lautet die allgemeine Lösung  

$\displaystyle {w}\left(x,t\right)=\left({{W}_{1}}\,{\sinh }\left({\frac {\alpha \,x}{\ell }}\right)+{{W}_{2}}\,\cosh {\left({\frac {\alpha \,x}{\ell }}\right)}+{{W}_{3}}\,\sin {\left({\frac {\alpha \,x}{\ell }}\right)}+{{W}_{4}}\,\cos {\left({\frac {\alpha \,x}{\ell }}\right)}\right)\cdot {{e}^{\displaystyle j\,{{\omega }_{0}}\,t}}$,

bei der wir davon ausgehen dürfen, dass alle W_i reellwertig sind.

Wie bei der partikularen Lösung w_p sind es die vier Integrationskonstanten W_i, die wir nun an die Randbedingungen

${\begin{array}{rl}w_{h}(0,t)&=0\\w_{h}'(0,t)&=0\\E\,I\,w_{h}''(\ell ,t)&=0\\E\,I\,w_{h}'''(\ell ,t)&=0\end{array}}$

anpassen müssen. Wir erhalten das lineare Gleichungssystem

${\begin{array}{ccc}\displaystyle {{W}_{2}}+{{W}_{4}}&=&0,\\\displaystyle {\frac {{{W}_{1}}\alpha }{\ell }}+{\frac {{{W}_{3}}\alpha }{\ell }}&=&0,\\\displaystyle {\frac {{{W}_{1}}\,{{\alpha }^{2}}\,{\sinh }\left(\alpha \right)}{\ell ^{2}}}+{\frac {{{W}_{2}}\,{{\alpha }^{2}}\,\cosh {\left(\alpha \right)}}{\ell ^{2}}}-{\frac {{{W}_{3}}\,{{\alpha }^{2}}\,\sin {\left(\alpha \right)}}{\ell ^{2}}}-{\frac {{{W}_{4}}\,{{\alpha }^{2}}\,\cos {\left(\alpha \right)}}{\ell ^{2}}}&=&0,\\\displaystyle {\frac {{{W}_{1}}\,{{\alpha }^{3}}\,\cosh {\left(\alpha \right)}}{\ell ^{3}}}+{\frac {{{W}_{2}}\,{{\alpha }^{3}}\,{\sinh }\left(\alpha \right)}{\ell ^{3}}}-{\frac {{{W}_{3}}\,{{\alpha }^{3}}\,\cos {\left(\alpha \right)}}{\ell ^{3}}}+{\frac {{{W}_{4}}\,{{\alpha }^{3}}\,\sin {\left(\alpha \right)}}{\ell ^{3}}}&=&0\end{array}}$,

vereinfacht und in Matrix-Schreibweise

$\underbrace {\begin{pmatrix}0&1&0&1\\1&0&1&0\\\operatorname {sinh} \left(\alpha \right)&\cosh {\left(\alpha \right)}&-\sin {\left(\alpha \right)}&-\cos {\left(\alpha \right)}\\\cosh {\left(\alpha \right)}&\operatorname {sinh} \left(\alpha \right)&-\cos {\left(\alpha \right)}&+\sin {\left(\alpha \right)}\end{pmatrix}} _{\displaystyle {\underline {\underline {D}}}(\kappa )}\cdot {\begin{pmatrix}W_{1}\\W_{2}\\W_{3}\\W_{4}\end{pmatrix}}={\underline {0}}$.

Für dieses lineare, homogenen Gleichungen gibt es nur dann nicht-triviale Lösungen, wenn die Matrix einen Rangabfall >= 1 hat, also

$\det({\underline {\underline {D}}}(\alpha )){\stackrel {!}{=}}0$.

Wir finden die charakteristische Gleichung

$-{{{\sinh }\left(\alpha \right)}^{2}}+{{\sin {\left(\alpha \right)}}^{2}}+{{\cosh {\left(\alpha \right)}}^{2}}+2\cos {\left(\alpha \right)}\,\cosh {\left(\alpha \right)}+{{\cos {\left(\alpha \right)}}^{2}}{\stackrel {!}{=}}0$.

Beim Plotten der Funktion sehen wir bereits, dass es mehr als eine Nullstelle gibt, hier bei α≈+/- 2 und α≈+/- 5.

Charakteristische Gleichung

Die charakteristische Gleichung können wir umformen zu

$\cos(\alpha )\cdot \cosh(\alpha )+1{\stackrel {!}{=}}0$

und beide Anteile auftragen. So wird klar, dass die beiden Funktionen sich wegen der Periodizität der cos-Funktion unendlich oft schneiden - es also unendlich viele Nullstellen gibt:

Nullstellensuche (Ausschnitt)

Nullstellensuche

Analytisch finden wir keine Lösungen dieser charakteristischen Gleichung - wir müssen sie numerisch bestimmen – hier mit der Maxima-Funktion "find_root". Die ersten (positiven) Lösungen sind

${\begin{array}{l}\alpha _{1}=1.875104068711961,\\\alpha _{2}=4.694091132974175,\\\alpha _{3}=7.854757438237613,\\\alpha _{4}=10.99554073487547,\\\alpha _{5}=14.13716839104647,\\\alpha _{6}=17.27875953208824\end{array}}$.

Sie zu finden ist einfach, weil die α_i für i>1 ungefähr π auseinander liegen - da können wir keine Nullstelle verpassen.

In diesen α_i stecken über $\omega _{0,i}=2\,\pi /T_{i}$ auch der Periodendauern, z.B. für die erste Mode

$\displaystyle T_{1}=0.5688\ldots \pi \ell ^{2}{\sqrt {\frac {\varrho \,A}{E\,I}}}$.

Wie bei einem normalen Eigenwertproblem berechnen wir nun die Eigenvektoren W aus der Beziehung

${\underline {\underline {D}}}(\alpha _{i})\cdot {\underline {W}}={\underline {0}}$.

Aus unendlich vielen Lösungen α_i folgen nun also auch unendlich viele W_j, die wir zu

$\left({\begin{array}{c}W_{1}\\W_{2}\\W_{3}\\W_{4}\end{array}}\right)_{j}=C_{j}\cdot \left({\begin{array}{r}\gamma _{j}\\-1\\-\gamma _{j}\\1\end{array}}\right){\text{ mit }}\displaystyle \gamma _{j}={\frac {\cosh {\left(\alpha _{j}\right)}+\cos {\left(\alpha _{j}\right)}}{\sinh \left(\alpha _{j}\right)+\sin {\left(\alpha _{j}\right)}}}$

bestimmen. Die homogene Lösung der Bewegungsgleichung ist also

$\displaystyle w_{h}(x,t)=\sum _{j=1}^{\infty }C_{j}\cdot \underbrace {\left(\gamma _{i},-1,-\gamma _{j},1\right)\cdot \left({\begin{array}{r}\sinh(\alpha _{j}\,\xi )\\\cosh(\alpha _{j}\,\xi )\\\sin(\alpha _{j}\,\xi )\\\cos(\alpha _{j}\,\xi )\end{array}}\right)} _{\displaystyle =:\phi _{j}(\xi )}\cdot e^{\displaystyle \omega _{0,j}\,t}{\text{ mit }}\xi ={\frac {x}{\ell }}$.

Denken Sie daran: die Eigenkreisfrequenz ω_0j ist dabei wieder eine Funktion von α_i - siehe oben. Die C_i sind wieder Konstanten, die wir brauchen, um die Lösung an die Anfangsbedingungen anzupassen. Zunächst schauen wir uns die Lösung jeweils zu einem α_i an. Diese Funktionen heißen Modalformen ϕ und deren Schwingungen können wir plotten:

Mode	Modalform ϕj	Mode	Modalform ϕj
#1 $\displaystyle \omega _{0,1}=1\cdot {\frac {2\pi }{T}}$		#2 $\displaystyle \omega _{0,2}=6.27\cdot {\frac {2\pi }{T}}$
#3 $\displaystyle \omega _{0,3}=17.5\cdot {\frac {2\pi }{T}}$		#4 $\displaystyle \omega _{0,4}=34.4\cdot {\frac {2\pi }{T}}$
#5 $\displaystyle \omega _{0,5}=56.8\cdot {\frac {2\pi }{T}}$		#6 $\displaystyle \omega _{0,6}=84.9\cdot {\frac {2\pi }{T}}$

Modalformen

Diese Moden ϕ_j können wir auch einzeln übereinander darstellen:

/***************************************************************************/
/* h o m o g e n e o u s   s o l u t i o n                                 */
/***************************************************************************/
ansatz: w(x,t) = W[i] * %e^(kappa*x) * %e^(%i*omega[0]*t);
heom : subst(ansatz,subst([g=0],eom));
heom : ratsimp(subst(ansatz,ev(heom,nouns)/w(x,t)));
sol[1] : subst(abbrev,solve(heom,kappa));
/* ansatz: w(x,t) = sum(subst(sol[1][i],subst(ansatz,w(x,t))),i,1,length(sol[1])); */
/* this is possible - but results in complex W[i] - too demanding */
ansatz : w(x,t) = subst(sol[1][1], W[1]*%e^(kappa*x)*%e^(%i*omega[0]*t))
                 +subst(sol[1][2], W[2]*%e^(kappa*x)*%e^(%i*omega[0]*t))
                 +subst(sol[1][3], W[3]*%e^(kappa*x)*%e^(%i*omega[0]*t))
                 +subst(sol[1][4], W[4]*%e^(kappa*x)*%e^(%i*omega[0]*t));
/* .. but we can rephrase w(x,t) with
             https://de.wikipedia.org/wiki/Sinus_hyperbolicus_und_Kosinus_hyperbolicus

      to have only real-values W[i]:            */
ansatz : w(x,t) =( W[1]*sinh(kappa*x)
                  +W[2]*cosh(kappa*x)
                  +W[3]*sin (kappa*x)
                  +W[4]*cos (kappa*x))*%e^(%i*omega[0]*t);

/* test:subst(ansatz,heom);ratsimp(ev(%,nouns));subst([x=xi*l[0]],%);subst(abbrev[1],%);*/

/* boundary conditions */
BouCon : [subst([x=  0 ],      subst(ansatz,w(x,t))     )/%e^(%i*omega[0]*t)=0,
          subst([x=  0 ], diff(subst(ansatz,w(x,t)),x  ))/%e^(%i*omega[0]*t)=0,
          subst([x=l[0]], diff(subst(ansatz,w(x,t)),x,2))/%e^(%i*omega[0]*t)=0,
          subst([x=l[0]], diff(subst(ansatz,w(x,t)),x,3))/%e^(%i*omega[0]*t)=0];
D : submatrix(augcoefmatrix(expand(BouCon),makelist(W[i],i,1,4)),5);
/* what happens to the six alpha=0 ???? */
D : diag_matrix (1, 1/kappa, 1/kappa^2,1/kappa^3).subst(abbrev,D);
D : subst(solve(sol[1][4],kappa),D);
CharEqua : expand(determinant(D))=0;
CharEqua : ratsimp(subst([sin(alpha)^2=1-cos(alpha)^2,sinh(alpha)^2=cosh(alpha)^2-1], CharEqua/2));

/* numerically find zeros for omega[0] (alpha) */
plot2d(lhs(CharEqua),[alpha,-6,6], [y,-6,6], [ylabel,"lhs(CharEqua)→"], [xlabel,"α→"], [style, [lines,3]])$
plot2d([cos(alpha),-1/cosh(alpha)],[alpha,0,60], [ylabel,"f(α)→"], [xlabel,"α→"], [style, [lines,3]])$
guessed: [2,5,8,11,14,17];
zeros : makelist(alpha = find_root(CharEqua/alpha^6, alpha, guessed[i]-1/2, guessed[i]+1/2),i,1,length(guessed));
sol[2] : makelist(omega[0] = subst(zeros[i],subst(abbrev[3],sqrt(omega[0]))^2),i,1,length(zeros));
params : append(params, [T[ref] = subst(zeros[1],subst(abbrev,2*%pi/omega[0]))]);

/* use upper triangular matrix instead of D */
[P,L,U] : get_lu_factors(lu_factor(D))$

/* get eigenvectors */
EigVec : subst(solve(makelist(sum(U[j][i]*W[i],i,1,4)=0,j,1,3),[W[1],W[2],W[3]])[1],
                                   matrix([W[1]],[W[2]],[W[3]],[W[4]]));
EigVec : subst(abbrev,subst([W[4]=1],EigVec));
Gammas : makelist(gamma = subst(zeros[i],subst(abbrev[2],gamma)),i,1,length(zeros));

modeshape : subst(abbrev[2],subst([x=xi*l[0]],subst(sol[1][4],subst(makelist(W[i]=EigVec[i][1],i,1,4),
                       subst(params,subst([t=T[ref]*tau],subst(abbrev,subst(ansatz,w(x,t)))))))));
modes : makelist(subst(zeros[i],modeshape),i,1,length(zeros));

/* plot mode snapshots for animation -> convert with imageMagick */
fpprintprec: 3$
titles : makelist(simplode(["ω/ω1 = ",subst(sol[2][1],subst(sol[2][i],omega[0])/omega[0])]),i,1,6);
NSteps : 19$
for mode: 1 thru 6 do
   for j:0 thru NSteps do
      (tstep : simplode([if j<10 then "0" else "",j]),
       plot2d(subst([tau=j/(NSteps+1)],subst([tau=0],modes[mode])*sin(2*%pi*tau)),[xi,0,1],
                        [title, titles[mode]], [legend, simplode(["t/T=",float(j/NSteps+0.01)])],
                        [style, [lines,3,1]], [xlabel, "ξ →"], [ylabel, "← w "],
                        [png_file,simplode([workdir,"mode-",mode,"-step-",tstep,".png"])],
                        [gnuplot_preamble, "set yrange [2.1:-2.1] reverse"]))$

/* plot all modeshapes (1...6)*/
legends : append([legend],makelist(simplode(["mode ",i]), i,1,length(zeros)));
plot2d(subst([tau=0],modes), [xi,0,1], legends,
             [xlabel, "ξ →"], [ylabel, "← w "],
             [gnuplot_preamble, "set yrange [2.1:-2.1] reverse"])$

Adapt to Initial Condition

Zum Zeitpunkt t=0 soll

$\underbrace {w_{h}(x,0)} _{\displaystyle =\sum _{j=1}^{\infty }C_{j}\cdot \phi _{j}(x)}+w_{p}(x)=0$

gelten. Das kriegen wir hin, indem wir die (unendlich vielen) C_i der homogenen Lösung passend wählend. Das klingt komplizierter als es ist - die Mathematik bereitet uns den Weg.

Die Beziehung oben multiplizieren wir dazu mit einer der Modalformen ϕ_k und bilden das Integral über die Balken-Lange ℓ, also

$\displaystyle \sum _{j=1}^{\infty }\left(C_{j}\cdot \int _{0}^{\ell }\phi _{j}(x)\cdot \phi _{k}(x)\;dx\right)+\int _{0}^{\ell }w_{p}(x)\cdot \phi _{k}(x)\;dx=0$

Wir projizieren also die Modalformen auf sich selbst und auf die partikulare Lösung. Und nutzen aus, dass diese Faltungs-Integrale der Modalformen eine besondere Eigenschaft haben:

$\displaystyle \int _{0}^{\ell }\phi _{j}\cdot \phi _{k}\;dx\left\{{\begin{array}{l}=0{\text{ für }}j\neq k\\\neq 0{\text{ für }}j=k\end{array}}\right.$

Man sagt: die Funktionen sind zueinander verallgemeinert orthogonal. Damit berechnen wir

${\begin{array}{l}C_{1}=0.506693895146775,\\C_{2}=0.007150056863941781,\\C_{3}=5.347110562156877\;10^{-4},\\C_{4}=9.955236981148623\;10^{-5},\\C_{5}=2.833383071646469\;10^{-5},\\C_{6}=1.038598312000442\;10^{-5}\\\;\;\;\;\;\vdots \end{array}}$,

Nur die erste Modalform liefert einen signifikanten Beitrag bei der Anpassung an die Anfangsbedingungen, wir nehmen bei der Auftragung der Gesamtlösung w_t(x,t) trotzdem die ersten sechs Modalformen mit:

Animation: Loslassen des Systems aus der unverformten Rugelage.

In der Animation sehen Sie das Loslassen des Stabes aus der unverformten, geraden Referenzlage für genau eine Periodenlänge T₁ der ersten Modalform. Allerdings ist die Schwingung, die entsteht, nicht T₁-periodisch - die anderen Modalformen leisten einen (kleinen) Beitrag, dessen Schwingungsperiode kein Vielfaches der Grundschwingung ist. Das sehen sie an dem kleinen "Ruck" in der Darstellung der Lösung, wenn nach der Dauer T₁ der Grundschwingung wieder die Anfangslage eingeblendet wird.

/***************************************************************************/
/* t o t a l   ( c o m p o s e d )   s o l u t i o n                       */
/***************************************************************************/

w[h] : sum(C[i]*modes[i],i,1,length(zeros));

/* adapt to inital condition w[h]+w[p]=0  */
/* Faltungsintegrale */
FI : [integrate(expand(subst([tau=0],modeshape^2)),xi,0,1),
      integrate(expand(subst([tau=0],modeshape*ratsimp(subst(params,subst(w[p],w(x,t))/w[s])))),xi,0,1)];
/* coefficients for mode shapes */
IVs : makelist(C[i] = -subst(zeros[i],FI[2]/FI[1]),i,1,length(zeros));

w[h] : subst(IVs,realpart(w[h]));

/* produce slides for animation */
fpprintprec: 3$
NSteps : 39$
for j:0 thru NSteps-1 do
   (tstep : simplode([if j<10 then "0" else "",j]),
    plot2d(subst([tau=j/NSteps],w[h])+ratsimp(subst(params,subst(w[p],w(x,t))/w[s])),[xi,0,1], 
               [legend, simplode(["t/T=",float(j/NSteps+0.01)])],
               [style, [lines,3,1]], [xlabel, "ξ →"], [ylabel, "← w"],
               [png_file,simplode([workdir,"initial-value-solution-step-",tstep,".png"])],
               [gnuplot_preamble, "set yrange [2.5:-0.5] reverse"]))$

🧨 Timoshenko-Balken:

Wieviele Moden man berücksichtigen muss, hängt von vielen praktischen Faktoren - z.B. der modalen Dämpfung - ab. Vergessen sollte man nicht, dass je nach Stab-Länge, der Einfluss der Schubverformung bei den höheren Moden zunimmt. Dann muss man ggf. das Modell wechseln ...

✔ Saiteninstruente:

Bei Saiteninstrumenten möchte man Biegeschwingungen der Saiten verhindern. Darum sind besonders die Saiten, die eine hohe Massebelegung brauchen, um in niedrigen Frequenzen zu schwingen, besondere Konstruktionen. Dabei wird oft ein Draht um die eigentliche Saite gewickelt. Können Sie aus diesen Untersuchungen erkennen, was das Ziel dieser Konstruktion ist?

---

Links

• Schwingungen von Kontinua eines Dehnstabes (SKER)

Literature

• ...