Gelöste Aufgaben/Bike

Aufgabenstellung

Ein Hersteller von Scheibenbremsbelägen erhält wiederholt Reklamationen, weil Fahrräder beim Bremsen erhebliche Geräusche entwickeln oder die Gabel starke Schwingungen in der Laufradebene ausführt.

Vorderradgabel mit Scheibenbremse.

Gesucht ist ein mathematisches Modell, das die Entstehung von Geräuschen und Schwingungen an einer Vorderradgabel zu erklären hilft.

Wir simulieren dazu das dynamische Verhalten eines Vorderrades beim Bremsen mit Scheibenbremse.

Strukturieren

Wir suchen einen Mechanismus, bei dem sich eine Schwingung ohne äußere Erregung einstellt. In der Mechanik spricht man dann von einer selbsterregten Schwingung. Bei rein mechanischen Systemen wie unserem Vorderrad sind dabei oft Stick-Slip-Schwingungen die Ursache, bei denen z.B. die Bremse kurzzeitig auf der Bremsscheibe blockiert und sich dann wieder losreißt. Das charakteristische Merkmal, das diesen Mechanismus treibt, ist dabei eine fallende Kennlinie: bei fast allen Materialpaarungen ist die Haftkraft der Bremse größer als die Reibkraft.

Aufgabenstellung konkretisieren

Wir erwarten, den gesuchten Mechanismus nur anhand der Dynamik des Vorderrades beschreiben zu können. Das Vorderrad wird deshalb am Steuersatz vom Rahmen getrennt und dort – gedacht - in einer starren Führung mit konstanter horizontaler Geschwindigkeit v₀ geführt. Um das Modell einfach zu halten, untersuchen wir Schwingungen nur in der Radebene. Das ist eine wesentliche Vereinfachung, weil wir bei Schwingungen unter Einwirkung der Bremse von Bewegungen des Rades um die Steuerrad-Achse rechnen müssen. Diese Vernachlässigen wir. Eine weitere Vereinfachung führen wir durch Setzen des Steuerwinkels von 90° ein. Damit vernachlässigen wir den Beitrag einer vertikalen Bewegung: beim Schwingen der Gabel in der Radebene bleibt der Steuersatz unverändert in der Ausgangshöhe. Diese Auf- und Ab-Bewegung berücksichtigen wir also nicht.

System strukturieren

Unser verbleibendes System besteht nun aus

• Federgabel (Länge $H$) und Bremssattel (Abmessungen $h$ und $b$),
• Nabe (Außenradius $r_1$)mit Bremsscheibe (Radius $r_2$)
• Felge (Innenradius $R_2$)und Mantel (Außenradius $R_1$) sowie den
• Speichen.

Vorderradgabel mit Scheibenbremse.

Dieses System wollen wir als dynamisches System mit Bewegungen nur in der Radebene beschreiben.

Im folgenden müssen wir Annahmen zu den Systemkomponenten treffen, die unseren Vorstellung von den relevanten physikalischen Prozessen entsprechen. Dabei gehen wir nach dem Prinzip

• "keep it dead simple"

vor - was allerdings interpretationsbedürtig ist. Ein erster Iterationsschritt ist:

• "so einfach wie möglich - so kompliziert wie nötig".

Wir versuchen also, die Anzahl der Freiheitsgrade des Systems so gering wie möglich zu halten, ohne dabei das Ziel aus den Augen zu verlieren: die Simulation eines komplexen Vorgangs.

So gehe ich davon aus, dass die selbsterregte Schwingung des Systems ihren Grund in Stick-Slip-Schwingungen an der Bremse hat. Das bedeutet: bei niedrigen Fahrtgeschwindigkeiten $v_0$ wäre die Relativgeschwindigkeit zwischen Bremsscheibe und Bremsbacken Null. Das kann nur aus einer elastischen Verformung der Gabel resultieren.

Der Mechanismus könnte ungefähr so aussehen:

Beim Bremsen verformt sich die Gabel elastisch, die Nabe wandert - relativ zum Steuersatz - in Fahrtrichtung nach hinten. Damit der Mantel dabei weiter in Punkt A abrollt, reduziert sich seine Drehgeschwindigkeit und es kommt zu einer geringeren Relativgeschwindigkeit an der Bremse. Je nach Betriebsbedingungen kann hier die Relativgeschwindigkeit Null werden, die Bremsebacken haften an der Bremsscheibe und ein Stick-Slip-Mechanismus nimmt seinen Anfang. Dabei wird offensichtlich auch ein veränderliches Drehmoment zwischen der Nabe (Bremsscheibe) und Felge übertragen.

Aus diesen Anfangsüberlegungen leiten wir folgende Annahmen für die Modellierung ab:

• Als starre Körper beschrieben wir die Teilkomponenten
  • Felge & Mantel sowie,
  • Nabe & Bremsscheibe.
• Als rein elastischen Körper beschreiben wir die Speichen, die eine Verdrehung zwischen Nabe und Felge erlauben.
• Als Kontinuum mit Masse und Elastizität beschrieben wir die Federgabel. Diese hat einen über ihre Länge veränderlichen Querschnitt.

Die Herausforderung hier liegt darin, die relevanten Mechanismen nicht durch unpassende Annahmen zu unterbinden. So würde eine als starr modellierte Gabel hier vermutlich nicht zu den erhofften Schwingungserscheinungen in einer Simulation führen.

Wenn man also nicht zu zufriedenstellenden Simulationsergebnissen kommt, muss man an dieser Stelle wieder neu einsteigen und das Modell ändern.

Für die Modellierung unseres System zeichnmen sich zwei zentrale Herausforderungen ab:

1. die Gleichgewichtsbedingungen für die Gabel und
2. die Kennlinie für die Bremse und den Kontakt Mantel/Asphalt.

Für die Bremskraft als Funktion der Relativgeschwindigkeit müssen wir dann passende Kennlinien ansetzen.

Und für die Gabel wählen wir einen Euler-Bernoulli-Balken mit veränderlichem Querschnitt und schreiben die Gleichgewichtsbedingungen mit der Methode der Finiten Elemente an. Mindestens brauchen wir zwei Finite Elemente - vom Steuersatz zur Bremse und von der Bremse zur Nabe.

