Oscilaciones acopladas de un péndulo y un muelle elástico

En la figura, se muestra los dos cuerpos y las fuerzas que actúan sobre los mismos.

Ecuaciones del movimiento

La posición del bloque de masa m es x. La posición de la masa puntual M es

xp=x+lsinθ
yp=lcosθ

Derivamos para obtener la velocidad de la masa puntual M

d x p dt = dx dt +lcosθ dθ dt d y p dt =lsinθ dθ dt

Volvemos a derivar para obtener la aceleración

d 2 x p d t 2 = d 2 x d t 2 +lcosθ d 2 θ d t 2 lsinθ ( dθ dt ) 2 d 2 y p d t 2 =lcosθ ( dθ dt ) 2 +lsinθ d 2 θ d t 2

Aplicamos la segunda ley de Newton para describir el movimiento del bloque de masa m a lo largo del eje X

m d 2 x d t 2 =Tsinθkx

El movimiento de la masa puntual M a lo largo del eje X y del eje Y es

M d 2 x p d t 2 =Tsinθ M d 2 y p d t 2 =TcosθMg

Las tres ecuaciones que describen el movimiento del sistema formado por el bloque y la masa puntual, son

m d 2 x d t 2 =Tsinθkx M( d 2 x d t 2 +lcosθ d 2 θ d t 2 lsinθ ( dθ dt ) 2 )=Tsinθ M( lcosθ ( dθ dt ) 2 +lsinθ d 2 θ d t 2 )=TcosθMg

Eliminamos la variable desconocida T, la tensión de la cuerda. Sumamos la primera y la segunda ecuación. Multiplicamos la segunda ecuación por cosθ y la tercera por sinθ y las sumamos. Obteniendo un sistema de dos ecuaciones diferenciales en x y θ.

( m+M ) d 2 x d t 2 +Ml( cosθ d 2 θ d t 2 sinθ ( dθ dt ) 2 )+kx=0 d 2 x d t 2 cosθ+l d 2 θ d t 2 +gsinθ=0

Formulación de Lagrange

Disponemos de otro modo alternativo de obtener las ecuaciones del movimiento. La energía cinética del sistema es

T= 1 2 m ( dx dt ) 2 + 1 2 M( ( d x p dt ) 2 + ( d y p dt ) 2 )= 1 2 ( m+M ) ( dx dt ) 2 +Mlcosθ dx dt dθ dt + 1 2 M l 2 ( dθ dt ) 2

La energía potencial es la suma de la energía potencial gravitatoria de la masa puntual M y la energía potencial elástica del muelle deformado x

V=Mglcosθ+ 1 2 k x 2

La Lagrangiana en términos de las coordenadas x y θ es

L=TV= 1 2 ( m+M ) x ˙ 2 +Mlcosθ( x ˙ θ ˙ +g )+ 1 2 M l 2 θ ˙ 2 1 2 k x 2

La ecuación de Euler-Lagrange para x es

t ( L x ˙ ) L x =0 L x ˙ =( m+M ) x ˙ +Mlcosθ θ ˙ t ( L x ˙ )=( m+M ) x ¨ +Ml( θ ¨ cosθ θ ˙ 2 sinθ ) L x =kx ( m+M ) d 2 x d t 2 +Ml( d 2 θ d t 2 cosθ ( dθ dt ) 2 sinθ )+kx=0

La ecuación de Euler-Lagrange para θ es

t ( L θ ˙ ) L θ =0 L θ ˙ =Ml( l θ ˙ + x ˙ cosθ ) t ( L θ ˙ )=Ml( l θ ¨ + x ¨ cosθ x ˙ θ ˙ sinθ ) L θ =Mlsinθ( x ˙ θ ˙ +g ) l d 2 θ d t 2 + d 2 x d t 2 cosθ+gsinθ=0

Llegamos por dos caminos distintos a las mismas ecuaciones

Solución numérica de las ecuaciones del movimiento

Para resolver este sistema de ecuaciones diferenciales, cada una de las dos ecuaciones solamente puede tener una única derivada segunda. Despejamos d2x/dt2 en la primera ecuación y la sustituimos en la segunda

d 2 θ d t 2 = Mlsinθcosθ ( dθ dt ) 2 kxcosθ+( m+M )gsinθ ( M sin 2 θ+m )l d 2 x d t 2 = Mlsinθ ( dθ dt ) 2 kx+Mgsinθcosθ ( M sin 2 θ+m )

Elaboramos un script con los siguientes valores de los parámetros que podremos modificar: masa del bloque m=1.5 kg, masa puntual M=0.5 kg, constante del muelle k=20 N/m, longitud de la cuerda L=0.7 m. A continuación, las condiciones iniciales: posición del bloque x=0.3, velocidad incial dx/dt=0. Desplazamiento angular inicial del péndulo θ=π/6, velocidad angular inicial dθ/dt=0.

