G-19: Sistema de masas y muelles
| Trabajo realizado por estudiantes | |
|---|---|
| Título | Modelo masa-resorte. Grupo 19-B |
| Asignatura | Ecuaciones Diferenciales |
| Curso | Curso 2013-14 |
| Autores | Javier Abad, Jesús Castaño, Ignacio Embid, Ángela Pozo, Javier Pérez, Cristino Pérez |
| Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura | |
En este trabajo vamos estudiar el comportamiento del sistema formado por 3 masas y 4 muelles de constantes k1, k2, k3 y k4, limitados por paredes verticales en ambos extremos.
Contenido
1 Ecuaciones del movimiento
Se van a considerar 3 grados de libertad x1, x2, x3, que van a ser los desplazamientos de las masas respecto de sus posiciones de equilibrio. Estableciendo [math]F=m*x''[/math] y aplicando la ley de Hooke obtenemos:
m1: [math]m_1x_1'' = - k_1x_1 + k_2(x_2 - x_1)[/math] m2: [math]m_2x_2'' = - k_2(x_2 - x_1) - k_3(x_2 - x_3)[/math] m3: [math]m_3x_3’’= -k_3(x_3-x_2)- k_4x_3[/math]
2 Resolución de sistemas
Si en el instante t= 0 las tres masas están desplazadas 0.5, 1 y 0.8 metros hacia la derecha de la posición de equilibrio y se sueltan repentinamente, sin velocidad, vamos a calcular la posición x1(t), x2(t), x3(t), con respecto a su estado de equilibrio. Teniendo en cuanta los siguiente datos [math]k_1 = 4N/m[/math]; [math]k_2 = 2N/m[/math]; [math]k_3 = 1N/m[/math]; [math]k_4 = 3N/m[/math]; [math]m_1 = 2kg[/math]; [math]m_2 = 1kg[/math];[math]m_3 = 3kg[/math]; que la distancia entre las paredes es de 12 metros y que en equilibrio las tres masas están en las posiciones [math]x_1 = 2.5m[/math] ; [math]x_2 = 4m[/math]; y [math]x_3 = 8m[/math]. Primero procedemos a reducir el orden del sistema, descomponiendolo en 6 ecuaciones de primer orden:
[math]x_1'=x_a[/math]: [math]x_2'=x_b[/math]: [math]x_3'=x_c[/math]: [math]x_a'=x_1''[/math]: [math]x_b'=x_2''[/math]: [math]x_c'=x_3''[/math]:
Resolución por el método de Euler modificado:
%euler modificado
%Datos
t0=0;
tf=10;
x0=[0,0,0,0.5,1,0.8]';
k1=4;
k2=2;
k3=1;
k4=3;
m1=2;
m2=1;
m3=3;
%Discretización
h=0.1;
N=(tf-t0)/h;
%vector tiempo y matriz posicion
t=t0:h:tf;
x=zeros(6,N+1);
%Inicializamos
x(:,1)=x0;
xx=x0;
% Nota: la funcion no depende del tiempo luego el tiempo no esta presente
%Interacciones
for n=1:N
K1=[-k1/m1*xx(4)+k2/m1*(xx(5)-xx(4));-k2/m2*(xx(5)-xx(4))+k3/m2*(xx(6)-xx(5));-k3/m3*(xx(6)-xx(5))-k4/m3*xx(6);xx(1);xx(2);xx(3)];
xp=(xx+h*K1/2);
K2=[-k1/m1*xp(4)+k2/m1*(xp(5)-xp(4));-k2/m2*(xp(5)-xp(4))+k3/m2*(xp(6)-xp(5));-k3/m3*(xp(6)-xp(5))-k4/m3*xp(6);xp(1);xp(2);xp(3)];
xx=xx+h/2*(K1+K2);
x(:,n+1)=xx;
end
plot(t,x);
Resolución por el método de Runge-kutta:
%runhekuter4
%Datos
t0=0;
tf=10;
x0=[0,0,0,0.5,1,0.8]';
k1=4;
k2=2;
k3=1;
k4=3;
m1=2;
m2=1;
m3=3;
%Discretización
h=0.025;
N=(tf-t0)/h;
%vector tiempo y matriz posicion
t=t0:h:tf;
x=zeros(6,N+1);
%Inicializamos
x(:,1)=x0;
xx=x0;
% Nota: la funcion no depende del tiempo luego el tiempo no esta presente
%Interacciones
for n=1:N
K1=[-k1/m1*xx(4)+k2/m1*(xx(5)-xx(4));-k2/m2*(xx(5)-xx(4))+k3/m2*(xx(6)-xx(5));-k3/m3*(xx(6)-xx(5))-k4/m3*xx(6);xx(1);xx(2);xx(3)];
xp=(xx+h*K1/2);
K2=[-k1/m1*xp(4)+k2/m1*(xp(5)-xp(4));-k2/m2*(xp(5)-xp(4))+k3/m2*(xp(6)-xp(5));-k3/m3*(xp(6)-xp(5))-k4/m3*xp(6);xp(1);xp(2);xp(3)];
xp=(xp+h*K2/2);
K3=[-k1/m1*xp(4)+k2/m1*(xp(5)-xp(4));-k2/m2*(xp(5)-xp(4))+k3/m2*(xp(6)-xp(5));-k3/m3*(xp(6)-xp(5))-k4/m3*xp(6);xp(1);xp(2);xp(3)];
xp=(xp+h*K3/2);
K4=[-k1/m1*xp(4)+k2/m1*(xp(5)-xp(4));-k2/m2*(xp(5)-xp(4))+k3/m2*(xp(6)-xp(5));-k3/m3*(xp(6)-xp(5))-k4/m3*xp(6);xp(1);xp(2);xp(3)];
xx=xx+h/6*(K1+2*K2+2*K3+K4);
x(:,n+1)=xx;
end
plot(t,x);
3 Interpretación Gráfica
4 Variaciones en la velocidad y la posición iniciales
4.1 Caso a)
Se desplazan las masas, 0.5m a la izquierda la primera, 1m a la derecha la segunda y 0.5m a la izquierda la tercera. Las tres masas se sueltan con una velocidad inicial de 1m/s hacia la derecha. Por el método del trapecio:
[math] y_0, \; t_0 [/math]
[math]y_{n+1} = y_{n} + \frac{h}{2}\left[ f(t_{n},y_{n}) + f(t_{n+1},y_{n+1}) \right][/math]
La ecuación a resolver sería
[math]y_{n+1} = y_{n} + \frac{h}{2} \left[ y_{n}\cdot(1 - y_{n}) + y_{n+1}\cdot(1 - y_{n+1}) \right][/math]
%trapecio
%Datos
t0=0;
tf=10;
x0=[1,1,1,-0.5,1,-0.5]';
k1=4;
k2=2;
k3=1;
k4=3;
m1=2;
m2=1;
m3=3;
%DiscretizaciÛn
h=0.2;
N=(tf-t0)/h;
%vector tiempo y matriz posicion
t=t0:h:tf;
x=zeros(6,N+1);
%A=matriz coeficientes
A=[0,0,0,-(k1+k2)/m1,k2/m1,0;0,0,0,k2/m2,-(k2+k3)/m2,k3/m2;0,0,0,0,k3/m3,-(k3+k4)/m3;1,0,0,0,0,0;0,1,0,0,0,0;0,0,1,0,0,0];
%Inicializamos
x(:,1)=x0;
%Interacciones
for n=1:N
x(:,n+1)=inv(eye(size(A))-h/2*A)*(x(:,n)+h/2*A*x(:,n));
end
hold on
plot(t,x(4,:));
plot(t,x(5,:));
plot(t,x(6,:));
hold off
Posiciones relativas de las masas por el método de euler modificado.
%euler modificado
%Datos
t0=0;
tf=10;
x0=[1,1,1,-0.5,1,-0.5]';
k1=4;
k2=2;
k3=1;
k4=3;
m1=2;
m2=1;
m3=3;
%DiscretizaciÛn
h=0.05;
N=(tf-t0)/h;
%vector tiempo y matriz posicion
t=t0:h:tf;
x=zeros(6,N+1);
%Inicializamos
x(:,1)=x0;
xx=x0;
% Nota: la funcion no depende del tiempo luego el tiempo no esta presente
%Interacciones
for n=1:N
K1=[-k1/m1*xx(4)+k2/m1*(xx(5)-xx(4));-k2/m2*(xx(5)-xx(4))+k3/m2*(xx(6)-xx(5));-k3/m3*(xx(6)-xx(5))-k4/m3*xx(6);xx(1);xx(2);xx(3)];
xp=(xx+h*K1/2);
K2=[-k1/m1*xp(4)+k2/m1*(xp(5)-xp(4));-k2/m2*(xp(5)-xp(4))+k3/m2*(xp(6)-xp(5));-k3/m3*(xp(6)-xp(5))-k4/m3*xp(6);xp(1);xp(2);xp(3)];
xx=xx+h/2*(K1+K2);
x(:,n+1)=xx;
end
subplot(2,1,1)
hold on
plot(t,x(4,:));
plot(t,x(6,:));
hold off
subplot(2,1,2)
hold on
plot(t,x(4,:));
plot(t,x(5,:));
plot(t,x(6,:));
hold off
4.2 Caso b)
En este caso se desplazan las masas, 0.5m a la derecha la primera, 1m a la derecha la segunda y 0.5m a la izquierda la tercera. Las primera masa se suelta con una velocidad inicial de 1m/s hacia la izquierda, la segunda sin velocidad inicial y la tercera con una velocidad de 0.5 m/s hacia la derecha.
En este caso el código de matlab será el mismo que en el apartado anterior, únicamente tendremos que modificar las condiciones iniciales.
5 Sistema de una masa con una fuerza exterior y que se desplaza en un medio viscoso
Vamos a trabajar ahora con una sistema formado por una sola masa y un solo muelle, dicho sistema esta sumergido en un medio viscoso que produce un amortiguamiento en el comportamiento del muelle y sobre ella actúa una fuerza exterior en la dirección del movimiento, con lo que la ecuación queda:
[math] mx''=-kx-cx'+2\ltsup\gt(-0.01t)\lt/sup\gtsen(2t) [/math]
%trapecio sistema no homogeneo
%Datos
t0=0;
tf=60;
x0=[0.3,0.3]';
k=4;
m=1;
nu=1;
%DiscretizaciÛn
h=0.2;
N=(tf-t0)/h;
%vector tiempo y matriz posicion
t=t0:h:tf;
x=zeros(2,N+1);
f=[];
for i=1:(length(t))
f=[f;2*exp(-0.01*t(i))*sin(2*t(i))];
end
%A=matriz coeficientes
A=[-nu/m,-k/m;1,0];
%Inicializamos
x(:,1)=x0;
%Interacciones
for n=1:N
x(:,n+1)=inv(eye(size(A))-h/2*A)*(x(:,n)+h/2*(A*x(:,n)+f(n)+f(n+1)));
end
subplot(2,1,1)
plot(t,x(2,:))
E=k/2*(x(2,:).^2)+m/2*(x(1,:).^2);
subplot(2,1,2)
semilogx(t,E)
