Gelöste Aufgaben/Kw24: Unterschied zwischen den Versionen

Aus numpedia
Zur Navigation springen Zur Suche springen
Keine Bearbeitungszusammenfassung
Zeile 19: Zeile 19:


== Lösung mit Maxima ==
== Lösung mit Maxima ==
Lorem Ipsum ....


==tmp==
<!-------------------------------------------------------------------------------->
Gesucht ist die Lösung zu dem nichtlinearen Anfangswertproblem.<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Header
{{MyCodeBlock|title=Header
|text=Text
|text=
Gesucht ist die Lösung zu dem nichtlinearen Anfangswertproblem.
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/* MAXIMA script                                        */
/* version: wxMaxima 15.08.2                            */
/* author: Andreas Baumgart                              */
/* last updated: 2017-09-30                              */
/* ref: Kw24 (TM-C, Labor 6)                            */
/* description: finds the solution for                  */
/*              the nonlinear IVP                        */
/*********************************************************/
</syntaxhighlight>
</syntaxhighlight>
}}
}}


==tmp==
 
<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Declarations
|text=
Wir wählen die Längenverhältnisse
Wir wählen die Längenverhältnisse


<math>\displaystyle {{\ell}_{2}}={{\ell}_{1}}, {{\ell}_{3}}=\sqrt{2}\cdot {{\ell}_{1}}</math>
::<math>\displaystyle {{\ell}_{2}}={{\ell}_{1}}, {{\ell}_{3}}=\sqrt{2}\cdot {{\ell}_{1}}</math>


so dass die Feder für φ''(t) = 45''° entspannt ist.<!-------------------------------------------------------------------------------->
so dass die Feder für φ''(t) = 45''° entspannt ist.
 
{{MyCodeBlock|title=Declarations
|text=Text
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/* declare variational variables - see 6.3 Identifiers  */
declare("δW", alphabetic); /* total virtual work        */
declare( "φ", alphabetic); /* coordinate of rotatin      */
declare("δφ", alphabetic); /* variation of coordinate phi*/
declare("δw", alphabetic); /* variation of coordinate w  */
 
/* parameter */
assume(l[1]>0, k>0, J>0);
params: [l[2] = l[1], l[3] = sqrt(2)*l[1]];
</syntaxhighlight>
</syntaxhighlight>
}}
}}


==tmp==
<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Kinematics
|text=
Der Mechanismus hat genau einen Freiheitsgrad, den wir entweder durch die Federlängung ''w'' oder den Drehwinkel ''φ'' erfassen. Den Zusammen zwischen den beiden stellen wir über
Der Mechanismus hat genau einen Freiheitsgrad, den wir entweder durch die Federlängung ''w'' oder den Drehwinkel ''φ'' erfassen. Den Zusammen zwischen den beiden stellen wir über


<math>{{\ell}_{1}^{2}}\cdot {{\mathrm{sin}\left( \mathrm{\varphi}\left( t\right) \right) }^{2}}+{{\left( {{\ell}_{3}}-{{\ell}_{1}}\cdot \mathrm{cos}\left( \mathrm{\varphi}\left( t\right) \right) \right) }^{2}}={{\left( {{\ell}_{2}}+w\right) }^{2}}</math>
::<math>{{\ell}_{1}^{2}}\cdot {{\mathrm{sin}\left( \mathrm{\varphi}\left( t\right) \right) }^{2}}+{{\left( {{\ell}_{3}}-{{\ell}_{1}}\cdot \mathrm{cos}\left( \mathrm{\varphi}\left( t\right) \right) \right) }^{2}}={{\left( {{\ell}_{2}}+w\right) }^{2}}</math>


her. Nun steht uns frei, ob wir diese Beziehung nach ''w'' oder ''φ'' auflösen. Wir wählen ''φ'' als Minimal-Koordinate und erhalten
her. Nun steht uns frei, ob wir diese Beziehung nach ''w'' oder ''φ'' auflösen. Wir wählen ''φ'' als Minimal-Koordinate und erhalten