Als Referenzpunkte wählen wir deshalb

• $A$: den Aufstandspunkt des Rades auf dem Asphalt,
• $B$: den Berührungspunkt zwischen Bremsscheibe und Bremsbelag,
• $C$: die Achse des Rades, die Naben und Gabel verbindet sowie
• $N_0,N_1,N_2$: als Punkte entlang der x-Achse der Federgabel.

Um die Auslenkung der Federgabel zu erfassen, führen wir materielle Koordinaten $x,z$ nach Vorgabe des Euler-Bernoulli-Balkens ein. Den Ursprung des Koordinatensystems legen wir passend in die Wurzel der Gabel bei $N_0$.

Der Einfachheit halber sollen $C$ und $N_2$ zusammen fallen. Die Gabel erfassen wir also mit zwei Finiten Elementen:

1. von $N_0$ nach $N_1$ und
2. von $N_1$ nach $N_2$.

Vorderradgabel mit Scheibenbremse.

Wir trennen nun zunächst unser System in Gabel und Vorderrad und führen Koordinaten der Bewegung für sie ein. Gleichzeitig führen dabei die Naben-Schnittkäfte $C_x,C_z$, die Bremskraft $B$ und die Kontaktkräfte $A_x,A_z$ ein.

Koordinaten der Gabel für $N_1$.

Die Bewegung der Gabel erfassen wir zunächst an den Punkten $N_0,N_1,N_2$ mit den Koordinaten der Bewegung $W_i(t)$ (horizontale Auslenkung) und ${\mathrm{\Phi }}_i(t)$ (Kippen des Querschnitts) für $i=0,1,2$. Das ist für Knoten 1 rechts dargestellt. Offensichtlich ist bei einer festen Einspannung im Steuersatz $W_0(t)\equiv 0$ und ${\mathrm{\Phi }}_0(t)\equiv 0$ - was wir wie üblich bei der Methode der Finiten Elemente erst bei der Einarbeitung der Randbeidngungen berücksichtigen.

Für die Drehbewegung wollen des Vorderrades wollen wir zwei relevante Koordinaten einführen: die Drehung der Komponente Nabe & Bremsschreibe sowie die Drehung der Modellkomponenten Mantel & Felge. Um die Beschleunigung des Steuersatzes (und des Rest-Fahrrades) nicht berücksichtigen zu müssen, gehen wir von konstantem $v_0$ aus. Dann können wir die Verdrehung des Vorderrades aus den einer großen Bewegung mit konstanter Drehgeschwindigkeit und einer Abweichung davon erfassen:

Freikörperbild Vorderrrad.

Als Koordinaten der Drehbewegung führen wir

• für die Rollbewegung $R_1\mathrm{\Omega }\cdot t=v_0\cdot t$ als Referenzgeschwindigkeit ein und relativ dazu
• den Winkel ${\mathrm{\Psi }}_1(t)$ für Felge/Mantel und
• den Winkel ${\mathrm{\Psi }}_2(t)$ für Nabe/Bremsscheibe.

Die horizontale Translationsbewegung des Rades erfassen wir über die Nabe - seine Auslenkung ist dann $W_2(t)$.

Wir sammeln die Koordinaten der Bewegung für unser System nun in der Spaltenmatrix

$\underset{\_}{\tilde{Q}}=\left(\begin{array}{c}{W}_0(t)\\ {\mathrm{\Phi }}_0(t)\\ W_1(t)\\ {\mathrm{\Phi }}_1(t)\\ W_2(t)\\ {\mathrm{\Phi }}_2(t)\\ {\mathrm{\Psi }}_1(t)\\ {\mathrm{\Psi }}_2(t)\end{array}\right)$

---

Modellieren

some text

Zunächst: ein kurzer Blick auf den roten Faden durch die kommenden Abschnitte: Unser mathematischen Modell wird wegen der festen Einspannung der Gabel im Steuersatz noch sechs Koordinaten $\underset{\_}{Q}=\left(\begin{array}{c}{W}_1\\ {\mathrm{\Phi }}_1\\ W_2\\ {\mathrm{\Phi }}_2,\\ {\mathrm{\Psi }}_1,\\ {\mathrm{\Psi }}_2\end{array}\right)$ haben. Wir werden ein System von sechs gewöhnlichen, nichtlinearen Differentialgleichungen zweiter Ordnung erhalten $\underset{\_}{\underset{\_}{M}}\cdot \ddot{\underset{\_}{Q}}+\underset{\_}{\underset{\_}{K}}\cdot \underset{\_}{Q}=\underset{\_}{P}(\underset{\_}{Q},\dot{\underset{\_}{Q}})$, wobei $\underset{\_}{\underset{\_}{M}},\underset{\_}{\underset{\_}{K}}$ die Massen- bzw. Steifigkeitsmatrix sind und wir in $\underset{\_}{P}$ die nichtlinearen Anteile aus den Brems-Kennlinien sammeln.

Schnittmomente statischer Fall

Bevor wir das mathematische Modell anschreiben, berechnen wir uns als Referenz die Schnittgrößen im Euler-Bernoulli-Balken für den statischen Fall. Anhand dieser Referenz können wir später unser FE-Modell auf Plausibilität prüfen.

Freikörperbild Gabel.

Als Ausgangspunkt geben wir uns eine Bremskraft $B$ in Punkt B vor und schrieben die Gleichgewichtsbedingungen für das Vorderrad an, aus denen wir z.B. $A_H,C_H$ berechnen.

Dafür benötigen wir die Bremskraft $B$ aufgeteilt in die Komponenten

$\underset{\_}{B}=\left(\begin{array}{c}{B}_x\\ B_z\end{array}\right)$,

wobei

$B_x=B\cdot \sin (\alpha ),B_z\cdot \cos (\alpha )$ mit $\cos (\alpha )=h/r_2,\sin (\alpha )=b/r_2$.

So kommen aus den Momentengleichgewichten

$\begin{array}{l}{R}_1\cdot A_H=r_2\cdot B\\ (h+R_1)\cdot B_z+b\cdot B_x+R_1\cdot C_H=0\end{array}$