Se resuelve el sistema de ecuaciones diferenciales mediante el procedimiento ode45, durante un tiempo de 10 s. Se representa la posición del bloque x y la posición angular del péndulo θ en función del tiempo t.

m=1.5;  %bloque
M=0.5; %masa puntual, péndulo
k=20; %constante del muelle
L=0.7; %longitud del péndulo
%x0=[0, 0, pi/6,0]; %ejemplos de condiciones iniciales
%x0=[0.1, 0.1, pi/6,pi/20];
x0=[0.3, 0, pi/4,0]; %condiciones iniciales
tspan=[0,10];
% x(1)=x, x(2)=dx/dt, x(3)=theta; x(4)=dtheta/dt
fg=@(t,x)[x(2);(M*L*sin(x(3))*x(4)^2-k*x(1)+M*9.8*sin(x(3))*cos(x(3)))
/(M*sin(x(3))^2+m); 
    x(4); -(M*sin(x(3))*cos(x(3))*x(4)^2-k*x(1)*cos(x(3))
/L+(m+M)*9.8*sin(x(3))/L)/(M*sin(x(3))^2+m)];
[t,x]=ode45(fg,tspan,x0);
hold on
plot(t,x(:,1)) 
plot(t,x(:,3))
%plot(t,x(:,2)) 
%plot(t,x(:,4))
hold off
grid on
legend('bloque','péndulo');
xlabel('t')
ylabel('x, \theta')
title('Bloque y péndulo')

Comprobamos que la energía total E=T+U es constante

E= 1 2 ( m+M ) ( dx dt ) 2 +MLcosθ( dx dt dθ dt g )+ 1 2 M L 2 ( dθ dt ) 2 + 1 2 k x 2

>> E=(m+M)*x(:,2).^2/2+M*L*cos(x(:,3)).*(x(:,2).*x(:,4)-9.8)+
M*L^2*x(:,4).^2/2+k*x(:,1).^2/2
E =-1.5254
   -1.5254
   -1.5254
   -1.5254
   -1.5254
   ...

Actividades

En primer lugar, especificamos los parámetros del sistema

En segundo lugar, especificamos las condiciones iniciales

Si queremos ver las gráficas de la posición x del bloque (en color rojo) y la posición angular θ del péndulo (en color azul), activamos la casilla titulada Gráfica

Se resuelve el sistema de dos ecuaciones diferenciales por el procedimiento de Runge-Kutta. El programa calcula la energía del sistema E en cada instante y lo compara con la energía inicial E0. Si la desviación es mayor que el 1% el programa se detiene.

| E E 0 E 0 |<0.01





Aproximaciones

Si nos restringimos a pequeños valores de los ángulos θ, las ecuaciones del movimiento se hacen mucho más simples. Aproximamos, sinθθ. Desarrollamos en serie de Taylor cosθ y tomamos los dos primeros términos no nulos del desarrollo en serie

cosθ1 1 2 θ 2

La energía cinética, potencial y la lagrangiana L=T-U con esta aproximación, se expresa

T 1 2 ( m+M ) ( dx dt ) 2 +Ml dx dt dθ dt + 1 2 M l 2 ( dθ dt ) 2 U 1 2 k x 2 Mgl( 1 1 2 θ 2 ) L= 1 2 ( m+M ) ( dx dt ) 2 +Ml dx dt dθ dt + 1 2 M l 2 ( dθ dt ) 2 1 2 k x 2 1 2 Mgl θ 2 +Mgl

Con esta aproximación, la ecuación de Euler-Lagrange para x es

t ( L x ˙ ) L x =0 ( m+M ) d 2 x d t 2 +Ml d 2 θ d t 2 +kx=0

y la ecuación de Euler-Lagrange para θ es

t ( L θ ˙ ) L θ =0 l d 2 θ d t 2 + d 2 x d t 2 +gθ=0

Despejamos las derivadas segundas de x y de θ en el sistema de las dos ecuaciones