<math>w = \displaystyle \frac{{{\ell}_{1}}\cdot \sqrt{12-{{2}^{\frac{7}{2}}}\cdot \cos(\varphi) }-2\cdot {{\ell}_{1}}}{2}</math>
::<math>w = \displaystyle \frac{{{\ell}_{1}}\cdot \sqrt{12-{{2}^{\frac{7}{2}}}\cdot \cos(\varphi) }-2\cdot {{\ell}_{1}}}{2}</math>


Nach den Regeln der [[Sources/Lexikon/Variationsmethoden|Variationsmethoden]] erhalten wir daraus:
Nach den Regeln der [[Sources/Lexikon/Variationsmethoden|Variationsmethoden]] erhalten wir daraus:


<math>\delta w = \displaystyle \frac{{2^{\frac{3}{2}}}\cdot {{\ell}_{1}}\cdot \sin( \varphi ) \cdot \delta \varphi}{\sqrt{12-{{2}^{\frac{7}{2}}}\cdot \cos(\varphi)}}</math>.<!-------------------------------------------------------------------------------->
::<math>\delta w = \displaystyle \frac{{2^{\frac{3}{2}}}\cdot {{\ell}_{1}}\cdot \sin( \varphi ) \cdot \delta \varphi}{\sqrt{12-{{2}^{\frac{7}{2}}}\cdot \cos(\varphi)}}</math>.
 
{{MyCodeBlock|title=Kinematics
|text=Text
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/**** kinematics ***/
 
/* coordinates */
coords : [φ, δφ];
 
kin : (l[1]*sin(φ(t)))^2+(l[3]-l[1]*cos(φ(t)))^2=(l[2]+w)^2;
kin : trigsimp(ratsimp(solve(subst(params, kin),w)[2]));
kin : [kin,
      δw = subst([epsilon=0],
                diff(subst([φ(t)=φ(t)+epsilon*δφ],
                                  subst(kin,w)),epsilon))];
</syntaxhighlight>
</syntaxhighlight>
}}
}}


==tmp==
 
<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Equilibrium Conditions
|text=
Mit dem Prinzip der virtuellen Verrückungen schreiben wir
Mit dem Prinzip der virtuellen Verrückungen schreiben wir


<math>\displaystyle \delta W=-\left( \frac{{{d}^{2}}}{d\,{{t}^{2}}}\cdot \mathrm{\varphi}\left( t\right) \right) \cdot J\cdot \delta \varphi-k\cdot w\cdot \delta w</math>
::<math>\displaystyle \delta W=-\left( \frac{{{d}^{2}}}{d\,{{t}^{2}}}\cdot \mathrm{\varphi}\left( t\right) \right) \cdot J\cdot \delta \varphi-k\cdot w\cdot \delta w</math>


und erhalten die Bewegungsgleichung
und erhalten die Bewegungsgleichung


<math>\displaystyle J \cdot \ddot{\varphi} + \sqrt{2}\cdot l_1^2 \cdot k \cdot \sin\varphi - \frac{{{2}^{\frac{3}{2}}}\cdot l_1^2 \cdot k \cdot \sin\varphi}{\sqrt{12-{{2}^{\frac{7}{2}}}\cdot \cos\varphi}} = 0</math>
::<math>\displaystyle J \cdot \ddot{\varphi} + \sqrt{2}\cdot l_1^2 \cdot k \cdot \sin\varphi - \frac{{{2}^{\frac{3}{2}}}\cdot l_1^2 \cdot k \cdot \sin\varphi}{\sqrt{12-{{2}^{\frac{7}{2}}}\cdot \cos\varphi}} = 0</math>
<!-------------------------------------------------------------------------------->
 
 
{{MyCodeBlock|title=Equilibrium Conditions
|text=Text
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/* equilibrium conditions */
PvV : δW = -J*'diff(φ(t),t,2)*δφ - k*w*δw;
 
/* dgl for φ ist the coefficient of δφ in δW */
dgl[1] : expand(coeff(subst([δW = 0], subst(kin, -PvV)),δφ));
</syntaxhighlight>
</syntaxhighlight>
}}
}}


