Application Center - Maplesoft

App Preview:

Péndulo Doble

You can switch back to the summary page by clicking here.

Learn about Maple
Download Application


 

Image 

 

P?ndulo Doble. 

 

 

Sebasti?n Varas Kittel, 

Ingenier?a Civil Industrial, 

Universidad Adolfo Ib??ez. 

 

sebastian.varask2006@alumnos.uai.cl 

 

 

Descripci?n: 

 

A continuaci?n analizaremos el comportamiento del p?ndulo doble. Utilizando el Principio de M?nima Acci?n derivaremos las ecuaciones del movimiento, las cuales linealizaremos entorno a sus puntos de equilibro, para as? poder analizar el comportamiento del p?ndulo doble, para peque?as oscilaciones, a trav?s de gr?ficos y animaciones que simular?n el movimiento de las masas. 

 

Image 

Package de Maple: 

 

restart; -1 

with(plots); -1 

with(LinearAlgebra); -1 

with(VariationalCalculus); -1 

 

Coordenadas: 

 

Definiremos las coordenadas de la posici?n de la masas del p?ndulo doble en funci?n de los ?ngulos (Coordenadas Generalizadas) θ1(t) y θ2(t).  

 

Masa 1: 

`:=`(X[1], proc (t) options operator, arrow; `*`(l[1], `*`(sin(theta[1](t)))) end proc); -1 

`:=`(Y[1], proc (t) options operator, arrow; `+`(`-`(`*`(l[1], `*`(cos(theta[1](t)))))) end proc); -1 

 

Masa 2: 

`:=`(X[2], proc (t) options operator, arrow; `+`(X[1](t), `*`(l[2], `*`(sin(theta[2](t))))) end proc); -1 

`:=`(Y[2], proc (t) options operator, arrow; `+`(Y[1](t), `-`(`*`(l[2], `*`(cos(theta[2](t)))))) end proc); -1 

 

 

Energ?a: 

 

Energ?a Cin?tica: 

 

Masa 1: 

 

 

`:=`(EKM1, proc (t) options operator, arrow; `+`(`*`(`/`(1, 2), `*`(m[1], `*`(`+`(`*`(`^`(diff(X[1](t), t), 2)), `*`(`^`(diff(Y[1](t), t), 2))))))) end proc) 

 

Masa 2: 

 

`:=`(EKM2, proc (t) options operator, arrow; `+`(`*`(`/`(1, 2), `*`(m[2], `*`(`+`(`*`(`^`(diff(X[2](t), t), 2)), `*`(`^`(diff(Y[2](t), t), 2))))))) end proc) 

 

Sistema: 

 

`:=`(E[k], proc (t) options operator, arrow; `+`(EKM1(t), EKM2(t)) end proc); -1 

 

 

Energ?a Potencial: 

 

Masa 1: 

`:=`(EPM1, proc (t) options operator, arrow; `*`(m[1], `*`(g, `*`(Y[1](t)))) end proc); -1 

 

Masa 2: 

`:=`(EPM2, proc (t) options operator, arrow; `*`(m[2], `*`(g, `*`(Y[2](t)))) end proc); -1 

 

Sistema: 

`:=`(E[p], proc (t) options operator, arrow; `+`(EPM1(t), EPM2(t)) end proc); -1 

 

Principio de M?nima Acci?n: 

 

Debemos minimizar el funcional que se forma al hacer la diferencia entre la Energ?a Cin?tica Total del Sistema y la Energ?a Potencial Total del Sistema. Para esto recurriremos a las ecuaciones de Euler-Lagrange. 

Image 

 

Funcional: 

`:=`(L, proc (t) options operator, arrow; `+`(E[k](t), `-`(E[p](t))) end proc); -1 

 

L(t) 

 

Euler-Lagrange: 

 

θ1(t): 

`:=`(EDO[1], simplify(EulerLagrange(L(t), t, theta[1](t))[1]) = 0) 

θ2(t): 

`:=`(EDO[2], simplify(EulerLagrange(L(t), t, theta[2](t))[1]) = 0) 

 

Puntos Cr?ticos: 

 

