Oscilador Duffing con un péndulo
En este ejemplo demostramos el uso de la funcionalidad proporcionada en
sympy.physics.mechanics para derivar las ecuaciones de movimiento para un sistema
que consiste en un oscilador con un péndulo. Este ejemplo está inspirado en el
Documento [P.Brzeskia2012] Sección 2.
El sistema se modelará con las ecuaciones de LaGrange. M, es masa del oscilador de duffing, m, es masa del péndulo, l, es longitud del péndulo. k1 y k2, son lineales y partes no lineales de la rigidez de la primavera, y c1, es un coeficiente de amortiguación viscoso del oscilador de duffing.
>>import sympy as sm
>>import sympy.physics.mechanics as me
>>me.init_vprinting()
Definir variables
>>M, m, l, k1, k2, c1, g, h, w, d, r = sm.symbols('M, m, l, k1, k2, c1, g, h, w, d, r')
>>q1, q2 = me.dynamicsymbols('q1 q2')
>>q1d = me.dynamicsymbols('q1', 1)
- h: altura del oscilador de Duffin
- w: ancho del oscilador de Duffin
- d: Profundidad del oscilador Duffing
- r: Radio del enorme bob del péndulo
- q1: Coordenada generalizada que representa la posición del oscilador de Duffing
- q2: Coordenada generalizada que representa el ángulo del péndulo
Definir cinemática
Defina todos los marcos y puntos de referencia.
>>N = me.ReferenceFrame('N')
>>B = N.orientnew('B', 'axis', (q2, N.z))
La velocidad angular del péndulo en el marco de referencia es:
>>B.ang_vel_in(N) >>q2'(t) n_z
Las ubicaciones y las velocidades del bloque de osciladores de Duffing y el péndulo son:
>>O = me.Point('O')
>>block_point = O.locatenew('block', q1 * N.y)
>>pendulum_point = block_point.locatenew('pendulum', l * B.y)
O es un punto fijo en el marco de referencia inercial.
>>O.set_vel(N, 0)
>>block_point.set_vel(N, q1d * N.y)
>>pendulum_point.v2pt_theory(block_point, N, B) >>q1'(t) n_y + -l*q2'(t) b_x
Definir inercia y cuerpos rígidos. Aquí, asumimos un péndulo simple que consiste en un bob de masa m colgando de una cadena de longitud sin masa l y fijado en un punto de pivote (bloqueo del oscilador de duffing).
>>I_block = M / 12 * me.inertia(N, h**2 + d**2, w**2 + d**2, w**2 + h**2)
>>I_pendulum = 2*m*r**2/5*me.inertia(B, 1, 0, 1)
>>block_body = me.RigidBody('block', block_point, N, M, (I_block, block_point))
>>pendulum_body = me.RigidBody('pendulum', pendulum_point, B, m, (I_pendulum, pendulum_point))
Definir fuerzas
Calculamos las fuerzas que actúan sobre el sistema. En este ejemplo, establecemos la energía potencial en cero en el lagrangiano e incluimos el Fuerzas conservadoras (gravedad y resorte de muescas) en las cargas.
>>path = me.LinearPathway(O, block_point)
>>spring = me.DuffingSpring(k1, k2, path, 0)
>>damper = me.LinearDamper(c1, path)
>>loads = spring.to_loads() + damper.to_loads()
>>bodies = [block_body, pendulum_body]
>>for body in bodies:
>> loads.append(me.Force(body, body.mass * g * N.y))
>>loads / _____ 3/2\ / _____ 3/2\ | / 2 / 2\ | | / 2 / 2\ | \k1*\/ q1 + k2*\q1 / /*q1 \- k1*\/ q1 - k2*\q1 / /*q1 [(O, ------------------------------ n_y), (block, -------------------------------- n_y), (O, c1*q1'(t) n_y), (block, -c1*q1'(t) n_y), (block, M*g n_y), (pendulum, g*m n_y)] _____ _____ / 2 / 2 \/ q1 \/ q1
Método de LaGrange
Con la configuración del problema, se puede calcular el Lagrangiano y se forman las ecuaciones de movimiento.
L = me.Lagrangian(N, block_body, pendulum_body)
L 2 2 2 / 2 2 2\ M*q1'(t) m*r *q2'(t) m*\l *q2'(t) - 2*l*sin(q2)*q1'(t)*q2'(t) + q1'(t) / --------- + ------------ + ---------------------------------------------------- 2 5 2
LM = me.LagrangesMethod(L, [q1, q2], bodies=bodies, forcelist=loads, frame=N)
sm.simplify(LM.form_lagranges_equations()) [ / 2 \ / 2\ ] [-M*g + M*q1''(t) + c1*q1'(t) - g*m - m*\l*sin(q2)*q2''(t) + l*cos(q2)*q2'(t) - q1''(t)/ + \k1 + k2*q1 /*q1] [ ] [ / 2 2 \ ] [ m*\5*g*l*sin(q2) + 5*l *q2''(t) - 5*l*sin(q2)*q1''(t) + 2*r *q2''(t)/ ] [ --------------------------------------------------------------------- ] [ 5
No hay comentarios.:
Publicar un comentario