==tmp==
 
<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Equilibrium Points
|text=
Die Bewegungsgleichung machen wir dimensionslos, indem wir charakteristische Größen aus der linearisierten Bewegungsgleichung verwenden.
Die Bewegungsgleichung machen wir dimensionslos, indem wir charakteristische Größen aus der linearisierten Bewegungsgleichung verwenden.


Die drei Gleichgewichts-Lagen bzgl. derer wir die Bewegungsgleichung linearisieren können, sind
Die drei Gleichgewichts-Lagen bzgl. derer wir die Bewegungsgleichung linearisieren können, sind


<math>\begin{array}{ll}\varphi_{0,1} &= \displaystyle -\frac{\pi}{4}\\\varphi_{0,2} &= 0\\\varphi_{0,3} &= \displaystyle +\frac{\pi}{4}\end{array}</math>.
::<math>\begin{array}{ll}\varphi_{0,1} &= \displaystyle -\frac{\pi}{4}\\\varphi_{0,2} &= 0\\\varphi_{0,3} &= \displaystyle +\frac{\pi}{4}\end{array}</math>.


Eine Taylor-Reihenentwicklung des Rückstellmoments bzgl. φ''<sub>0,3</sub>'' ,also
Eine Taylor-Reihenentwicklung des Rückstellmoments bzgl. φ''<sub>0,3</sub>'' ,also


<math>\varphi = \varphi_{0,3} + \phi</math>
::<math>\varphi = \varphi_{0,3} + \phi</math>


liefert die linearisierte Bewegungsgleichung  
liefert die linearisierte Bewegungsgleichung  


<math>\displaystyle\ddot{\phi}+\omega_0^2\cdot\phi = 0 \text{ mit } \omega_0^2 = \frac{{{\ell}_{1}^{2}}\cdot k}{J}</math>.
::<math>\displaystyle\ddot{\phi}+\omega_0^2\cdot\phi = 0 \text{ mit } \omega_0^2 = \frac{{{\ell}_{1}^{2}}\cdot k}{J}</math>.


Dabei ist ''ω<sub>0</sub>'' die Eigenkreisfreuquenz der linearisierten Bewegungsgleichung. <!-------------------------------------------------------------------------------->
Dabei ist ''ω<sub>0</sub>'' die Eigenkreisfreuquenz der linearisierten Bewegungsgleichung.  
 
