Masse-Feder-Dämpfer animieren TikZ

Masse-Feder-Dämpfer animieren TikZ

Ich habe ein TikZ eines Masse-Feder-Dämpfer-Systems erstellt, das durch die folgenden Differentialgleichungen beschrieben wird

dx1/dt = x2
dx2/dt = -w0^2*x1 -2*zeta*w0*x2 + u
y = K*w0*x1 + 2*zeta*w0*x2

wobei udie Eingabe und ydie Ausgabe ist. Nehmen wir an w0=K=1, und zeta=0.1und die Anfangsbedingungen x1(0)=x2(0)=0.

Das TikZ ist

\documentclass{article}
\usepackage{tikz}
\usetikzlibrary{arrows, shapes, patterns, calc, decorations.pathmorphing, decorations.markings, positioning, animations}

\begin{document}
    \begin{tikzpicture}
        \tikzstyle{wheel} = [draw, circle, minimum size = 0.75cm]
        \tikzstyle{mass} = [draw, rectangle, minimum height = 1cm, minimum width = 2cm]
        \tikzstyle{spring} = [decorate, decoration = {zigzag, pre length = 0.3cm, post length = 0.3cm, segment length = 6}]
        \tikzstyle{damper} = [decoration = {markings, mark connection node = dmp, mark = at position 0.5 with 
        {
            \node (dmp) [inner sep = 0pt, transform shape, rotate = -90, minimum width = 5pt, minimum height = 10pt, draw=none] {};
            \draw ($(dmp.north east)+(1.5pt,0)$) -- (dmp.south east) -- (dmp.south west) -- ($(dmp.north west)+(1.5pt,0)$);
            \draw ($(dmp.north)+(0,-1.5pt)$) -- ($(dmp.north)+(0,1.5pt)$);
        }
        }, decorate]
        \tikzstyle{groundflat} = [fill, pattern = north east lines, draw = none, minimum width = 0.75cm, minimum height = 0.3cm]
        \draw[fill, pattern = north east lines] (-.5,-.25) cos (0,0) sin (.5,.25) coordinate (top) cos (1,0) sin (1.5,-.25) cos (2,0) sin (2.5,.25) cos (3,0) sin (3.5,-.25) |- (-.5,-1)--cycle;
        \node[wheel, above = 0pt of top] (u) {$u$};
        \node[mass, above = 2cm of u] (m) {$y$};
        \draw [-] (u.north) |- ++(0.5,0.25cm)coordinate (uright);
        \draw [-] (u.north) |- ++(-0.5,0.25cm)coordinate (uleft);
        \draw [spring] (uleft) -- ($(m.south) +(-0.5,0)$);
        \draw [damper] (uright) -- ($(m.south) +(0.5,0)$);
    \end{tikzpicture}
\end{document}

Wie kann ich nun dieses Rad u=sin(a*t)mit unterschiedlichen Geschwindigkeiten aüber eine bestimmte Zeit animieren t? yEs sollte natürlich der Lösung der Differentialgleichungen folgen. Und das uRad sollte entlang der Bodensinuskurve „fahren“.

Antwort1

Aktualisieren:Rad „fahren“ auf der Straße und beweglicher Beobachter, wie im OP und den Kommentaren gefordert.

(Klicken Sie auf das Bild, um das interaktive SVG auszuführen [4,2 MiB] im Browser.)

Das zugrunde liegende ODE-System wurde gelöst mit \pstODEsolvedem Paketpst-ode.

Damit das Rad auf dem Straßenprofil „fährt“Sie(X) bei konstanter Geschwindigkeit, eine dritte Differentialgleichung, dX/DT, muss gelöst werden. Es beschreibt die Bewegung derXKoordinate des Kontaktpunktes des Rades mit der Straße. Das System der Differentialgleichungen lautet nun

