G-19: Sistema de masas y muelles

De MateWiki
Revisión del 20:43 9 mar 2014 de Grupo19 (Discusión | contribuciones) (Resolución de sistemas)

(dif) ← Revisión anterior | Revisión actual (dif) | Revisión siguiente → (dif)
Saltar a: navegación, buscar
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.

Sistema formado por 3 masas y 4 muelles, limitado por dos paredes verticales

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

Aproximación por el método de Euler modificado
Aproximación por el método de Runge-kutta

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
Aproximación por el método del trapecio

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
Aproximación por el método de Euler modificado

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.

Aproximación por el método del trapecio
Aproximación por el método de Euler modificado

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)
Aproximación por el método del trapecio