{{MyCodeBlock|title=Equilibrium Points and Stability
|text=Text
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/* Zeros of restoring function */
equilib : solve(subst(['diff(φ(t),t,2)=0],dgl[1]),φ(t));
 
/*********************************************************/
/* derive the eigenfrequency                            */
/* of linearized ODE at stable point of equilob.        */
LinODE : taylor (subst([φ(t)=psi],rhs(dgl[1]/J)), psi, %pi/4,1);
ref : omega[0]^2 = coeff(expand(LinODE),psi);
ref : [t[B] = 1/(subst(solve(ref,omega[0])[2],omega[0])/(2*%pi))];
</syntaxhighlight>
</syntaxhighlight>
}}
}}


==tmp==
 
<!-------------------------------------------------------------------------------->{{MyCodeBlock|title=Make Equations of Motion Dimensionless
|text=
Eine Bezugszeit für unsere nichtlineare Bewegungsgleichung kommt nun aus  
Eine Bezugszeit für unsere nichtlineare Bewegungsgleichung kommt nun aus  


<math>\displaystyle \omega_0 = \frac{2\pi}{t_{Bez}}</math>
::<math>\displaystyle \omega_0 = \frac{2\pi}{t_{Bez}}</math>


Damit lautet die dimensionslose Bewegungsgleichung
Damit lautet die dimensionslose Bewegungsgleichung


<math>\displaystyle \frac{d^2}{d\tau^2}\Phi + \frac{{{2}^{\frac{5}{2}}}\cdot {{\pi }^{2}}\cdot \sqrt{12-{{2}^{\frac{7}{2}}}\cdot \mathrm{cos}\left( \Phi\right) }\cdot \mathrm{sin}\left( \Phi\right) -{{2}^{\frac{7}{2}}}\cdot {{\pi }^{2}}\cdot \mathrm{sin}\left( \Phi\right) }{\sqrt{12-{{2}^{\frac{7}{2}}}\cdot \mathrm{cos}\left( \Phi\right) }} = 0 \text{ mit } \tau = t/t_{Bez}</math>[[Datei:Kw24-11.png|mini|Rückstellmoment der Feder|alternativtext=]]Den nichtlinearen (zweiten) Summand können wir als dimensionsloses "Rückstellmoment auffassen und über dem Drehwinkel auftragen:
::<math>\displaystyle \frac{d^2}{d\tau^2}\Phi + \frac{{{2}^{\frac{5}{2}}}\cdot {{\pi }^{2}}\cdot \sqrt{12-{{2}^{\frac{7}{2}}}\cdot \mathrm{cos}\left( \Phi\right) }\cdot \mathrm{sin}\left( \Phi\right) -{{2}^{\frac{7}{2}}}\cdot {{\pi }^{2}}\cdot \mathrm{sin}\left( \Phi\right) }{\sqrt{12-{{2}^{\frac{7}{2}}}\cdot \mathrm{cos}\left( \Phi\right) }} = 0 \text{ mit } \tau = t/t_{Bez}</math>


Nullstellen des Rückstellmoments zeigen die Gleichgewichtslagen des Systems an. Die Steigung dieses Rückstellmoments in den Gleichgewichtslagen deutet dabei auf die Stabilität der linearisierten Bewegung (oder Bewegung im Kleinen) hin. So ist die Lösung ''Φ'' = 0 eine instabile Gleichgewichtslage, die beiden anderen sind stabil.<!-------------------------------------------------------------------------------->{{MyCodeBlock|title=Make Equations of Motion Dimensionless
[[Datei:Kw24-11.png|mini|Rückstellmoment der Feder|alternativtext=]]Den nichtlinearen (zweiten) Summand können wir als dimensionsloses "Rückstellmoment auffassen und über dem Drehwinkel auftragen:
|text=Text
 
Nullstellen des Rückstellmoments zeigen die Gleichgewichtslagen des Systems an. Die Steigung dieses Rückstellmoments in den Gleichgewichtslagen deutet dabei auf die Stabilität der linearisierten Bewegung (oder Bewegung im Kleinen) hin. So ist die Lösung ''Φ'' = 0 eine instabile Gleichgewichtslage, die beiden anderen sind stabil.
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/* make dgl dimensionless employing the eigenfrequency  */
 
dgl[2] : subst([diff(φ(t),t,2) = 'diff(Phi,t,2)/t[B]^2,
                    φ(t)      =      Phi ] , dgl[1]);
dgl[2] : solve(subst(ref,dgl[2]), 'diff(Phi,t,2));
 
/* equilibrium angles für Φ */
plot2d(-rhs(dgl[2][1]), [Phi,-1,1],
                        [legend, "Rückstell'moment'"],
                        [xlabel, "Φ/1->"], [ylabel, "M/1->"]);
</syntaxhighlight>
</syntaxhighlight>
}}
}}


==tmp==
Für die numerische Lösung müssen wir die Bewegungsgleichung


<math>\ddot{\Phi} = f(\Phi)</math>
<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Numerical Solution
|text=Für die numerische Lösung müssen wir die Bewegungsgleichung
 
::<math>\ddot{\Phi} = f(\Phi)</math>


umschreiben auf eine Differentialgleichung erster Ordnung:
umschreiben auf eine Differentialgleichung erster Ordnung:


<math>\left(\begin{array}{c}\Phi\\\Omega\end{array}\right) = \left(\begin{array}{c}\Omega\\f(\Phi)\end{array}\right)</math>
::<math>\left(\begin{array}{c}\Phi\\\Omega\end{array}\right) = \left(\begin{array}{c}\Omega\\f(\Phi)\end{array}\right)</math>


Diese lösen wir mit einer fertigen Routine - hier einem Runge-Kutta-Verfahren - als Anfangswertproblem. Anders als bei linearisierten Bewegungsgleichungen hängt die Charakteristik der Lösung von den Anfangsbedingungen ab.
Diese lösen wir mit einer fertigen Routine - hier einem Runge-Kutta-Verfahren - als Anfangswertproblem. Anders als bei linearisierten Bewegungsgleichungen hängt die Charakteristik der Lösung von den Anfangsbedingungen ab.
<!-------------------------------------------------------------------------------->
{{MyCodeBlock|title=Numerical Solution
|text=Text
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
/*********************************************************/
/* numerical solution of IVP */
fpprintprec:3$
lhs(dgl[2][1]) = float(expand(rhs(dgl[2][1])));
times : subst([t0 = 0, tmax = 6, dt = 0.01],
                    [t, t0, tmax, dt]);
dgl1stOrder : [Omega,float(expand(rhs(dgl[2][1])))];
stateVabs : [Phi,Omega];
initiVals : [[1,0],[-1,0],
            [0,0.001],[0,3],
            [0,8], [0,10]];
/* solution of IVP (ivs) */
for i: 1 thru length(initiVals) do
  ivs[i] : rk(dgl1stOrder, stateVabs, initiVals[i], times)$
</syntaxhighlight>
</syntaxhighlight>
}}
}}
Zeile 153: Zeile 215:
Das schauen wir uns am Verlauf von Lösungen im Zeitbereich an:
Das schauen wir uns am Verlauf von Lösungen im Zeitbereich an:


 
<table class="wikitable" style="background-color:white; float: left;  margin-right:14px;">
<tr><td>a)[[Datei:K24-21.png|rahmenlos]]
</td><td>[[Datei:Kw24-24.png|rahmenlos]]
b)
</td></tr><tr><td>[[Datei:Kw24-22.png|rahmenlos]]
c)
</td><td>[[Datei:Kw24-25.png|rahmenlos]]
d)
</td></tr><tr><td>[[Datei:Kw24-23.png|rahmenlos]]
e)
</td><td>[[Datei:Kw24-26.png|rahmenlos]]
f)
</td></tr></table>