und $\underset{\_}{B}$

$\begin{array}{l}{A}_H=B\frac{\left(r_2\right)}{R_1},\\ C_H=-B\frac{h^2+R_{1}h+b^2}{R_1{r}_2}\end{array}$

Freikörperbild für die Schnittlasten.

Mit den Schnittbildern für die Bereiche $N_0-N_1$ und $N_1-N_2$ rechts und den zugehörigen Gleichgewichtsbedingungen können wir die Schnittmomente $M_1(x_1),M_2(x_2)$ berechnen. Mit ${\xi }_i=x_i/{\ell }_i$ ist

$\begin{array}{lll}{M}_1(x_1)& =r_{2}B\cdot & \frac{(R_1{h}^2-R_1{\ell }_{2}h-r_2^2{\ell }_2+(({\xi }_1-1){\ell }_1-R_1)r_2^2)}{(R_1{r}_2^2)}\\ M_2(x_2)& =r_{2}B\cdot & \frac{((R_1{\xi }_2-R_1){\ell }_{2}h+(r_2^2{\xi }_2-r_2^2){\ell }_2)}{(R_1{r}_2^2)}\end{array}$

Die Verläufe der Schnittlasten tragen wir dann auf:

Verlauf des Schnittmoments.

Mit dieser Rechnung haben wir eine Referenzlösung geschaffen, die wir später mit den FE-Ergebnissen vergleichen können.

/*********************************************************/
/* MAXIMA script                                         */
/* version: wxMaxima 21.05.2                             */
/* author: Andreas Baumgart                              */
/* last updated: 2025-03-27                              */
/* ref: self-excitaation of bike-breakes                 */
/* description: analytische Lösung für die Schnittgrößen */
/*********************************************************/

assume(mm>0);
/* paramteres */
params: [ℓ[1]=H-h,
         ℓ[2]=h,
         r[2] = 100*mm,
         R[1] = 350*mm,
         H = 375*mm,
         h =  90*mm];

/* equilibrium conditions components */
GGW: [A[H]*R[1] = B*r[2],
      R[1]*C[H] + (R[1]+h)*B[z] + b*B[x] = 0,
      B[x] = B*sin(α), B[z] = B*cos(α),
      cos(α) = h/r[2], sin(α) = b/r[2]];

/* equilibrium conditions for internal loads Q, M */
ggw: [Q[2] = -C[H],
      M[2] = (ℓ[2]-x[2])*C[H],
      Q[1] = -C[H]-B[z],
      M[1] = (ℓ[2]+ℓ[1]-x[1])*C[H] + (ℓ[1]-x[1])*B[z]-b*B[x]];

rep: [h^2 = r[2]^2 - b^2, h = ℓ[2]]; 

sol: solve(GGW, [A[H],C[H],B[x],B[z],cos(α),sin(α)])[1];
sol: append(sol,ratsimp(subst(sol, ggw)));

/* use dimensionless coordiantes */
dimless:  [x[1] = ℓ[1]*ξ[1],
           x[2] = ℓ[2]*ξ[2]];

/* cross-sectional moments */
mom: ratsimp(subst([b^2 = r[2]^2-h^2],expand(subst(dimless,subst(sol,[M[1],M[2]])))));

/* plot section-wise */
mmm: ratsimp(mom/(r[2]*B));
mmm: ratsimp(subst([ξ[1]=ξ,ξ[2]=ξ],subst(params, mmm)));

plot2d([[parametric, subst(params,     ℓ[1]*ξ)/mm, mmm[1], [ξ,0,1]],
        [parametric, subst(params,ℓ[1]+ℓ[2]*ξ)/mm, mmm[2], [ξ,0,1]]],
        [legend,"sec- I","sec-II"], [xlabel, "x/mm ->"], [ylabel, "M/r_2*B ->"])$

---

Nun liegt die Idee nahe, aus dem Momentenverlauf mit der Gleichung

$E\cdot I\cdot w^{″}(x)=-M(x)$

die Verschiebung $w(x)$ und damit eine Nachgibigkeit zu berechnen, also

$w(x)=-\int \int \frac{M(x)}{EI}\text{d}x\text{d}x$

zu lösen. Das funktioniert nicht gut, weil für die Gabel mit Hohlprofil

$I(x)=\frac{1}{64}(D^4(x)-d^4(x))$

gilt. Dabei sind $D$ und $d$ Funktionen von $x$, so dass wir diesen Weg aufgeben müssen und eine Näherungslösung - hier die Methode der Finiten Elemente - verfolgen.

Mathematisches Modell formulieren

Die Bewegungsgleichungen des Systems leiten wir mit dem Prinzip der virtuellen Verrückungen her:

$\begin{array}{lcc}\delta W& =& \delta W^a-\delta \mathrm{\Pi }\\ & \stackrel{!}{=}& 0\end{array}$

Dabei ist

• $\delta W^a$ ... die virtuelle Arbeit äußerer, eingeprägter Kräfte - und dazu zählen nach dem Prinzip von d'Alembert auch die Trägheitskräfte der Systemkomponenten;
• $\delta \mathrm{\Pi }$ ... die virtuelle Formänderungsenergie.

Die Aussage ist: die Summe der virtuellen Arbeiten äußerer Kräfte und der virtuellen Formänderungsenergie ist Null.

Nach dem Prinzip der virtuellen Verrückungen dürfen wir jede Minimalkoordinate aus $\tilde{Q}$ variieren - an ihr wackeln - und wir definieren als Komplement zu $\tilde{Q}$ die Spaltenmatrix

$\underset{\_}{\delta \tilde{Q}}=\left(\begin{array}{c}\delta W_0\\ \delta {\mathrm{\Phi }}_0\\ \delta W_1\\ \delta {\mathrm{\Phi }}_1\\ \delta W_2\\ \delta {\mathrm{\Phi }}_2,\\ \delta {\mathrm{\Psi }}_1,\\ \delta {\mathrm{\Psi }}_2\end{array}\right)$.