DX/DT= ω Rad RRad / Quadrat(1+(Sie'(X)) 2 )

DjMasse /dT=gegenMasse

DgegenMasse /dT= −ω 0 2 (jMassejVersatzjRad ) − 2ζω 0 (gegenMassegegenRad )

Der IntegrationsparameterTdie Zeit, ω Rad die Winkelgeschwindigkeit des Rades undRRad seinen Radius. (Daher ist ω Rad RRad Tist die zurückgelegte Strecke auf dem kurvenreichen Straßenprofil.) ω 0 ist die ungedämpfte Eigenfrequenz und ζ ist das Dämpfungsverhältnis. Basierend auf dX/DT,X(T), RadradiusRRad , lokale StraßenhöheSie(X), NeigungSie'(X) und 2. AbleitungSie''(X) die RadachsenkoordinatenXRad ,jRad und seine vertikale GeschwindigkeitgegenRad berechnet.jRad undgegenRad dienen nun als Eingang für das Masse-Feder-Dämpfer-System.

Das StraßenprofilSie(X) wird als Wellenformburst modelliert:

Sie(X) = cos(XXVersatz )tλ(XXVersatz ) 2

Beachten Sie die nicht ganz so harmonische Flugbahn des Rades dort, wo es in die „Schlaglöcher“ stößt, sowie die daraus resultierende Reaktionskurve.

Das bewegliche Fenster der animierten Handlung wird durch dynamisches Einstellen der xminund xmaxOptionen der axisUmgebung erreicht.

Der ODE-Solver (RKF45Methode) in pkg pst-odewurde in PostScript kodiert. Kürzlich implementierte Marcel Krüger einePostScript-Interpreterin Lua ist dies ps2pdfnicht mehr erforderlich.

Setzen Sie den Satz als PDF, indem Sie ihn lualatexdreimal ausführen. Für die GIF- und SVG-Ausgabe heben Sie die Kommentierung einer der Zeilen am Anfang der Quelle auf und folgen Sie den darin enthaltenen Anweisungen.

%\PassOptionsToPackage{export}{animate} % lualatex <file> ; convert -density 300 -alpha remove -delay 4 <file>.pdf <file>.gif
%\PassOptionsToPackage{dvisvgm,autoplay}{animate} % dvilualatex <file> ; dvisvgm --zoom=-1 --font-format=woff2 --bbox=papersize <file>.dvi
\documentclass[margin=3pt]{standalone}

\usepackage{pst-ode}
\usepackage{pgfplots}
\pgfplotsset{compat=newest}
\usepgfplotslibrary{fillbetween}
\usetikzlibrary{arrows.meta, calc, decorations}
\usepackage{listofitems} % read space separated items
\usepackage{animate}
\usepackage{xsavebox} %\xsbox

% adjustable parameters & definitions
\def\tEnd{80} % max. time (integration parameter)
\def\animationframes{201} % number of output points for animated objects
\def\outputpoints{1001} % number of output points for curves
\def\omegaWheel{1.0} % wheel angular velocity
\def\rWheel{1.0} % wheel radius
\def\rCenterOfMass{0.16em} % center-of-mass symbol size
\def\springlength{6} % length of spring in neutral state (and of damper cylinder)
\def\springturns{10} % (integer!) number of spring windings
\def\windowsizeright{5} % moving window x-size around moving object
\def\windowsizeleft{25}
\def\ymin{-1} % y axis limits
\def\ymax{17}
\pstVerb{
  tx@Dict begin
  % mass-spring-damper properties
  /w0 1.0 def
  /zeta 0.1 def
  % initial conditions
  /x0     0.0 def % x(t0)=0 starting x coordinate
  /yMass0 0.0 def % y coord of mass centre (/wo vert. offset)
  /vMass0 0.0 def % y velocity of mass centre
  % road profile U(x) with wave-form burst, and its first and second derivatives U'(x), U''(x)
  /A 1 def        % amplitude
  /a 1 def        % mode
  /lambda (1/40) AlgParser cvx exec def % burst growth/decay rate
  /xOff 20.0 def  % burst offset
  /U (A*Euler^(-lambda*(x[0]-xOff)^2)*cos(a*(x[0]-xOff))) AlgParser cvx bind def
  /dUdx (-A*Euler^(-lambda*(x[0]-xOff)^2)*a*sin(a*(x[0]-xOff))-2*lambda*(x[0]-xOff)*U) AlgParser cvx bind def
  /d2Udx2 (-(a^2+2*lambda+4*lambda^2*(x[0]-xOff)^2)*U-4*lambda*(x[0]-xOff)*dUdx) AlgParser cvx bind def
  % do not change anything below
  /rWheel \rWheel\space def
  /omegaWheel \omegaWheel\space def
  /massWheelOff (rWheel+1.5*\springlength+2)  AlgParser cvx exec def
  /dxdt (rWheel*omegaWheel/sqrt(1+dUdx^2)) AlgParser cvx bind def
  % local wheel hub coordinates
  /xWheel (x[0]-rWheel*dUdx/sqrt(1+dUdx^2)) AlgParser cvx bind def
  /yWheel (U+rWheel/sqrt(1+dUdx^2)) AlgParser cvx bind def
  /dyWheeldx (dUdx*(1-d2Udx2/(1+dUdx^2)*(yWheel-U))) AlgParser cvx bind def
  % angular position of wheel in terms of t
  /thetaWheel {omegaWheel t mul 2 Pi mul div dup cvi sub neg 360 mul} bind def
  end
}

% solve equations of motion with many output points for smooth curves
\pstODEsolve[algebraicAll,saveData]{curves}{% PS variable and file that take result table
                           % table format to be saved in `curves' and `curves.dat':
  x[0] |                   %    x coordinate of contact point
  U |                      %    y coordinate U(x) of contact point
  xWheel | yWheel |        %    wheel hub x (same as that of mass) and y coordinates
  x[1]+massWheelOff+rWheel %    y coordinate of mass (response with const. vertical offset)
}{0}{\tEnd}{\outputpoints}{% t_0, t_end, number of t steps plus 1
  x0 | yMass0 | vMass0     % initial conditions
}{                         % ODE system's RHS:
  dxdt |                   %   dx/dt, movement of contact point (x coordinate)
  x[2] |                   %   centre of mass vert. velocity
  -w0^2*(x[1]-yWheel+rWheel)%  centre of mass vert. acceleration
    -2*zeta*w0*(x[2]-dyWheeldx*dxdt)
}

% solve equations with fewer output points for mass-spring-damper animation; save solution
% table (wheel and mass centre coordinates, valve angular position) in `objects' and `objects.dat'
\pstODEsolve[algebraicAll,saveData]{objects}{
  xWheel | yWheel | x[1]+massWheelOff+rWheel | thetaWheel-90
}{0}{\tEnd}{\animationframes}{x0 | yMass0 | vMass0}{
  dxdt | x[2] | -w0^2*(x[1]-yWheel+rWheel)-2*zeta*w0*(x[2]-dyWheeldx*dxdt)
}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% \fileopenr{<file stream>}{<file name>}, opens file for reading
\newcommand\fileopenr[2]{%
  \newread#1%
  \immediate\openin#1=#2%
}
% \readtolist[<sep char>]{<file stream>}{\list}
% reads a line from file stream and splits at <sep char> into \list[1], \list[2], ...
\newcommand\readtolist[3][,]{{%
  \setsepchar{#1}%
  \immediate\read#2 to \inputline%
    \ifeof#2
      \immediate\closein#2%
      \ifdefined\multiframebreak\multiframebreak\fi%
    \else%
      \greadlist*#3\inputline%
    \fi%
}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\tikzset{spring/.style 2 args= {decorate, decoration = {zigzag, segment length = \dimexpr(#1-1pt)/\springturns\relax, amplitude=#2}}}%

\begin{document}
\IfFileExists{curves.dat}{}{dummy text\end{document}}%
\xsbox{centreofmass}{%
  \tikz[radius=\rCenterOfMass] {% from https://tex.stackexchange.com/questions/23453
    \fill (0,0) -- ++(\rCenterOfMass,0) arc [start angle=0,end angle=90] -- ++(0,{-2*\rCenterOfMass}) arc [start angle=270, end angle=180];%
    \draw [thin] (0,0) circle;%
  }%
}%
\begin{animateinline}[controls]{25}
  \fileopenr{\data}{objects.dat}
  \readtolist[ ]{\data}{\table}
  \multiframe{10000}{}{
    \begin{tikzpicture}[line cap=rect, line join=bevel]
      \begin{axis}[
        unit vector ratio = 1, axis on top,
        xmin={\table[1]-\windowsizeleft}, xmax={\table[1]+\windowsizeright},
        ymin=\ymin, ymax=\ymax,
        xtick distance=5,
        ytick distance=2,
        xlabel={$x$}, ylabel={$y$}, ylabel style={rotate=-90},
        extra x tick labels={\phantom{\strut00}},extra x ticks={\pgfkeysvalueof{/pgfplots/xmax}}, % prevents wobbling of axis
      ]
        % road profile; off-window coordinates filtered away for smaller file size
        \addplot [name path=road,no markers,unbounded coords=jump,x filter/.expression={x<\pgfkeysvalueof{/pgfplots/xmin}-0.1 || x>\pgfkeysvalueof{/pgfplots/xmax}+0.1 ? nan : x}]
          table [x index=0, y index=1] {curves.dat} -- (\pgfkeysvalueof{/pgfplots/xmax},0);
        \path [name path=bottom] (0,\pgfkeysvalueof{/pgfplots/ymin}) -- (\pgfkeysvalueof{/pgfplots/xmax},\pgfkeysvalueof{/pgfplots/ymin});
        \addplot [fill=black!20!white] fill between[of=road and bottom];
        % wheel hub trajectory
        \addplot [no markers,blue,unbounded coords=jump,x filter/.expression={x<\pgfkeysvalueof{/pgfplots/xmin}-0.1 || x>\table[1] ? nan : x}]
          table [x index=2, y index=3] {curves.dat};
        % centre of mass trajectory
        \addplot [no markers,red,unbounded coords=jump,x filter/.expression={x<\pgfkeysvalueof{/pgfplots/xmin}-0.1 || x>\table[1] ? nan : x}]
          table [x index=2, y index=4] {curves.dat};
        \coordinate (w) at (\table[1],\table[2]);  % wheel hub
        \coordinate (m) at (\table[1],\table[3]);  % centre of mass
        % wheel with valve
        \fill[even odd rule] (w) circle [radius={0.7*\rWheel}] (w) circle [radius=\rWheel];
        \draw[thick] ($(w)+(axis direction cs:{0.6*\rWheel*cos(\table[4])},{0.6*\rWheel*sin(\table[4])})$)
          -- ++(axis direction cs:{0.1*\rWheel*cos(\table[4])},{0.1*\rWheel*sin(\table[4])});
        % mass
        \node at (m) {\thecentreofmass};
        \draw (m) ++(axis direction cs:-1,-1) rectangle ++(axis direction cs:2,2);
        % spring/damper support
        \draw (m) ++(axis direction cs:-1,-1) -- ++(axis direction cs:0,{-0.25-0.25*\springlength}) coordinate (sprtop);
        \path (m) ++(axis direction cs:1,-1) coordinate (dmptop);
        \draw (w) ++(axis direction cs:0,{\rWheel+0.5}) -| ++(axis direction cs:-1,{0.25+0.25*\springlength}) coordinate (sprbot);
        \draw let \p1=(axis direction cs:0.35,0) in [{Circle[open, length=\x1]}-, shorten <=-0.5*\x1] (w) -- ++(axis direction cs:0,{\rWheel+0.5}) -| ++(axis direction cs:1,0.5) coordinate (dmpbot);
        % spring
        \path let \p1=($(sprtop)-(sprbot)$), \p2=(axis direction cs:0.6,0) % length and amplitude (half width)
              in [draw,spring={\y1}{\x2}] (sprbot) -- (sprtop);
        % damper lower part
        \draw (dmpbot) -| +(axis direction cs:0.4,\springlength) (dmpbot) -| +(axis direction cs:-0.4,\springlength);
        % damper upper part
        \draw (dmptop) -- ++(axis direction cs:0,-\springlength) ++(axis direction cs:-0.3,0) -- ++(axis direction cs:0.6,0);
      \end{axis}
    \end{tikzpicture}%
    \readtolist[ ]{\data}{\table}
  }
\end{animateinline}
\end{document}

verwandte Informationen