Podemos ver que las ecuaciones obtenidas son altamente no lineales, por lo que es necesario linealizarla entorno a los puntos cr?ticos estables del sistema. 

 

Para encontrar puntos cr?ticos, es necesario encontrar los valores de θ1(t) y θ2(t) para los que se cumplen las siguientes igualdades: 

Image  

solve(eval({EDO[1], EDO[2]}, [theta[1](t) = theta[1], theta[2](t) = theta[2], (D[1, 1](theta[1]))(t) = 0, (D[1, 1](theta[2]))(t) = 0, (D[1](theta[1]))(t) = 0, (D[1](theta[2]))(t) = 0]), {theta[1], the...
solve(eval({EDO[1], EDO[2]}, [theta[1](t) = theta[1], theta[2](t) = theta[2], (D[1, 1](theta[1]))(t) = 0, (D[1, 1](theta[2]))(t) = 0, (D[1](theta[1]))(t) = 0, (D[1](theta[2]))(t) = 0]), {theta[1], the...
 

 

 

Por lo que tenemos que el punto cr?tico estable se obtiene cuando: 

Image 

 

Linealizaci?n: 

 

La no linealidad es producida por las siguientes funciones de funciones del tiempo:
Image

Si linealizamos estas funciones a trav?s de una expansi?n de taylor de grado 1 centrada en los puntos cr?ticos, tendremos que:
Image
 

 

De esta forma, las ecuaciones del movimiento de las masas del p?ndulo doble, linealizadas entorno al punto de equilibrio estable, son: 

 

`:=`(EDO[3], eval(EDO[1], [sin(theta[1](t)) = theta[1](t), sin(theta[2](t)) = theta[2](t), cos(theta[1](t)) = 1, cos(theta[2](t)) = 1, `*`(`^`((D[1](theta[1]))(t), 2)) = 0, `*`(`^`((D[1](theta[2]))(t)...
`:=`(EDO[3], eval(EDO[1], [sin(theta[1](t)) = theta[1](t), sin(theta[2](t)) = theta[2](t), cos(theta[1](t)) = 1, cos(theta[2](t)) = 1, `*`(`^`((D[1](theta[1]))(t), 2)) = 0, `*`(`^`((D[1](theta[2]))(t)...
 

 

`:=`(EDO[4], eval(EDO[2], [sin(theta[1](t)) = theta[1](t), sin(theta[2](t)) = theta[2](t), cos(theta[1](t)) = 1, cos(theta[2](t)) = 1, `*`(`^`((D[1](theta[1]))(t), 2)) = 0, `*`(`^`((D[1](theta[2]))(t)...
`:=`(EDO[4], eval(EDO[2], [sin(theta[1](t)) = theta[1](t), sin(theta[2](t)) = theta[2](t), cos(theta[1](t)) = 1, cos(theta[2](t)) = 1, `*`(`^`((D[1](theta[1]))(t), 2)) = 0, `*`(`^`((D[1](theta[2]))(t)...
 

 

Condiciones Iniciales y Asignaci?n de Valores: 

 

Para poder realizar an?lisis gr?ficos del P?ndulo Doble asignaremos valores a las masas, a los largos de cuerda y a las condiciones iniciales del problema. 

 

Es importante tener en cuenta que las condiciones iniciales deben ser cercanas al punto de equilibrio en torno al cual hemos linealizado el sistema, ya que la linealizaci?n es una buena aproximaci?n siempre y cuando θ1(t) y θ2(t) se muevan en valores cercanos al (0,0). Asumamos los siguientes valores: 

 

`:=`(m[1], 3); -1; `:=`(m[2], 1); -1; `:=`(l[1], 16); -1; `:=`(l[2], 16); -1; `:=`(g, 10); -1 

 

`:=`(CI, theta[1](0) = `/`(1, 10), (D(theta[1]))(0) = 0, theta[2](0) = 0, (D(theta[2]))(0) = 0); -1 

 

Resoluci?n del Sistema de Ecuaciones Diferenciales: 

 

Para resolver el sistema de ecuaciones diferenciales usaremos el comando dsolve. 

 

dsolve([EDO[3], EDO[4], CI], [theta[1](t), theta[2](t)]) 

 

`:=`(theta[1], proc (t) options operator, arrow; `+`(`*`(`/`(1, 20), `*`(cos(`+`(`*`(`/`(1, 2), `*`(sqrt(5), `*`(t))))))), `*`(`/`(1, 20), `*`(cos(`+`(`*`(`/`(1, 6), `*`(sqrt(15), `*`(t)))))))) end pr... 

`:=`(theta[2], proc (t) options operator, arrow; `+`(`-`(`*`(`/`(1, 10), `*`(cos(`+`(`*`(`/`(1, 2), `*`(sqrt(5), `*`(t)))))))), `*`(`/`(1, 10), `*`(cos(`+`(`*`(`/`(1, 6), `*`(sqrt(15), `*`(t)))))))) e... 

 

Frecuencias Naturales de Oscilaci?n: 

 

Viendo las funciones que representan el movimiento del p?ndulo doble, podemos ver que las frecuencias naturales de oscilaci?n del sistema para peque?as oscilaciones son:
 

`:=`(w[1], `+`(`*`(`/`(1, 2), `*`(sqrt(5), `*`(t))))) 

`:=`(w[2], `+`(`*`(`/`(1, 6), `*`(sqrt(15), `*`(t))))) 

 

Gr?ficos de Energ?a: 

 

Energ?a Cin?tica: 

 

Masa 1: 

 

plot(EKM1(t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = green) 

Masa 2: 

 

plot(EKM2(t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = red) 

Masa 1 v/s Masa 2: 

 

plot([EKM1(t), EKM2(t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [green, red]) 

 

plot([EKM1(t), EKM2(t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([EKM1[1](t), EKM2(t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

Sistema: 

 

plot(E[k](t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = blue) 

plot([E[k](t), EKM1(t), EKM2(t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [blue, green, red]) 

 

Energ?a Potencial:
 

Masa 1: 

 

plot(EPM1(t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = green) 

 

Masa 2: 

plot(EPM2(t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = red) 

Masa 1 v/s Masa 2: 

 

plot([EPM1(t), EPM2(t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [green, red]) 

 

plot([EPM1(t), EPM2(t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([EPM1(t), EPM2(t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

 

Sistema: 

 

plot(E[p](t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = blue) 

 

plot([E[p](t), EPM1(t), EPM2(t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [blue, green, red]) 

 

 

Energ?a Cin?tica v/s Energ?a Potencial: 

 

Masa 1: 

 

plot([EKM1(t), EPM1(t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [blue, green]) 

 

plot([EKM1(t), EPM1(t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([EKM1[1](t), EPM1(t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

Masa 2: 

 

plot([EKM2(t), EPM2(t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [blue, green]) 

 

plot([EKM2(t), EPM2(t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([EKM2(t), EPM2(t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

Sistema: 

 

plot([E[k](t), E[p](t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [blue, green]) 

 

plot([E[k](t), E[p](t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([E[k](t), E[p](t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

 

Funcional: 

 

plot(L(t), t = 0 .. `+`(`*`(10, `*`(Pi)))) 

 

Gr?ficos de Movimiento: 

 

Masa 1: 

 

plot(theta[1](t), t = 0 .. `+`(`*`(10, `*`(Pi))), color = green) 

 

plot([theta[1](t), (D(theta[1]))(t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([theta[1](t), (D(theta[1]))(t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

Masa 2: 

 

plot(theta[2](t), t = 0 .. `+`(`*`(10, `*`(Pi)))) 

 

plot([theta[2](t), (D(theta[2]))(t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([theta[2](t), (D(theta[2]))(t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

 

Masa 1 v/s Masa 2: 

 

plot([theta[1](t), theta[2](t), t = 0 .. `+`(`*`(10, `*`(Pi)))]) 

 

animatecurve([theta[1](t), theta[2](t), t = 0 .. `+`(`*`(10, `*`(Pi)))], numpoints = 200, frames = 100, color = red) 

 

plot([theta[1](t), theta[2](t)], t = 0 .. `+`(`*`(10, `*`(Pi))), color = [green, red]) 

 

Simulaci?n: 

 

A continuaci?n realizaremos una animaci?n que representar? el movimiento del p?ndulo doble para peque?as oscilaciones. 

 

`:=`(soporte, plot([[-10, 0], [10, 0]], thickness = 2, color = black)); -1 

`:=`(animacion[1], animate(implicitplot, [`+`(`*`(`^`(`+`(x, `-`(evalf(X[1](t)))), 2)), `*`(`^`(`+`(y, `-`(evalf(Y[1](t)))), 2))) = 4, x = -35 .. 35, y = -35 .. 0, grid = [100, 100], view = [-35 .. 35...
`:=`(animacion[1], animate(implicitplot, [`+`(`*`(`^`(`+`(x, `-`(evalf(X[1](t)))), 2)), `*`(`^`(`+`(y, `-`(evalf(Y[1](t)))), 2))) = 4, x = -35 .. 35, y = -35 .. 0, grid = [100, 100], view = [-35 .. 35...
`:=`(animacion[1], animate(implicitplot, [`+`(`*`(`^`(`+`(x, `-`(evalf(X[1](t)))), 2)), `*`(`^`(`+`(y, `-`(evalf(Y[1](t)))), 2))) = 4, x = -35 .. 35, y = -35 .. 0, grid = [100, 100], view = [-35 .. 35...`:=`(animacion[2], animate(implicitplot, [`+`(`*`(`^`(`+`(x, `-`(evalf(X[2](t)))), 2)), `*`(`^`(`+`(y, `-`(evalf(Y[2](t)))), 2))) = 4, x = -35 .. 35, y = -35 .. 0, grid = [100, 100], view = [-35 .. 35...
`:=`(animacion[2], animate(implicitplot, [`+`(`*`(`^`(`+`(x, `-`(evalf(X[2](t)))), 2)), `*`(`^`(`+`(y, `-`(evalf(Y[2](t)))), 2))) = 4, x = -35 .. 35, y = -35 .. 0, grid = [100, 100], view = [-35 .. 35...
`:=`(animacion[2], animate(implicitplot, [`+`(`*`(`^`(`+`(x, `-`(evalf(X[2](t)))), 2)), `*`(`^`(`+`(y, `-`(evalf(Y[2](t)))), 2))) = 4, x = -35 .. 35, y = -35 .. 0, grid = [100, 100], view = [-35 .. 35...`:=`(animacion[3], animate(plot, [[[0, 0], [`+`(X[1](t), `-`(`*`(2, `*`(cos(`+`(`*`(`/`(1, 2), `*`(Pi)), `-`(theta[1](t)))))))), `+`(Y[1](t), `*`(2, `*`(sin(`+`(`*`(`/`(1, 2), `*`(Pi)), `-`(theta[1](t...
`:=`(animacion[3], animate(plot, [[[0, 0], [`+`(X[1](t), `-`(`*`(2, `*`(cos(`+`(`*`(`/`(1, 2), `*`(Pi)), `-`(theta[1](t)))))))), `+`(Y[1](t), `*`(2, `*`(sin(`+`(`*`(`/`(1, 2), `*`(Pi)), `-`(theta[1](t...`:=`(animacion[4], animate(plot, [[[`+`(X[1](t), `*`(2, `*`(sin(theta[2](t))))), `+`(Y[1](t), `-`(`*`(2, `*`(cos(theta[2](t))))))], [`+`(X[2](t), `-`(`*`(2, `*`(cos(`+`(`*`(`/`(1, 2), `*`(Pi)), `-`(th...
`:=`(animacion[4], animate(plot, [[[`+`(X[1](t), `*`(2, `*`(sin(theta[2](t))))), `+`(Y[1](t), `-`(`*`(2, `*`(cos(theta[2](t))))))], [`+`(X[2](t), `-`(`*`(2, `*`(cos(`+`(`*`(`/`(1, 2), `*`(Pi)), `-`(th...
 

display(soporte, animacion[1], animacion[2], animacion[3], animacion[4], insequence = false, axes = none) 

 

 

Legal Notice: The copyright for this application is owned by the author(s). Neither Maplesoft nor the author are responsible for any errors contained within and are not liable for any damages resulting from the use of this material. This application is intended for non-commercial, non-profit use only. Contact the author for permission if you wish to use this application in for-profit activities.