Und ja: am Steuersatz ist $W_0\equiv 0,{\mathrm{\Phi }}_0\equiv$, also auch $\delta W_0\equiv 0,{\mathrm{\Phi }}_0\equiv 0$ - aber das berücksichtigen erst mit der Einarbeitung der Randbedingungen.

Die Kräfte und Momente des Systems leisten nun virtuelle Arbeiten auf ihren virtuellen Verrückungen.

Virtuelle Formänderungsenergie

Spannungen leisten virtuelle Formänderungsarbeit auf vituellen Dehungen, also

$\delta \mathrm{\Pi }=\underset{V}{\int }{\underset{\_}{\sigma }}^T\cdot \delta \underset{\_}{ϵ}\text{d}x$.

Wir verwenden hier nur Bauteile, bei denen wir aus Formelwerken $\delta {\mathrm{\Pi }}_j$ als Funktion der Koorinaten kennen.

Die virtuelle Formänderungsenergie setzt sich additiv zusammen aus den Anteilen $\delta {\mathrm{\Pi }}_G$ der Gabel und $\delta {\mathrm{\Pi }}_S$ für die Drehfeder der Speichen:

$\delta \mathrm{\Pi }=\delta {\mathrm{\Pi }}_G+\delta {\mathrm{\Pi }}_S$.

Weil wir nur zwei Finite Element berücksichtigen, ist

$\delta {\mathrm{\Pi }}_G=\delta {\mathrm{\Pi }}_{G1}+\delta {\mathrm{\Pi }}_{G2}$

Für die Speichen schreiben wir als diskrete Feder direkt

$\delta {\mathrm{\Pi }}_S=K_S\cdot ({\mathrm{\Psi }}_2-{\mathrm{\Psi }}_1)\cdot (\delta {\mathrm{\Psi }}_2-\delta {\mathrm{\Psi }}_1)$

mit der Drehfedersteifigkeit $K_S$an.

Die virtuelle Formänderungsenergie eines Euler-Bernoulli-Balkens kommt aus der virtuellen Arbeit des Biegemoments $M(x)$ auf der virtuellen Krümmung $\delta w^{″}(x)$, also für Element i

$\delta {\mathrm{\Pi }}_{Gi}=\underset{0}{\overset{{\ell }_i}{\int }}{M}_i(x)\cdot \delta {w_i}^{″}(x)dx$

mit

$M_i(x)=E\cdot I(x)\cdot {w_i}^{″}(x)$ und

$(.):=\frac{d(.)}{dx}$.

In dieser Gleichung stehen nun das Flächenmoment 2ten Grades $I_i(x)$ und die Auslenkung $w_i(x)$ als Funktionen von $x$. Diese beiden Funktoinen müssen wir Element-weise definieren, die Ableitungen und Integrale lösen und nach den Koordinaten auseinandersortieren - fertig.

Querschnitt.

Wir starten mit

$I(x)=\frac{1}{64}(D^4(x)-d^4(x))$

und wählen für $D_i(x),d_i(x)$ lineare Ansatzfunktionen mit den Durchmessern ${\mathbf{𝐃}}_0,{\mathbf{𝐃}}_1,{\mathbf{𝐃}}_2$ sowie ${\mathbf{𝐝}}_0,{\mathbf{𝐝}}_1,{\mathbf{𝐝}}_2$ in den Knoten $N_0,N_1,N_2$. Dann ist in Element i

$D_i(x)={\mathbf{𝐃}}_{i-1}\cdot (1-{\xi }_i)+{\mathbf{𝐃}}_i\cdot {\xi }_i$

und

$d_i(x)={\mathbf{𝐝}}_{i-1}\cdot (1-{\xi }_i)+{\mathbf{𝐝}}_i\cdot {\xi }_i$

Ähnlich verfahren wir mit der Ansatzfunktionen für die Auslenkung und setzen

$w_i(x)={\underset{\_}{Q}}_i^T\cdot \underset{\_}{\phi }(x)$

an. Hier stehen in

${\underset{\_}{\tilde{Q}}}_i=\left(\begin{array}{c}{W}_{i-1}\\ {\mathrm{\Phi }}_{i-1}\\ W_i\\ {\mathrm{\Phi }}_i\end{array}\right)$

die Koordinaten der Knoten eines Elements.

Die Trial-Functions sind Hermite-Polynome wie in der FEM-Formulierung für den Euler-Bernoulli-Balken beschreiben:

$\underset{\_}{\phi }(x)=\left(\begin{array}{c}(\xi -1{)}^2\cdot (2\cdot \xi +1)\\ {\ell }_i\cdot \xi \cdot (\xi -1{)}^2\\ -{\xi }^2\cdot (2\xi -3)\\ {\ell }_i\cdot {\xi }^2\cdot (\xi -1)\end{array}\right)$

Einsetzen und differenzieren / integrieren der Terme liefert diesen langen Ausdruck

$\delta {\mathrm{\Pi }}_{Gi}=\underset{0}{\overset{{\ell }_i}{\int }}E\cdot I_i(x)\cdot (W_{i-1}{{\phi }_1}^{″}\cdot \delta W_{i-1}{{\phi }_1}^{″}+{\mathrm{\Phi }}_{i-1}{{\phi }_2}^{″}\cdot \delta W_{i-1}{{\phi }_1}^{″}+W_i{{\phi }_3}^{″}\cdot \delta W_{i-1}{{\phi }_1}^{″}+\dots +{\mathrm{\Phi }}_i{{\phi }_4}^{″}\cdot \delta {\mathrm{\Phi }}_i{{\phi }_4}^{″})$

den wir jedoch zweckmäßiger Weise als

$\delta {\mathrm{\Pi }}_{Gi}=\delta {\underset{\_}{\tilde{Q}}}_i^T\cdot {\underset{\_}{\underset{\_}{\tilde{K}}}}_i\cdot {\underset{\_}{\tilde{Q}}}_i$

hinschrieben. Dabei sind die Matrix-Koeffizienen von ${\underset{\_}{\underset{\_}{K}}}_i$

$k_{i,jk}=\underset{0}{\overset{{\ell }_i}{\int }}E\cdot I_i(x){{\phi }_j}^{″}\cdot {{\phi }_k}^{″}dx$

