Sistema de dos masas y tres muelles (grupo 13)

Sistemas resorte-masa. Muchos modelos de vibraciones se modelizan con sistemas de muelles y masas.Consideramos aquí un ejemplo simple formado por dos muelles y dos masas unidos a dos paredes que se deslizan libremente sobre un plano horizontal, como se muestra en la figura.

1 Apartado Primero

Escribir el sistema de ecuaciones diferenciales para el desplazamiento de ambas masas desde la posición de equilibrio, si suponemos que no hay rozamiento. Interpretar las condiciones iniciales.

Primero designamos las variables. Tomamos [math]x_1(t)[/math] y [math]x_2(t)[/math] como la distancia de las masas a su posición de equilibrio cuando han pasado t segundos.

Como no hay fierzas de rozamiento, suponemos que no existe el amortiguamiento debido a la resistencia del aire. Por lo tanto, las únicas fuerzas que actúan en nuestro sistema serán las fuerzas restauradoras ejercidas por cada uno de los muelles. Estas fuerzas serán proporcionales al acortamiento o alargamiento de cada muelle, que depende de lo que se hayan desplazado las partículas.

Entonces, obtenemos el sistema:

[math]m_1x_1'' = - k_1x_1 + k_2(x_2 - x_1)[/math]

[math]m_2x_2'' = - k_2(x_2 - x_1) - k_3x_2[/math]

2 Apartado Segundo

Suponiendo que en el instante [math]t = 0[/math] las dos masas están desplazadas 1 y 1,5 metros hacia la derecha de la posición de equilibrio y se sueltan repentinamente, calcular la posición de cada masa en cada instante de tiempo y dibujar la gráfica. Suponer [math]k_1 = 4N/m[/math]; [math]k_2 = 2N/m[/math]; [math]k_3 = 1N/m[/math]; [math]m_1 = 2kg[/math]; [math]m_1 = 1kg[/math]; que la distancia entre las paredes es de 9 metros y que en equilibrio las dos masas están en las posiciones [math]x_A = 2m[/math] y [math]x_B = 4m[/math] respectivamente. Usar el método de Newmark para resolver sistemas de orden 2 y un método Runge-Kutta de cuarto orden, reduciendo el sistema a un sistema de primer orden. Comparar los resultados. Como las masas están desplazadas 1 y 1,5m hacia la derecha de la posición de equilibrio, [math]x_1(t_0) = 1[/math], [math]x_2(t) = 1,5[/math]. Además, como las masas se sueltan repentinamente, ambas velocidades iniciales valdrán 0.

Por otro lado nos dan los valores de m1, m2, k1, k2 y k3. Introduciendo estos valores en nuestro sistema obtenemos: [math]x_1'' = - 3x_1 + x_2[/math] [math]x_2'' = 2x_1 - 3x_2[/math]

Método Newmark:

Siendo [math]y' = z[/math], este método establece que: [math]y_n+1 = y_n + hz_n + h^2(betaf_n+1 + (1/2 - beta)f_n)[/math] [math] z_n+1 = z_n + h(gammaf_n+1 + (1 - gamma)f_n); [/math]

Tomamos beta=1/4 y gamma=1/2, porque al tomar estos valores el método es estable y de segundo orden.

Introducimos el método en Matlab:

x1(1)=1;x2(1)=1.5;z1(1)=0;z2(1)=0;
t0=0;tN=10;x=[x1 x2]';
z=[z1 z2]';
N=100; h=(tN-t0)/N;
t=t0:h:tN;
I=[1 0;0 1];
D=[-3 1;2 -3];
A=I-((h^2)/4*D);
C=I+((h^2)/4*D);
for n=1:N    
    B=C*x + h*z;
    xx=A\B;
    zz=z + ((h/2)*D)*(xx +x);
    x1(n+1)=xx(1);
    x2(n+1)=xx(2);
    z1(n+1)=zz(1);
    z2(n+1)=zz(2);
    x=xx; z=zz;
    if x1(n+1)<-2
        x1(n+1)=-2;
    end
    if x2(n+1)<-4
        x2(n+1)=-4;
    end
    if x1(n+1)>5
        x1(n+1)=5;
    end
    if x2(n+1)>7
        x2(n+1)=7;
    end
end
w1=x1+2;w2=x2+4;
hold on
plot(t,w1,'g')
plot(t,w2)
title('Posicion respecto a la pared izquierda (Apdo 2 Newmark)')
xlabel('tiempo (s)')
ylabel('Posicion')
legend('masa 1','masa 2')
hold off

Método Runge-Kutta:

Descomponemos nuestro sistema de dos ecuaciones de segundo orden en uno de cuatro ecuaciones de primer orden. Con lo que el sistema que teníamos anteriormente se transformará en: [math]y_1'=y_2[/math] [math]y_2'=-3y_1+y_3[/math] [math]y_3'=y_4[/math] [math]y_4'=2y_1+y_3[/math]

Siendo [math]y_1=x_1[/math] y [math]y_3=x_2[/math].

Almacenamos en los vectores [math]y_1[/math] e [math]y_3[/math] las posiciones de las partículas a lo largo del tiempo, y en los vectores [math]y_2[/math] e [math]y_4[/math] las velocidades. Con estos vectores formamos la matriz Y de cuatro filas en función de t.

Aplicando el método de Runge-Kutta, crearemos un bucle para obtener el valor de cada columna de Y partiendo del valor anterior. Para eso, necesitamos definir una serie de matrices [math]k_1,k_2,k_3[/math] y [math]k_4[/math].

Además, añadimos dos condiciones para que en todo momento las masas se encuentren entre las dos paredes; es decir, que la posición de la partícula 1 esté entre -2 y 7 y que la posición de la partícula 2 esté entre -4 y 5.

t0=0;tN=10;y1(1)=1;y3(1)=1.5;y2(1)=0;y4(1)=0;
N=4000;h=(tN-t0)/N;
t=t0:h:tN;
Y=[y1;y2;y3;y4];
for n=1:N
    k0=[0 1 0 0;-3 0 1 0;0 0 0 1;2 0 -3 0];
    k1=k0*Y;
    k2=k0*(Y+(1/2)*k1*h);
    k3=k0*(Y+(1/2)*k2*h);
    k4=k0*(Y+k3*h);
    Y=Y+(h/6)*(k1+2*k2+2*k3+k4);
    y1(n+1)=Y(1);
    y2(n+1)=Y(2);
    y3(n+1)=Y(3);
    y4(n+1)=Y(4);
    if y1(n+1)<-2
       y1(n+1)=-2;
    end
    if y3(n+1)<-4
       y3(n+1)=-4;
    end
    if y1(n+1)>5
       y1(n+1)=5;
    end
    if y3(n+1)>7
       y3(n+1)=7;
    end
end
x1=y1+2;
x2=y3+4;
hold on
plot(t,x1,'g')
plot(t,x2)
title('Posicion respecto a la pared izquierda (Apdo 2 Runge-Kutta)')
xlabel('tiempo (s)')
ylabel('Posicion')
legend('masa 1','masa 2')
hold off