Lösung f) fällt aus dem Raster - hier gibt es keine periodische Lösung mehr: Bei einer anderen Skalierung der Φ-Achse erkennen wir, dass der Stab AC um den Punkt A rotiert:[[Datei:Kw24-30.png|mini|Lösung f) im Zeitbereich.|alternativtext=|ohne]]
Lösung f) fällt aus dem Raster - hier gibt es keine periodische Lösung mehr: Bei einer anderen Skalierung der Φ-Achse erkennen wir, dass der Stab AC um den Punkt A rotiert:[[Datei:Kw24-30.png|mini|Lösung f) im Zeitbereich.|alternativtext=|ohne]]
Zeile 164: Zeile 238:
Analog sehen Sie in der Lösung c) (grün) eine Separatrix zwischen der periodischen Lösung der kleinen Schwingung um die Gleichgewichtslage und solchen periodischen Lösungen, bei denen der Stab die Horizontale durchbricht: [[Datei:Kw24-32.png|mini|Phasendiagramme|alternativtext=|ohne]]
Analog sehen Sie in der Lösung c) (grün) eine Separatrix zwischen der periodischen Lösung der kleinen Schwingung um die Gleichgewichtslage und solchen periodischen Lösungen, bei denen der Stab die Horizontale durchbricht: [[Datei:Kw24-32.png|mini|Phasendiagramme|alternativtext=|ohne]]
<!-------------------------------------------------------------------------------->
<!-------------------------------------------------------------------------------->