Hier ist es zweckmäßig, zu dimensionslosen Koordinaten

${\xi }_1=\frac{x}{{\ell }_1}$

${\xi }_2=\frac{(x-H)}{{\ell }_2}$

zu wechseln, wir müssen aber beim Ableiten und Integrieren

$\frac{d(.)}{dx}=\frac{d(.)}{d{\xi }_i}\cdot \frac{1}{{\ell }_i}$

berücksichtigen. Die virtuelle Formänderungsarbeit des Gesamtsystems erhalten wir durch Aufsummieren aller Anteile zu

$\delta {\mathrm{\Pi }}_G=\delta {\underset{\_}{\tilde{Q}}}^T\cdot \underset{\_}{\underset{\_}{\tilde{K}}}\cdot \underset{\_}{\tilde{Q}}$

mit

$\underset{\_}{\underset{\_}{\tilde{K}}}=\left(\begin{array}{cccccccc}{k}_{1,11}& k_{1,12}& k_{1,13}& k_{1,14}& 0& 0& 0& 0\\ k_{1,12}& k_{1,22}& k_{1,23}& k_{1,24}& 0& 0& 0& 0\\ k_{1,13}& k_{1,23}& k_{1,33}+k_{2,11}& k_{1,34}+k_{2,12}& k_{2,13}& k_{2,14}& 0& 0\\ k_{1,14}& k_{1,24}& k_{1,34}+k_{2,12}& k_{1,44}+k_{2,22}& k_{2,23}& k_{2,24}& 0& 0\\ 0& 0& k_{2,13}& k_{2,23}& k_{2,33}& k_{2,34}& 0& 0& \\ 0& 0& k_{2,14}& k_{2,24}& k_{2,34}& k_{2,44}& 0& 0& \\ 0& 0& 0& 0& 0& 0& K_S& -K_S\\ 0& 0& 0& 0& 0& 0& -K_S& K_S\\ \end{array}\right)$

Virtuelle Arbeiten der d'Alembertschen Trägheitskräfte

Analog zur virtuellen Formänderungsenergie berechnen wir nun die virtuelle Arbeit äußerer Kräfte aus den Beiträgen $\delta W_A^a$ von d'Alembertschen Trägheitskräften und der Kontaktkräfte $\delta W_F^a$:

$\delta W^a=\delta W_A^a+\delta W_F^a$

Äußere, eingeprägte Kräfte - hier $\underset{\_}{A},\underset{\_}{C},\underset{\_}{B}$ - leisten virtuelle Arbeit auf vituellen Verschiebungen, also

$\delta W_F^a=\sum _i{\underset{\_}{F}}_i^T(t)\cdot \delta {\underset{\_}{r}}_i$,

wobei $\delta {\underset{\_}{r}}_i$ die Variation der Verschiebung des Kraftangriffspunktes ist. Und äußere, eingeprägte Volumenkräfte - hier die d'Alembert'schen Trägheitskräfte - leisten virtuelle Arbeit auf den vituellen Verschiebungen, also mit der Dichte $\varrho$

$\delta W_A^a=\underset{V}{\int }-\varrho {\underset{\_}{\ddot{r}}}_i^T(\underset{\_}{x},t)\cdot \delta {\underset{\_}{r}}_i(\underset{\_}{x})\text{d}\underset{\_}{x}$.

Aus diesen Ansätzen erhalten wir die Koeffizienten der Massenmatrix $\underset{\_}{\underset{\_}{\tilde{M}}}$ und der "rechten Seite" $\underset{\_}{\tilde{P}}(\underset{\_}{\dot{\tilde{Q}}},\underset{\_}{\tilde{Q}})$.

Wir starten mit

$\begin{array}{llll}\delta W_A^a& =& -& (M_f+M_n)\cdot {\underset{\_}{\ddot{r}}}_C^T\cdot \delta {\underset{\_}{r}}_C\\ & =& -& J_f{\ddot{\mathrm{\Psi }}}_1\cdot \delta {\mathrm{\Psi }}_1\\ & =& -& J_n{\ddot{\mathrm{\Psi }}}_2\cdot \delta {\mathrm{\Psi }}_2\\ & =& -& \underset{x}{\int }\varrho A(x)\ddot{w}(x,t)\cdot \delta w(x)\text{d}x\\ \end{array}$,

wobei

• $M_f,M_n$ die Massen von Felge & Mantel bzw. Nabe & Bremsscheibe sowie,
• $J_f,J_n$ die Massenmomente von Felge & Mantel bzw. Nabe & Bremsscheibe sind und
• $\varrho A(x)$ die Massenbelegung der Gabel ist.

Bei den Anteilen der Gabel steht dabei nur das $\ddot{w}(x,t)$, weil sich der Steuersatz mit konstanter Geschwindigkeit bewegt. Das untersuchen wir etwas genauer für die Beschleunigung ${\underset{\_}{\ddot{r}}}_C$.

Die Koordinaten des Ortsvektors der Nabe C in $x,z$-Richtung schrieben wir als

${\underset{\_}{r}}_C(t)=\left(\begin{array}{c}0\\ -v_0\cdot t\end{array}\right)+\left(\begin{array}{c}H\\ W_2(t)\end{array}\right)$

und wir sehen, dass

${\underset{\_}{\ddot{r}}}_C(t)=\left(\begin{array}{c}0\\ {\ddot{W}}_2(t)\end{array}\right)$ und

$\delta {\underset{\_}{r}}_C=\left(\begin{array}{c}0\\ \delta W_2\end{array}\right)$.

Demnach ist

$(M_f+M_n)\cdot {\underset{\_}{\ddot{r}}}_C^T\cdot \delta {\underset{\_}{r}}_C=(M_f+M_n)\cdot {\underset{\_}{\ddot{W}}}_2\cdot \delta {\underset{\_}{W}}_2$.

Wir setzen - wie für $\delta \mathrm{\Pi }$ - die Ansatzfunktionen für $w(x,t),\delta w(x),D(x),d(x)$ in die Gleichung von $\delta W_A^a$ ein, integrieren über die Gabel-Länge und erhalten

