Sistema de dos masas y tres muelles (grupo 13)

De MateWiki
Saltar a: navegación, buscar
Trabajo realizado por estudiantes
Título Sistemas resorte-masa
Asignatura Ecuaciones Diferenciales
Curso 2012-13
Autores {{{4}}}
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura


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.

Figuradosmasas.png

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:

Masa 1: [math]m_1x_1'' = - k_1x_1 + k_2(x_2 - x_1)[/math] Masa 2: [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:

 1 x1(1)=1;x2(1)=1.5;z1(1)=0;z2(1)=0;
 2 t0=0;tN=10;x=[x1 x2]';
 3 z=[z1 z2]';
 4 N=100; h=(tN-t0)/N;
 5 t=t0:h:tN;
 6 I=[1 0;0 1];
 7 D=[-3 1;2 -3];
 8 A=I-((h^2)/4*D);
 9 C=I+((h^2)/4*D);
10 for n=1:N    
11     B=C*x + h*z;
12     xx=A\B;
13     zz=z + ((h/2)*D)*(xx +x);
14     x1(n+1)=xx(1);
15     x2(n+1)=xx(2);
16     z1(n+1)=zz(1);
17     z2(n+1)=zz(2);
18     x=xx; z=zz;
19     if x1(n+1)<-2
20         x1(n+1)=-2;
21     end
22     if x2(n+1)<-4
23         x2(n+1)=-4;
24     end
25     if x1(n+1)>5
26         x1(n+1)=5;
27     end
28     if x2(n+1)>7
29         x2(n+1)=7;
30     end
31 end
32 w1=x1+2;w2=x2+4;
33 hold on
34 plot(t,w1,'g')
35 plot(t,w2)
36 title('Posicion respecto a la pared izquierda (Apdo 2 Newmark)')
37 xlabel('tiempo (s)')
38 ylabel('Posicion')
39 legend('masa 1','masa 2')
40 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.

 1 t0=0;tN=10;y1(1)=1;y3(1)=1.5;y2(1)=0;y4(1)=0;
 2 N=4000;h=(tN-t0)/N;
 3 t=t0:h:tN;
 4 Y=[y1;y2;y3;y4];
 5 for n=1:N
 6     k0=[0 1 0 0;-3 0 1 0;0 0 0 1;2 0 -3 0];
 7     k1=k0*Y;
 8     k2=k0*(Y+(1/2)*k1*h);
 9     k3=k0*(Y+(1/2)*k2*h);
10     k4=k0*(Y+k3*h);
11     Y=Y+(h/6)*(k1+2*k2+2*k3+k4);
12     y1(n+1)=Y(1);
13     y2(n+1)=Y(2);
14     y3(n+1)=Y(3);
15     y4(n+1)=Y(4);
16     if y1(n+1)<-2
17        y1(n+1)=-2;
18     end
19     if y3(n+1)<-4
20        y3(n+1)=-4;
21     end
22     if y1(n+1)>5
23        y1(n+1)=5;
24     end
25     if y3(n+1)>7
26        y3(n+1)=7;
27     end
28 end
29 x1=y1+2;
30 x2=y3+4;
31 hold on
32 plot(t,x1,'g')
33 plot(t,x2)
34 title('Posicion respecto a la pared izquierda (Apdo 2 Runge-Kutta)')
35 xlabel('tiempo (s)')
36 ylabel('Posicion')
37 legend('masa 1','masa 2')
38 hold off


Apdo2newmarkgraph.1.jpg

Apdo2rungegraph.1.jpg

3 Apartado Tercero

Realizar gráficas de las posiciones de las masas a lo largo del tiempo en los casos particulares siguientes: A) El sistema parte del equilibrio con velocidades iniciales de 1m/s para las dos masas, en el mismo sentido. B) El sistema parte del equilibrio con velocidades iniciales de 1m/s para las dos masas, pero en sentidos opuestos. Dibujar en una gráfica las posiciones de las masas, es decir los puntos [math](x_A(t); x_B(t))[/math] para valores de t entre 0 y 10.


Tendremos que definir nuevamente las condiciones iniciales. Por un lado, las masas parten de sus posiciones de equilibrio. Por otro lado, las velocidades serán de 1m/s en el mismo sentido para ambas masas.

Método Newmark para mismo sentido.

Utilizamos el mismo método del apartado dos variando las condiciones iniciales.

 1 t0=0;tN=10;x1(1)=0;x2(1)=0;z1(1)=1;z2(1)=1;
 2 x=[x1 x2]';z=[z1 z2]';
 3 N=400; h=(tN-t0)/N;
 4 t=t0:h:tN;
 5 I=[1 0;0 1];
 6 D=[-3 1;2 -3];
 7 A=I-((h^2)/4*D);
 8 C=I+((h^2)/4*D);
 9 for n=1:N
10     B=C*x + h*z;
11     xx=A\B;
12     zz=z + ((h/2)*D)*(xx +x);
13     x1(n+1)=xx(1);
14     x2(n+1)=xx(2);
15     z1(n+1)=zz(1);
16     z2(n+1)=zz(2);
17     x=xx; z=zz;
18 end