{{MyCodeBlock|title=Post-Processing
{{MyCodeBlock|title=Post-Processing
Zeile 170: Zeile 243:
|code=
|code=
<syntaxhighlight lang="lisp" line start=1>
<syntaxhighlight lang="lisp" line start=1>
1+1
 
/*********************************************************/
/* plot results */
for i: 1 thru length(initiVals) do
(T[i] : makelist(ivs[i][j][1],j,1,length(ivs[i])),
P[i] : makelist(ivs[i][j][2],j,1,length(ivs[i])),
O[i] : makelist(ivs[i][j][3],j,1,length(ivs[i])),
plot2d([discrete, T[i], P[i]], [y,-%pi,%pi],
      [title, sconcat("start @: ",string(initiVals[i]))],
                  [xlabel,"τ/1->"], [ylabel,"Φ/1->"]));
 
/* phase plot */
curves : makelist([discrete,P[i],O[i]],i,1,length(initiVals))$
plot2d(curves, [legend, false], [x,-2,10],
                  [ylabel,"Ω/1->"], [xlabel,"Φ/1->"]);
plot2d(curves, [legend, false], [x,-2,2], [y,-5,5],
                  [ylabel,"Ω/1->"], [xlabel,"Φ/1->"]);
plot2d(curves, [legend, false], [x,-2,2], [y,-5,5], [svg_file, "Kw24-20.svg"],
                  [ylabel,"Ω/1->"], [xlabel,"Φ/1->"]);
</syntaxhighlight>
</syntaxhighlight>
}}
}}
<table class="wikitable" style="background-color:white; float: left;  margin-right:14px;">
<tr><td>[[Datei:K24-21.png|rahmenlos]]
a)
</td><td>[[Datei:Kw24-24.png|rahmenlos]]
b)
</td></tr><tr><td>[[Datei:Kw24-22.png|rahmenlos]]
c)
</td><td>[[Datei:Kw24-25.png|rahmenlos]]
d)
</td></tr><tr><td>[[Datei:Kw24-23.png|rahmenlos]]
e)
</td><td>[[Datei:Kw24-26.png|rahmenlos]]
f)
</td></tr></table>


<hr/>
<hr/>

Version vom 29. März 2021, 06:02 Uhr


Aufgabenstellung

Ein Mechanismus besteht aus einem starren Stab AC (Länge 1, Massenmoment J) und einer Feder BC (Länge 2 im entspannten Zustand, Steifigkeit k).

Die Körper sind drehbar gelagert, der Abstand AB ist 3.


Lageplan

Hier geht es um die verschiedenen Lösungstypen des Systems - in Abhängigkeit von den Anfangsbedingungen. Gesucht sind Lösungen als Anfangswertproblem, nichtlineare Schwingungen und Aussagen über die Stabilität der Bewegungen.

Lösung mit Maxima

Header

Gesucht ist die Lösung zu dem nichtlinearen Anfangswertproblem.


/*********************************************************/
/* MAXIMA script                                         */
/* version: wxMaxima 15.08.2                             */
/* author: Andreas Baumgart                              */
/* last updated: 2017-09-30                              */
/* ref: Kw24 (TM-C, Labor 6)                             */
/* description: finds the solution for                   */
/*              the nonlinear IVP                        */
/*********************************************************/




Declarations

Wir wählen die Längenverhältnisse

so dass die Feder für φ(t) = 45° entspannt ist.


/*********************************************************/
/* declare variational variables - see 6.3 Identifiers   */
declare("δW", alphabetic); /* total virtual work         */
declare( "φ", alphabetic); /* coordinate of rotatin      */
declare("δφ", alphabetic); /* variation of coordinate phi*/
declare("δw", alphabetic); /* variation of coordinate w  */