$\delta W_A^a=-\delta \underset{\_}{\tilde{Q}}\cdot \underset{\_}{\underset{\_}{\tilde{M}}}\cdot \underset{\_}{\ddot{\tilde{Q}}}$

mit

$\underset{\_}{\underset{\_}{\tilde{M}}}=\left(\begin{array}{cccccccc}{m}_{1,11}& m_{1,12}& k_{m,13}& m_{1,14}& 0& 0& 0& 0\\ m_{1,12}& m_{1,22}& k_{m,23}& m_{1,24}& 0& 0& 0& 0\\ m_{1,13}& m_{1,23}& k_{m,33}+m_{2,11}& m_{1,34}+m_{2,12}& m_{2,13}& m_{2,14}& 0& 0\\ m_{1,14}& m_{1,24}& k_{m,34}+m_{2,12}& m_{1,44}+m_{2,22}& m_{2,23}& m_{2,24}& 0& 0\\ 0& 0& m_{2,13}& m_{2,23}& m_{2,33}+M_f+M_n& m_{2,34}& 0& 0& \\ 0& 0& m_{2,14}& m_{2,24}& m_{2,34}& m_{2,44}& 0& 0& \\ 0& 0& 0& 0& 0& 0& J_f& 0\\ 0& 0& 0& 0& 0& 0& 0& J_n\\ \end{array}\right)$

Die Koeffizienten $m_{i,jk}$ kommen dabei aus

$m_{i,jk}=\varrho {\ell }_i\underset{0}{\overset{1}{\int }}{A}_i({\xi }_i){\phi }_j({\xi }_i)\cdot {\phi }_k({\xi }_i)\text{d}{\xi }_i$

Virtuelle Arbeiten der Reibkräfte in A und B

Die äußeren, eingeprägten Reibkräfte - die Bremskraft in B und die Kontaktkraft in A - leisten virtuelle Arbeit $\delta W_F^a$ auf den virtuellen Verrückungen in A und B. Wir benätigen die Ortsvektoren zu den beiden Punkten, um daraus die Relativgeschwindigkeit der beiden Reibpartner und die variation der Ortsvektoren zu berechnen.

Wir kommen von zu Punkt A, indem wir von Punkt C aus den materiellen Punkt des Mantels auswählen, der gerade im Kontakt mit dem Asphalt ist. Wir führen dazu einen Ortsvektor C senkrect nach oben ein und rotieren diesen

• um $\mathrm{\Omega }\cdot t$ für die reguläre Drehbewegung,
• um ${\psi }_1(t)$, um einen Schlupf zwischen Mantel und Asphalt zu erfassen und
• um ${\mathrm{\Theta }}_1$, um den Kontaktpunkt in A zu erreichen.

Damit ist muss also

$\mathrm{\Omega }\cdot t+{\mathrm{\Psi }}_1(t)+{\mathrm{\Theta }}_1=\pi$

gelten. Die Koordinaten des Ortsvektors zu A auf dem Mantel schreiben wir demnach als

${\underset{\_}{r}}_A(t)={\underset{\_}{r}}_C(t)+\underset{\_}{\underset{\_}{D}}(\mathrm{\Omega }\cdot t+{\mathrm{\Psi }}_1(t)+{\mathrm{\Theta }}_1)\cdot \left(\begin{array}{c}-R_1\\ 0\end{array}\right)$

mit der Euler-Drehmatrix

$\underset{\_}{\underset{\_}{D}}(\alpha ):=\left(\begin{array}{cc}\cos (\alpha )& \sin (\alpha )\\ -\sin (\alpha )& \cos (\alpha )\end{array}\right).$

Wir finden

${\underset{\_}{r}}_A(t)=\left(\begin{array}{c}H-R_1\cos (\mathrm{\Omega }\cdot t+{\mathrm{\Psi }}_1(t)+{\mathrm{\Theta }}_1)\\ -v_0\cdot t+W_2(t)-R_1\sin (\mathrm{\Omega }\cdot t+{\mathrm{\Psi }}_1(t)+{\mathrm{\Theta }}_1)\end{array}\right)$

und

$\delta {\underset{\_}{r}}_A=\left(\begin{array}{c}0\\ \delta W_2+R_1\delta {\mathrm{\Psi }}_1\end{array}\right)$.

Und damit ist die virtuelle Arbeit der Kontaktkraft $A_H$

$\begin{array}{ccc}\delta W_{F,A}^a& =& {\left(\begin{array}{r}-A_V\\ A_H\end{array}\right)}^T\cdot \delta {\underset{\_}{r}}_A\\ & =& A_H\cdot (\delta W_2+R_1\delta {\mathrm{\Psi }}_1)\end{array}$,

wobei

$A_H=A_H({\underset{\_}{\dot{r}}}_A(t))$

Analog dazu finden wir die Koordinaten des Ortsvektors zu B - einmal ausgehend von der Gabel und einmal ausgehend vom Vorderrad.

Wir sehen, dass

$\begin{array}{l}{\underset{\_}{r}}_{B1}(t)=\left(\begin{array}{c}H-h\\ -v_0\cdot t+W_1(t)\end{array}\right)+\underset{\_}{\underset{\_}{D}}({\mathrm{\Phi }}_1(t))\cdot \left(\begin{array}{c}0\\ b\end{array}\right)\\ {\underset{\_}{r}}_{B2}(t)={\underset{\_}{r}}_C(t)+\underset{\_}{\underset{\_}{D}}(\mathrm{\Omega }\cdot t+{\mathrm{\Psi }}_2(t)+{\mathrm{\Theta }}_2)\cdot \left(\begin{array}{c}-r_2\\ 0\end{array}\right)\end{array}$

wobei

$\mathrm{\Omega }\cdot t+{\mathrm{\Psi }}_2(t)+{\mathrm{\Theta }}_2=\psi$, mit $\tan (\psi )=-\frac{b}{h}$

ist. Und die Variation der Koordinaten liefert

