
He creado un TikZ de un sistema de amortiguador de resorte y masa que se describe en las siguientes EDO
dx1/dt = x2
dx2/dt = -w0^2*x1 -2*zeta*w0*x2 + u
y = K*w0*x1 + 2*zeta*w0*x2
donde u
es la entrada y y
es la salida. Supongamos w0=K=1
y zeta=0.1
y las condiciones iniciales x1(0)=x2(0)=0
.
El TikZ es
\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}
Ahora, ¿cómo puedo animar este st u=sin(a*t)
con diferentes a
a lo largo de un tiempo determinado t
? y
Por supuesto, deberíamos seguir la solución de las EDO. Y la u
rueda debe "pasar" a lo largo de la sinusoide del suelo.
Respuesta1
Actualizar:Rueda "montando" en la carretera y observador en movimiento, según lo solicitado en OP y comentarios.
(Haga clic en la imagen para ejecutar el SVG interactivo [4,2 MB] en el navegador.)
El sistema ODE que lo rige se ha resuelto con \pstODEsolve
el paquete frompst-ode
.
Para que la rueda "se desplace" sobre el perfil de la carreteratu(X) a velocidad constante, una tercera ecuación diferencial, dX/dt, necesita ser solucionado. Describe el movimiento delXCoordenada del punto de contacto de la rueda con la carretera. El sistema de ecuaciones diferenciales ahora es
dX/dt= ω rueda rrueda / sqrt(1+(tu'(X)) 2 )
dymasa /dt=vmasa
dvmasa /dt= −ω 0 2 (ymasa -ycompensación -yrueda ) − 2ζω 0 (vmasa -vrueda )
El parámetro de integracióntrepresenta el tiempo, ω rueda la velocidad angular de la rueda yrrueda su radio. (Por lo tanto, rueda ω rrueda tes la distancia recorrida a lo largo del perfil de la carretera con curvas.) ω 0 es la frecuencia natural no amortiguada y ζ es la relación de amortiguación. Basado en dX/dt,X(t), radio de la ruedarrueda , elevación de la carretera localtu(X), pendientetu'(X) y segunda derivadatu''(X) las coordenadas del eje de la ruedaXrueda ,yrueda y su vert. velocidadvSe calculan las ruedas .yrueda yvLa rueda sirve ahora como entrada para el sistema masa-resorte-amortiguador.
El perfil de la carreteratu(X) se modela como una ráfaga de forma de onda:
tu(X) = porque(X−Xcompensar )mi−λ(X−Xcompensación ) 2
Observe la trayectoria no tan armónica de la rueda donde choca con los "baches", así como la curva de respuesta resultante.
La ventana móvil de la trama animada se logra configurando dinámicamente las opciones xmin
y xmax
del axis
entorno.
El solucionador de ODE (RKF45método) en el paquete pst-ode
fue codificado en PostScript. Recientemente, Marcel Krüger implementó unIntérprete de PostScripten Lua tal que ps2pdf
ya no es necesario.
Componga como PDF ejecutándolo lualatex
3 veces. Para salida GIF y SVG, descomente una de las líneas al principio de la fuente y siga las instrucciones que contiene.
%\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}