{ m d 2 x d t 2 Mgθ+kx=0 lm d 2 θ d t 2 kx+( m+M )gθ=0

Buscamos una solución de la forma

x=Xsin(ωt+φ), θ=Θsin(ωt+φ)

que representa MAS de amplitud X, Θ y frecuencia angular ω.

( km ω 2 Mg k ( m+M )glm ω 2 )( X Θ )=( 0 0 )

Tenemos un sistema homogéneo, los cuadrados de las frecuencias de los modos normales de vibración se calculan haciendo que el determinante de los coeficientes sea igual a cero

| km ω 2 Mg k ( m+M )glm ω 2 |=0

Obtenemos una ecuación de segundo grado en ω2 con dos raíces positivas, ω1 y ω2.

m 2 l ω 4 ( lmk+m(m+M)g ) ω 2 +kmg=0

La forma general de la posición del bloque x en función del tiempo t es una combinación de los dos modos normales de vibración

x(t)=A1sin(ω1t)+ B1cos(ω1t)+ C1sin(ω2t)+ D1cos(ω2t)

Lo mismo para θ

θ(t)=A2sin(ω1t)+ B2cos(ω1t)+ C2sin(ω2t)+ D2cos(ω2t)

Las velocidades son

d x dt = A 1 ω 1 cos( ω 1 t) B 1 ω 1 sin( ω 1 t)+ C 1 ω 2 cos( ω 2 t) D 1 ω 2 sin( ω 2 t) d θ dt = A 2 ω 1 cos( ω 1 t) B 2 ω 1 sin( ω 1 t)+ C 2 ω 2 cos( ω 2 t) D 2 ω 2 sin( ω 2 t)

Relacionamos los coeficientes A1 y A2, B1 y B2, C1 y C2, D1 y D2, a través de una de las dos ecuaciones diferenciales del movimiento

m d 2 x d t 2 Mgθ+kx=0

Tenemos cuatro relaciones

{ ( km ω 1 2 ) A 1 Mg A 2 =0 ( km ω 1 2 ) B 1 Mg B 2 =0 ( km ω 2 2 ) C 1 Mg C 2 =0 ( km ω 2 2 ) D 1 Mg D 2 =0

Las condiciones iniciales determinan los valores de los cuatro coeficientes A1, B1, C1 y D1 y por tanto, de los coeficientes A2, B2, C2 y D2, a través de las relaciones establecidas anteriormente.

Desplazamos el bloque x0 y el péndulo un ángulo θ0 a continuación, los soltamos. La velocidad inicial es  dx/dt=0, dθ/dt=0 en el instante t=0.

Resolvemos el sistema de cuatro ecuaciones junto a las relaciones establecidas anteriormente.

x0=B1+D1
0=ω1·A1+ ω2·C1
θ
0=B2+D2
0=ω1·A2+ ω2·C2

La solución de este sistema es

A1=C1=A2=C2=0

B 1 = Mg θ 0 ( km ω 2 2 ) x 0 m( ω 2 2 ω 1 2 ) B 2 = km ω 1 2 Mg B 1 D 1 = Mg θ 0 +( km ω 1 2 ) x 0 m( ω 2 2 ω 1 2 ) D 2 = km ω 2 2 Mg D 1

La posición angular en función del tiempo t de cada una de las partículas es

x(t)=B1cos(ω1t)+D1cos(ω2t)
θ
(t)=B2cos(ω1t)+D2cos(ω2t)

Modos normales de vibración

Para comprobar los modos normales en la simulación, vamos a establecer los siguientes valores para los parámetros: Masa del bloque m=1.5 kg, masa puntual del péndulo M=0.5 kg, constante del muelle k=20 N/m, longitud del péndulo l=0.7 m. El desplazamiento angular inicial del péndulo es de 30°=π/6. El script nos calcula el desplazamiento inicial del bloque x0 que introducimos con los datos anteriores en el programa interactivo

m=1.5;  %bloque
M=0.5; %péndulo
k=20; %constante del muelle
L=0.7; %longitud del péndulo
%condiciones iniciales
fi0=pi/6; %desviación inicialdel péndulo

%cuadrado d elas frecuencias w_1 y w_2
a=m^2*L;
b=(m+M)*9.8*m+k*L*m;
c=m*9.8*k;
w2_1=(b+sqrt(b^2-4*a*c))/(2*a);
w2_2=(b-sqrt(b^2-4*a*c))/(2*a);

%primer modo de vibración
x0=M*9.8*fi0/(k-m*w2_1);

B1=(M*9.8*fi0-(k-m*w2_2)*x0)/(m*(w2_2-w2_1));
B2=(k-m*w2_1)*B1/(M*9.8);
D1=(-M*9.8*fi0+(k-m*w2_1)*x0)/(m*(w2_2-w2_1));
D2=(k-m*w2_2)*D1/(M*9.8);

x=@(t) B1*cos(sqrt(w2_1)*t)+D1*cos(sqrt(w2_2)*t);
theta=@(t) B2*cos(sqrt(w2_1)*t)+D2*cos(sqrt(w2_2)*t);
hold on
fplot(x,[0,10]) 
fplot(theta,[0,10])
hold off
grid on
legend('bloque','péndulo');
xlabel('t')
ylabel('x, \theta')
title('Bloque y péndulo')