/* parameter */
assume(l[1]>0, k>0, J>0);
params: [l[2] = l[1], l[3] = sqrt(2)*l[1]];




Kinematics

Der Mechanismus hat genau einen Freiheitsgrad, den wir entweder durch die Federlängung w oder den Drehwinkel φ erfassen. Den Zusammen zwischen den beiden stellen wir über

her. Nun steht uns frei, ob wir diese Beziehung nach w oder φ auflösen. Wir wählen φ als Minimal-Koordinate und erhalten

Nach den Regeln der Variationsmethoden erhalten wir daraus:

.

/*********************************************************/
/**** kinematics ***/

/* coordinates */
coords : [φ, δφ];

kin : (l[1]*sin(φ(t)))^2+(l[3]-l[1]*cos(φ(t)))^2=(l[2]+w)^2;
kin : trigsimp(ratsimp(solve(subst(params, kin),w)[2]));
kin : [kin,
       δw = subst([epsilon=0],
                 diff(subst([φ(t)=φ(t)+epsilon*δφ],
                                  subst(kin,w)),epsilon))];




Equilibrium Conditions

Mit dem Prinzip der virtuellen Verrückungen schreiben wir

und erhalten die Bewegungsgleichung


/*********************************************************/
/* equilibrium conditions */
PvV : δW = -J*'diff(φ(t),t,2)*δφ - k*w*δw;

/* dgl for φ ist the coefficient of δφ in δW */
dgl[1] : expand(coeff(subst([δW = 0], subst(kin, -PvV)),δφ));




Equilibrium Points

Die Bewegungsgleichung machen wir dimensionslos, indem wir charakteristische Größen aus der linearisierten Bewegungsgleichung verwenden.

Die drei Gleichgewichts-Lagen bzgl. derer wir die Bewegungsgleichung linearisieren können, sind

.

Eine Taylor-Reihenentwicklung des Rückstellmoments bzgl. φ0,3 ,also

liefert die linearisierte Bewegungsgleichung

.

Dabei ist ω0 die Eigenkreisfreuquenz der linearisierten Bewegungsgleichung.


/*********************************************************/
/* Zeros of restoring function */
equilib : solve(subst(['diff(φ(t),t,2)=0],dgl[1]),φ(t));

/*********************************************************/
/* derive the eigenfrequency                             */
/* of linearized ODE at stable point of equilob.         */
LinODE : taylor (subst([φ(t)=psi],rhs(dgl[1]/J)), psi, %pi/4,1);
ref : omega[0]^2 = coeff(expand(LinODE),psi);
ref : [t[B] = 1/(subst(solve(ref,omega[0])[2],omega[0])/(2*%pi))];




Make Equations of Motion Dimensionless

Eine Bezugszeit für unsere nichtlineare Bewegungsgleichung kommt nun aus

Damit lautet die dimensionslose Bewegungsgleichung

Rückstellmoment der Feder

Den nichtlinearen (zweiten) Summand können wir als dimensionsloses "Rückstellmoment auffassen und über dem Drehwinkel auftragen:

Nullstellen des Rückstellmoments zeigen die Gleichgewichtslagen des Systems an. Die Steigung dieses Rückstellmoments in den Gleichgewichtslagen deutet dabei auf die Stabilität der linearisierten Bewegung (oder Bewegung im Kleinen) hin. So ist die Lösung Φ = 0 eine instabile Gleichgewichtslage, die beiden anderen sind stabil.


/*********************************************************/
/* make dgl dimensionless employing the eigenfrequency   */

dgl[2] : subst([diff(φ(t),t,2) = 'diff(Phi,t,2)/t[B]^2,
                     φ(t)      =      Phi ] , dgl[1]);
dgl[2] : solve(subst(ref,dgl[2]), 'diff(Phi,t,2));

/* equilibrium angles für Φ */
plot2d(-rhs(dgl[2][1]), [Phi,-1,1],
                        [legend, "Rückstell'moment'"],
                        [xlabel, "Φ/1->"], [ylabel, "M/1->"]);




Numerical Solution

Für die numerische Lösung müssen wir die Bewegungsgleichung

umschreiben auf eine Differentialgleichung erster Ordnung:

Diese lösen wir mit einer fertigen Routine - hier einem Runge-Kutta-Verfahren - als Anfangswertproblem. Anders als bei linearisierten Bewegungsgleichungen hängt die Charakteristik der Lösung von den Anfangsbedingungen ab.


/*********************************************************/
/* numerical solution of IVP */
fpprintprec:3$
lhs(dgl[2][1]) = float(expand(rhs(dgl[2][1])));
		
times : subst([t0 = 0, tmax = 6, dt = 0.01],
                    [t, t0, tmax, dt]);
dgl1stOrder : [Omega,float(expand(rhs(dgl[2][1])))];
stateVabs : [Phi,Omega];
initiVals : [[1,0],[-1,0],
             [0,0.001],[0,3],
             [0,8], [0,10]];
/* solution of IVP (ivs) */
for i: 1 thru length(initiVals) do
   ivs[i] : rk(dgl1stOrder, stateVabs, initiVals[i], times)$




tmp

Das schauen wir uns am Verlauf von Lösungen im Zeitbereich an:

a)

b)

c)

d)