Método Runge-Kutta para mismo sentido.

Utilizamos el mismo método del apartado dos variando las condiciones iniciales.

 1 t0=0;tN=10;
 2 N=4000;h=(tN-t0)/N;
 3 t=t0:h:tN;
 4 y1(1)=0;y2(1)=1;y3(1)=0;y4(1)=1;
 5 Y=[y1;y2;y3;y4];
 6 for n=1:N
 7     k0=[0 1 0 0;-3 0 1 0;0 0 0 1;2 0 -3 0];
 8     k1=k0*Y;
 9     k2=k0*(Y+(1/2)*k1*h);
10     k3=k0*(Y+(1/2)*k2*h);
11     k4=k0*(Y+k3*h);
12     Y=Y+(h/6)*(k1+2*k2+2*k3+k4);
13     y1(n+1)=Y(1);
14     y2(n+1)=Y(2);
15     y3(n+1)=Y(3);
16     y4(n+1)=Y(4);
17 end

Apto3anewmark.1.jpg Apto3arungek.1.jpg

Método Newmark y Runge-Kutta para sentidos opuestos.

Se realizan igual que en los apartados anteriores, cambiando la velocidad inicial de la segunda masa por -1m/s.

Apto3bnewmark.1.jpg Apto3brungek.1.jpg

4 Apartado cuarto

Supondremos ahora que el sistema tiene sólo una masa. Es decir, desenganchamos el muelle de la segunda masa. Supondremos también que el sistema está sumergido en un medio viscoso que provoca un amortiguamiento en el comportamiento del muelle proporcional a la velocidad de la masa, con coeficiente [math]μ = 1[/math]. Escribir el nuevo sistema de ecuaciones y resolverlo para los parámetros anteriores. Usar tanto el método de Runge-Kutta de cuarto orden como el de Newmark. Dibujar en una gráfica la energía a lo largo del tiempo, en escala logarítmica. Es cierto que al incrementar μ aumenta la tasa de decaimiento de la energía, es decir, decrece más rápidamente? Explicar si esto es razonable.

Llamamos [math]x(t)[/math] a la distancia de la masa a su posición de equilibrio cuando han pasado t segundos.

Esta vez, sí existe fuerza de amortiguamiento debida a la resistencia del medio, que será proporcional a la velocidad. Además, actuará la fuerza restauradora del muelle, que será proporcional al acortamiento o alargamiento del muelle (distancia de la masa a la posición de equilibrio).

Entonces, obtendremos la ecuación de segundo orden:

[math]mx''=-kx-mux'[/math]

Donde, al sustituir los parámetros, obtenemos:

[math]2x''=-4x-mux'[/math] Es decir: [math]x''=-2x-(mu/2)x'[/math]

Método de Newmark

 1 t0=0;tN=20;
 2 y(1)=input('posicion inicial de la masa (mayor que -2)')
 3 z(1)=input('velocidad inicial de la masa')
 4 %beta=1/4 gamma=1/2;
 5 N=4000; h=(tN-t0)/N; 
 6 t=t0:h:tN;
 7 for n=1:N
 8     y(n+1)=((1-((h^3)/(8+2*h)))*y(n)+(h-((h^2)/(8+2*h))-((h^3)/(4*(8+2*h))))*z(n))/(1-((h^2)/2)-((h^3)/(8+2*h)));
 9     z(n+1)=(z(n)-h*y(n+1)-h*y(n)-(h/4)*z(n))/(1+(h/4));   
10 end
11 E=2*(y.^2)+(z.^2);
12 semilogx(t,E,'g')
13 title('Variacion energia sistema (Newmark)')
14 xlabel('tiempo (s)')
15 ylabel('Energia(J)')


Método de Runge-Kutta

 1 t0=0;tN=20;
 2 y1(1)=input('posicion inicial de la masa (mayor que -2)')
 3 y2(1)=input('velocidad inicial de la masa')
 4 N=4000;h=(tN-t0)/N;
 5 t=t0:h:tN;mu=1;
 6 Y=[y1;y2];
 7 for n=1:N
 8     k0=[0 1;-2 -2*(1/2)*mu];
 9     k1=k0*Y;
10     k2=k0*(Y+(1/2)*k1*h);
11     k3=k0*(Y+(1/2)*k2*h);
12     k4=k0*(Y+k3*h);
13     Y=Y+(h/6)*(k1+2*k2+2*k3+k4);
14     y1(n+1)=Y(1);
15     y2(n+1)=Y(2);
16 end
17 E=2*(y1.^2)+(y2.^2);
18 semilogx(t,E,'g')
19 title('Variacion energia sistema (Runge-Kutta)')
20 xlabel('tiempo (s)')
21 ylabel('Energia(J)')

Apdo4energianewmark.1.jpg Apdo4energiarunge.1.jpg

5 La energía en función de μ

Apdo4mu10.1.jpg Para μ=10

Apdo4mu100.1.jpg Para μ=100

Observamos que para valores mayores de μ hay una mayor pérdida inicial de energía, pero tarda más tiempo en anularse completamente.