$\begin{array}{l}\delta {\underset{\_}{r}}_{B,1}=\left(\begin{array}{c}-b\cdot \delta {\mathrm{\Phi }}_1\\ \delta W_1\end{array}\right)\\ \delta {\underset{\_}{r}}_{B2}=\left(\begin{array}{c}{r}_2\cdot \delta {\mathrm{\Psi }}_2\sin (\psi )\\ \delta W_2-r_2\cdot \delta {\mathrm{\Psi }}_2\cos (\psi )\end{array}\right)\end{array}$

Und damit sind die virtuellen Arbeiten der Kontaktkräfte $\underset{\_}{B}$

$\delta W_{F,B}^a=-{\underset{\_}{B}}^T\cdot \delta {\underset{\_}{r}}_{B1}+{\underset{\_}{B}}^T\cdot \delta {\underset{\_}{r}}_{B2}$, wobei $\underset{\_}{B}=\underset{\_}{B}(\mathrm{\Delta }{\underset{\_}{v}}_{B,r})$,

mit der Relativgeschwindigkeit

$\begin{array}{ccc}\mathrm{\Delta }{\underset{\_}{v}}_B& =& {\underset{\_}{\dot{r}}}_{B,1}-{\underset{\_}{\dot{r}}}_{B,2}\\ & =& \left(\begin{array}{c}-r_2\sin \psi \cdot \mathrm{\Omega }-r_2\sin \psi \cdot {\dot{\mathrm{\Psi }}}_2-b\cdot \cos ({\mathrm{\Phi }}_1)\cdot {\dot{\mathrm{\Phi }}}_1\\ r_2\cos \psi \cdot \mathrm{\Omega }+r_2\cos \psi \cdot {\dot{\mathrm{\Psi }}}_2-b\cdot \sin ({\mathrm{\Phi }}_1)\cdot {\dot{\mathrm{\Phi }}}_1+{\dot{W}}_1-{\dot{W}}_2\end{array}\right)\end{array}$,

was für den stationären Fall

${\underset{\_}{v}}_{B,rel}=\left(\begin{array}{r}-r_2\sin (\psi )\mathrm{\Omega }\\ r_2\cos (\psi )\mathrm{\Omega }\end{array}\right)$

liefert. Und Achtung: hier ist $\psi$ negaiv!

Die Kraftkomponenten $B_x,B_y$ sind also in Richtung der Relativgeschwindigkeit positiv. Wir setzen deshalb mit der Reibkennlinie $\mu (v)$ und der vorzugebenden Anpresskraft $N$ an:

$\underset{\_}{B}=2\cdot \mu (||{\underset{\_}{v}}_{B,rel}||)\cdot N\cdot {\underset{\_}{e}}_B$

mit dem Einheitsvektor der Relativgeschwindigkeit

${\underset{\_}{e}}_B=\frac{{\underset{\_}{v}}_{B,rel}}{||{\underset{\_}{v}}_{B,rel}||}$

Der Faktor "2" kommt in dieser Gleichung übrigens daher, dass wir zwei Bremsbacken gegen die Bremsscheibe drücken.

Einsetzten in das Prinzip der virtuellen Verrückungen liefert nun

$\begin{array}{lll}\delta W_F^a& =& \delta W_{F,A}^a+\delta W_{F,B}^a\\ & =& -B_z\cdot \delta W_1+b\cdot B_x\cdot \delta {\mathrm{\Phi }}_1+(A_H+B_z)\dot{\delta }{W}_2+R_1\cdot A_H\cdot \delta {\mathrm{\Psi }}_1+(r_2\cdot B_x\cdot \sin (\psi )-r_2\cdot B_z\cdot \cos (\psi ))\cdot \delta {\mathrm{\Psi }}_2\\ & =& {\underset{\_}{\tilde{Q}}}^T\cdot \underset{\_}{P}\end{array}$

mit

$\underset{\_}{P}=\left(\begin{array}{c}0\\ 0\\ -B_z\\ +b{B}_x\\ A_H+B_z\\ 0\\ R_1{A}_H\\ r_2{B}_x\sin (\psi )-r_2{B}_z\cos (\psi )\\ 0\end{array}\right)$.

Zusammengefasst

Wir fassen alle Anteile zusammen und erhalten

$\begin{array}{lll}\delta W& =& \delta {\underset{\_}{\tilde{Q}}}^T\cdot (-\underset{\_}{\underset{\_}{\tilde{K}}}\cdot \underset{\_}{\tilde{Q}}-\underset{\_}{\underset{\_}{\tilde{M}}}\cdot \underset{\_}{\tilde{\ddot{Q}}}+\underset{\_}{P})\\ & \stackrel{!}{=}& 0\end{array}$

$\delta W$ kann aber nur Null sein, wenn der Klammerausdruck verschwindet, also

$\underset{\_}{\underset{\_}{\tilde{M}}}\cdot \underset{\_}{\tilde{\ddot{Q}}}+\underset{\_}{\underset{\_}{\tilde{K}}}\cdot \underset{\_}{\tilde{Q}}=\underset{\_}{P}$

Die Randbedingungen müssen wir jetzt noch einarbeiten! Das geht einfach, weil wir

$\begin{array}{ccc}{W}_0(t)& =& 0\\ \delta W_0& =& 0\\ {\mathrm{\Phi }}_0(t)& =& 0\\ \delta {\mathrm{\Phi }}_0& =& 0\end{array}$

durch Streichen der ersten beiden Zeilen des Gleichungssystems und dann jeweils der ersten beiden Spalten in den Systemmatrizen einarbeiten können.

Wir behalten also die Koordinaten $3\dots 8$, also

$\underset{\_}{Q}=\left(\begin{array}{c}{W}_1(t)\\ {\mathrm{\Phi }}_1(t)\\ W_2(t)\\ {\mathrm{\Phi }}_2(t)\\ {\mathrm{\Psi }}_1(t)\\ {\mathrm{\Psi }}_2(t)\end{array}\right)$

und die Rest-Systemmatrizen