e)

f)

Lösung f) fällt aus dem Raster - hier gibt es keine periodische Lösung mehr: Bei einer anderen Skalierung der Φ-Achse erkennen wir, dass der Stab AC um den Punkt A rotiert:

Lösung f) im Zeitbereich.

Oft erkennt man das Charakteristische einer Lösung im Phasendiagramm. Hier wird die Geschwindigkeit Ω über dem Drehwinkel Φ aufgetragen:

Phasendiagramme der Lösungen

Periodische Lösungen - die oft gesucht sind - erkennt man an einem geschlossenen Verlauf. Die hellblaue Lösung gehört zu den Anfangswerten von f), bei dem der Drehwinkel immer größer wird - sich also keine periodische Lösung ergibt.

Wir sehen, dass die Lösung e) (schwarz) eine "gerade noch" periodische Lösung ist - etwas mehr Anfangsgeschwindigkeit im System und die Lösung sieht wie in f) aus. Diese Lösung nennt man auch Separatrix - sie separiert zwei Charakteristiken von Lösungen.

Analog sehen Sie in der Lösung c) (grün) eine Separatrix zwischen der periodischen Lösung der kleinen Schwingung um die Gleichgewichtslage und solchen periodischen Lösungen, bei denen der Stab die Horizontale durchbricht:

Phasendiagramme

Post-Processing

Text


/*********************************************************/
/* plot results */
for i: 1 thru length(initiVals) do
	(T[i] : makelist(ivs[i][j][1],j,1,length(ivs[i])),
	 P[i] : makelist(ivs[i][j][2],j,1,length(ivs[i])),
	 O[i] : makelist(ivs[i][j][3],j,1,length(ivs[i])),
	 plot2d([discrete, T[i], P[i]], [y,-%pi,%pi],
	      [title, sconcat("start @: ",string(initiVals[i]))],
	                   [xlabel,"τ/1->"], [ylabel,"Φ/1->"]));

/* phase plot */
curves : makelist([discrete,P[i],O[i]],i,1,length(initiVals))$
plot2d(curves, [legend, false], [x,-2,10],
                   [ylabel,"Ω/1->"], [xlabel,"Φ/1->"]);
plot2d(curves, [legend, false], [x,-2,2], [y,-5,5],
                   [ylabel,"Ω/1->"], [xlabel,"Φ/1->"]);
plot2d(curves, [legend, false], [x,-2,2], [y,-5,5], [svg_file, "Kw24-20.svg"],
                   [ylabel,"Ω/1->"], [xlabel,"Φ/1->"]);





Links

  • ...

Literature

  • ...