$\underset{\_}{\underset{\_}{K}}=\left(\begin{array}{cccccc}{k}_{1,33}+k_{2,11}& k_{1,34}+k_{2,12}& k_{2,13}& k_{2,14}& 0& 0\\ k_{1,34}+k_{2,12}& k_{1,44}+k_{2,22}& k_{2,23}& k_{2,24}& 0& 0\\ k_{2,13}& k_{2,23}& k_{2,33}& k_{2,34}& 0& 0& \\ k_{2,14}& k_{2,24}& k_{2,34}& k_{2,44}& 0& 0& \\ 0& 0& 0& 0& K_S& -K_S\\ 0& 0& 0& 0& -K_S& K_S\\ \end{array}\right)$,

$\underset{\_}{\underset{\_}{M}}=\left(\begin{array}{cccccc}{k}_{m,33}+m_{2,11}& m_{1,34}+m_{2,12}& m_{2,13}& m_{2,14}& 0& 0\\ k_{m,34}+m_{2,12}& m_{1,44}+m_{2,22}& m_{2,23}& m_{2,24}& 0& 0\\ m_{2,13}& m_{2,23}& m_{2,33}+M_f+M_n& m_{2,34}& 0& 0& \\ m_{2,14}& m_{2,24}& m_{2,34}& m_{2,44}& 0& 0& \\ 0& 0& 0& 0& J_f& 0\\ 0& 0& 0& 0& 0& J_n\\ \end{array}\right)$

und

$\underset{\_}{P}=A_H\cdot \left(\begin{array}{c}0\\ 0\\ 1\\ 0\\ R_1\\ 0\\ 0\end{array}\right)+\left(\begin{array}{c}-B_z\\ b{B}_x\\ B_z\\ 0\\ 0\\ r_2{B}_x\sin (\psi )-r_2{B}_z\cos (\psi )\\ 0\end{array}\right)$.

Übrigens: in $\underset{\_}{P}$ stehen nur die Komponenten der Brems- und Kontaktkraft - und diese sind nur noch von den Zeitableitungen der Koordinaten abhängig, nachdem wir kleine Terme vernachlässigt haben. Also ist $\underset{\_}{P}=\underset{\_}{P}\left(\underset{\_}{\dot{Q}}\right)$.

Modell für die numerische Behandlung anpassen

Matlab® - und viele andere Numerik-Pakete - lösen das Anfangswertproblem mit Routinen, für die Bewegungsgleichungen als Differentialgleichungssystem erster Ordnung vorliegen müssen. Wir machen aus sechs Differentialgleichungen zweiter Ordnung deshalb 12 Differentialgleichungen erster Ordnung:

$\begin{array}{lll}\dot{\underset{\_}{Q}}& =& \underset{\_}{Q}\\ \dot{\underset{\_}{R}}& =& {\underset{\_}{\underset{\_}{M}}}^{-1}\cdot (\underset{\_}{P}-\underset{\_}{\underset{\_}{K}}\cdot \underset{\_}{Q})\end{array}$,

so dass wir mit

$\underset{\_}{Y}=\left(\begin{array}{l}\underset{\_}{Q}\\ \underset{\_}{R}\\ \end{array}\right)$

unser mathematisches Modell als

$\begin{array}{lll}\dot{\underset{\_}{Y}}& =& \left(\begin{array}{l}\underset{\_}{R}\\ {\underset{\_}{\underset{\_}{M}}}^{-1}\cdot (\underset{\_}{P}-\underset{\_}{\underset{\_}{K}}\cdot \underset{\_}{Q})\end{array}\right)\\ & =& \underset{\_}{f}\left(\underset{\_}{Y}\right)\\ \end{array}$

schreiben.

Wir führen also die Geschwindigkeiten hier explizit als Zustandsgrößen in die numerische Lösung des Probles ein.

Simulieren

Computerprogramm schreiben

stationäre Lösung berechnen

$\underset{\_}{\underset{\_}{A}}=\left(\begin{array}{ccccccc}& & & & & & 0\\ & & & & & & 0\\ & & & \underset{\_}{\underset{\_}{K}}& & & -1\\ & & & & & & 0\\ & & & & & & -R_1\\ & & & & & & 0\\ 0& 0& 1& 0& R_1& 0& 0\\ \end{array}\right)$

$\underset{\_}{b}=\left(\begin{array}{c}\\ \underset{\_}{P}\\ \\ 0\end{array}\right)$

$\underset{\_}{x}=\left(\begin{array}{c}\\ \underset{\_}{Q}\\ \\ A_H\end{array}\right)$

$\begin{array}{lll}{\underset{\_}{v}}_{A,rel}& =& {\underset{\_}{\dot{r}}}_A\\ & =& \left(\begin{array}{c}0\\ -v_0+R_1\cdot (\mathrm{\Omega }+{\dot{\mathrm{\Psi }}}_1(t))+{\dot{W}}_2(t)\end{array}\right)\end{array}$

$\dot{(.)}:=\frac{d(.)}{dt}$

$\delta W^a=\delta W_A^a+\delta W_B^a+\delta W_{d^{\prime }Alembert}^a$

$\begin{array}{lll}{\mathrm{\Theta }}_1& =& -\mathrm{\Omega }t-{\mathrm{\Psi }}_1(t)+\pi \\ {\mathrm{\Theta }}_2& =& -\mathrm{\Omega }t-{\mathrm{\Psi }}_1(t)+\psi \end{array}$

$\delta {\underset{\_}{r}}_A=\left(\begin{array}{c}0\\ \delta W_2+R_1\cdot \delta {\mathrm{\Psi }}_1\end{array}\right)$

$\delta W_A^a={\underset{\_}{F}}_A^T\cdot \delta {\underset{\_}{r}}_A$

$\delta W_B^a={\underset{\_}{F}}_B^T\cdot \delta {\underset{\_}{r}}_{B,1}-{\underset{\_}{F}}_B^T\cdot \delta {\underset{\_}{r}}_{B,2}$

${\underset{\_}{F}}_B=\left(\begin{array}{c}{B}_x\\ B_z\end{array}\right)$

${\underset{\_}{F}}_A=\left(\begin{array}{c}{A}_V\\ A_H\end{array}\right)$

$\begin{array}{lll}h& =& r_2\cos (\psi )\\ b& =-& r_2\sin (\psi )\end{array}$

Title

Text

1+1

---

Links

• ...

Literature

• ...