Ecos en la caverna – Forzando al oscilador


Autor:  Alonso Fernández

Hoy Alonso nos trae un fantástico resumen técnico de los primeros detalles de las oscilaciones forzadas y amortiguadas.  Estos elementos, además de ser muy interesantes en física, tienen una importancia mayúscula en cuestiones ingenieriles.  Así que nada mejor que ver cómo se las gasta un ingeniero con estos temas.  Si tenías pensado repasar todo esto pues esta es tu entrada.  

En esta entrada vamos a hablar un poco del oscilador armónico… .  No os preocupéis, en el blog ya se ha hablado de eso en muchas ocasiones. Sobrarían las razones para dedicar otra entrada al tema ya que el oscilador armónico es, al humilde entender de quien esto escribe, uno de los sistemas más ricos que se pueden analizar con las herramientas de la física. Dicho esto,hoy nos vamos a enfrentar a un tema que no se había tratado en la web. Vamos a presentar al oscilador forzado.

Espero que os sea de interés.

SISTEMA MECÁNICO Y ECUACIONES DEL MOVIMIENTO

Y como siempre que presentamos un sistema mecánico, lo hacemos como con el puchero, con todos sus avíos. En física, los avíos suelen ser un dibujito y sus ecuaciones del movimiento. El dibujito es el siguiente:

01

Como se ve en el dibujo, nuestro sistema mecánico es un péndulo simple, que consiste en un punto material de masa m suspendido de un hilo inextensible de longitud l. La diferencia con el caso normal es que en nuestro sistema al péndulo se le aplica un momento conocido M(t).  Es decir, hay un par de fuerzas actuando sobre el sistema.

En principio el sistema queda descrito por las coordenadas X e Y. Sin embargo, podemos usar el hecho de que el hilo es inextensible para emplear una única coordenada en la descripción del movimiento: el ángulo θ.

Una vez explicado el artefacto, vamos a plantear las ecuaciones del movimiento. Para ello debemos hacer dos cosas:

  1. Expresar la posición del punto en función del ángulo θ.

  2. Escribir las ecuaciones en función de dicho ángulo, que llamaremos a partir de ahora coordenada generalizada.

Es fácil ver que, en el sistema de referencia que hemos dibujado, las relaciones entre coordenadas son:

x=l sin \theta

y = -l cos\theta

Hecho esto, sólo nos queda plantear las ecuaciones. Para ello usaremos una vieja conocida: La Lagrangiana. Para los sistemas mecánicos, la Lagrangiana se escribe como la diferencia entre la energía cinética y la energía potencial del sistema:

L=T-V

La energía cinética en este caso vale:

T=\frac{1}{2}m(\dot{x}^2+\dot{y}^2)

Donde el puntito sobre las coordenadas indica derivada respecto al tiempo. Si usamos el cambio de coordenadas antes descrito vemos que la expresión de la energía cinética es:

T=\frac{ml^2}{2}\dot{\theta}^2

La energía potencial es de tipo gravitatorio en este caso, con \vec{g}=-g\vec{j} que usando la coordenada de interés queda:

V=mgy=-mglcos\theta

Con lo que nuestra Lagrangiana, que contiene toda la información del sistema mecánico, es:

L=\frac{ml^2}{2}\dot{\theta}^2+mglcos\theta

Llegados a este punto, sólo nos queda plantear la ecuación de Euler Lagrange para la única coordenada generalizada que tenemos:

\frac{d}{dt}(\frac{\partial L}{\partial \dot{\theta}})-\frac{\partial L}{\partial \theta}=Q_\theta

Donde el término del lado derecho se denomina fuerza generalizada. Cuando un sistema es puramente conservativo, esa fuerza vale cero. Pero en nuestro caso tenemos el momento externo M(t) y, además, en aras de representar un sistema más realista, vamos a tener un mecanismo disipativo del que hablaremos en breve.

Bueno, dicho lo anterior, ya sólo queda tomar la Lagrangiana, calcularle las derivadas y darle forma a la ecuación del movimiento. Después de algunas manipulaciones llegamos a la siguiente ecuación:

ml^2\ddot{\theta}+mglsin\theta=Q_\theta

Vamos a darle ahora forma a la fuerza generalizada. En primer lugar vamos a tratar de aclarar lo de ese mecanismo disipativo. Como sabemos, en los movimientos reales normalmente existe una pérdida de energía. Un péndulo normal, por ejemplo, terminará parándose debido a la fricción con el aire; si pensamos en un sólido deformable vibrando, como una cuerda de guitarra, por ejemplo, además del rozamiento externo con el aire también tendrá lugar un rozamiento interno a nivel molecular (estructural, como decimos los ingenieros) que contribuirá a amortiguar el movimiento. Una forma bastante típica de modelar estos fenómenos de pérdida de energía es mediante la función de disipación de Rayleigh, que es una especie de potencial, con la salvedad de que en vez de la posición, es una función de la velocidad. La forma de la función de Rayleigh es:

U(\dot{\theta})=\frac{c\dot{\theta}^2}{2}

A partir de este potencial, la fuerza generalizada asociada al mismo se calcula como siempre se hace con los potenciales: derivándolo respecto de la coordenada de interés y cambiándolo de signo. Vamos, hacer el gradiente cambiado de signo. Si lo hacemos obtenemos el valor de la fuerza generalizada asociada al amortiguamiento:

Q_{\theta, disip}=-c\dot{\theta}

Para terminar, nos falta calcular la fuerza generalizada asociada al momento aplicado sobre el péndulo. Para calcularla hacemos uso del principio de los trabajos virtuales, que nos dice en pocas palabras que el trabajo virtual asociado a la fuerza generalizada es igual al trabajo de la fuerza (o momento) externo:

\delta W_q=Q\delta q=\delta W_m=M\delta \theta

Nota para los puristas: Soy consciente de que el principio de los trabajos virtuales tiene mucha más pelusilla de la que aquí acabo de exponer, sobre todo a poco que uno se enfrente a un sistema con varios grados de libertad (discreto o contínuo). Sin embargo, dado que esta entrada se centra en el oscilador armónico forzado, tomamos el mínimo necesario para nuestra tarea.

De donde se sigue, ya que la coordenada generalizada es el propio ángulo:

Q=M

Y nuestra ecuación del movimiento es:

ml^2\dot{\theta}+c\dot{\theta}+mgl sin\theta=M

Que, dado que es una ecuación de segundo orden en el tiempo, debemos resolver usando dos condiciones iniciales: Una en posición y otra en velocidad.

MEJOR SEGUIR LA LÍNEA RECTA. AL MENOS AL PRINCIPIO

Ahora vamos a ver la ecuación con un poco de detenimiento. ¿Notáis algo? ¡Exacto! La ecuación es no lineal. Hay un término en el que el ángulo está dentro de una función seno. Y eso es, hablando mal y pronto, jodido. Todos los teoremas de existencia y unicidad y demás cosas maravillosas como el principio de superposición son válidos para casos lineales. El comportamiento no lineal es mucho más oscuro (y rico en matices y en fenómenos asombrosos como el caos). y requiere un tratamiento mucho más personalizado. Para esta entrada nos vamos a centrar en el comportamiento lineal del sistema, dejando para más adelante el análisis no lineal.

Para convertir el sistema en lineal tomamos el término problemático y lo desarrollamos en serie de Taylor:

sin\theta=\theta-\frac{\theta^3}{6}+O(\theta^5)\approx \theta-\frac{\theta^3}{6}

El término O(θ5) significa que si truncamos ahí cometemos errores del orden de la quinta potencia del ángulo. Para ver más claro a dónde queremos llegar veamos las gráficas siguientes:

02

En azul vemos la función seno. En rojo vemos la aproximación lineal, y en verde la aproximación no lineal cúbica. Vemos que, cerca de cero, el seno se aproxima bastante bien por la función lineal, y aun mejor por la función cúbica. Cuando nos alejamos un poco, llega un momento en el que la función lineal deja de ser una buena aproximación, pero la cúbica aun se mantiene en el ajo. Si nos alejamos más, la función cúbica deja de funcionar, y ya necesitaríamos la potencia quinta, y así sucesivamente. Para aclarar esto representemos el error (%) cometido con cada aproximación:

im3-error_taylor

Vemos que, hasta aproximadamente 0.3 rad podemos usar tranquilamente la aproximación lineal, ya que cometeríamos errores de menos del 1 %. La cúbica por su parte alcanza ese valor del error en 1 rad. Por supuesto, el error admisible es un parámetro que debe evaluarse de forma crítica en cada caso. Con esto en mente vamos a suponer que, en primera aproximación, vamos a describir el movimiento en el caso de pequeñas oscilaciones de forma que podemos linealizar tranquilamente.

LIBRE, LIBRE QUIERO SER: RESPUESTA EN VIBRACIÓN LIBRE

Ahora sí, la ecuación nos queda como sigue:

ml^2\dot{\theta}+c\dot{\theta}+mgl\theta=M

Los coeficientes de la ecuación tienen nombre propio:

El término que multiplica a la aceleración se llama masa equivalente y es el responsable de los efectos de inercia del sistema.

El término que multiplica a la velocidad se llama amortiguamiento equivalente y es el responsable de los efectos disipativos.

El término que multiplica a la posición se llama rigidez equivalente y es el que introduce los efectos elásticos (recuperadores) en el sistema.

Por último, en el miembro derecho tenemos la excitación.

Nota: Como algún avispado habrá notado, en este caso de movimiento angular la masa equivalente tiene dimensiones de momento de inercia y la rigidez equivalente de momento de fuerzas.

Los nombres de los términos se eligieron históricamente por analogía con el sistema MASA+MUELLE+AMORTIGUADOR, que es la base para el análisis de vibraciones mecánicas.

Ahora bien, los términos tal cual están ya nos dan la primera información interesante: La frecuencia natural de vibración. En el caso del oscilador armónico usual la frecuencia natural es:

\omega^2_n=\frac{k}{m}

Si usamos los conceptos de masa y rigidez equivalente, para nuestro péndulo, tenemos:

\omega^2_n=\frac{k_{eq}}{m_{eq}}=\dfrac{mgl}{ml^2}=\dfrac{g}{l}

Que es la archifamosa ecuación que demuestra que para pequeñas oscilaciones la frecuencia del péndulo es independiente de la amplitud de la oscilación.

Si re-escribimos la ecuación usando la frecuencia natural tenemos que:

\ddot{\theta}+2\zeta\omega_n\dot{\theta}=\dfrac{M_0}{m_{eq}}=\hat{M}

Donde, ya que estábamos, hemos introducido la relación de amortiguamiento:

\zeta=\dfrac{c}{C_{crit}}=\dfrac{c}{2\sqrt{k_{eq}m_{eq}}}

Aquí el C crítico es el valor del amortiguamiento necesario para que el sistema no oscile.  Esto podemos verlo si resolvemos el sistema sin excitación:

im4-respuesta_amortiguamiento

FORZANDO LA MÁQUINA: VIBRACIONES FORZADAS NO AMORTIGUADAS

Bueno, pues ahora vamos a lo interesante. Vamos a ver qué pasa si la excitación, en vez de ser nula, es una función armónica del tipo:

\hat{M}=\hat{M}_0 exp(i\Omega t)

Nota matemática: Aquí hemos usado la fórmula de Euler, la cual nos enseña:

exp(ix)=cos x + i sin x

Y que de paso sirve para enunciar la, para mí, ecuación más bonita de las matemáticas:

e^{i\pi}+1=0

¿Que por qué es tan bonita? Veréis, contiene las tres operaciones básicas (suma, multiplicación y exponenciación), los elementos neutros de la suma y de la multiplicación (el cero y el uno), y los números más importantes del álgebra, de la geometría y del análisis (el imaginario i, el número π y el número e). Lo que yo os diga, preciosa.

En primer lugar consideremos el caso de amortiguamiento nulo:

\ddot{\theta}+\omega^2_n=\hat{M}_0 exp(i \Omega t)

La solución de este problema, aprovechando que estamos en el régimen lineal, puede escribirse como una solución de la ecuación homogénea, que contendrá la información del transitorio inicial, sumada a una solución particular de la ecuación completa, que contendrá la información del régimen permanente. Resolver eso es sencillo, dado un valor de la frecuencia de excitación y su magnitud. Sin embargo, es más interesante estudiar el comportamiento del sistema al hacer un barrido paramétrico en las frecuencias de excitación. Como no vamos a resolver infinitos problemas, uno para cada frecuencia, vamos a hacer un análisis paramétrico.

Por ello, vamos a centrarnos en la solución en régimen permanente (imaginemos que inicialmente la posición y la velocidad valen cero) y ensayamos una solución del tipo:

\theta=A_\theta exp(i \Omega t)

Donde A es la amplitud del movimiento. Después de derivar y sustituir obtenemos la siguiente expresión para la amplitud de la respuesta:

A_\theta (\omega^2_n-\Omega^2)exp(i \Omega t)=\hat{m} exp(i \Omega t)

A_\theta \omega_n^2(1-\dfrac{\Omega^2}{\omega^2_n})=\dfrac{M_0}{m_eq}

A_\theta=\dfrac{M_0/k_{eq}}{(1-\dfrac{\Omega^2}{\omega^2_n})}

A_\theta=\dfrac{M_0/k_{eq}}{(1-\tau^2)}

Donde hemos definido la relación, sí otra más, de frecuencias \tau=\Omega/\omega_n.

Antes de pintar la función fijémonos en lo que pasa si la frecuencia de excitación, y con ella τ, es cero (es decir, el caso de un problema estático). En ese caso, simplemente tenemos una deformación estática dada por M0/keq, que es exactamente lo que dice la ley de Hooke, base de la teoría del sólido deformable. A ese valor lo llamaremos deflexión estática, con lo que nos queda:

\dfrac{A_\theta}{A_{\theta,est}=\dfrac{1}{1-\tau^2}=\beta}

Donde β puede llamarse factor de amplificación o factor dinámico de carga. El por qué del nombre está más o menos claro, es el valor que nos indica cómo afecta al sistema el cambio de una carga estática por una dinámica. Con todo lo anterior en la cabeza, vamos a representar el factor dinámico de carga en función de la frecuencia:

im5-fdc_zeta0

Nota: Si notáis algo extraño, tranquilos. Efectivamente hay algo extraño, volveremos a ello en breve.

De un vistazo la curva nos enseña la respuesta del sistema ante cualquier carga dinámica. Y eso es maravilloso a poco que escarbemos, porque nos informa de que:

  • A frecuencias muy bajas la respuesta del sistema es prácticamente la estática, esto ya lo hemos comentado antes al hablar de la deflexión estática.

  • A frecuencias muy altas, la respuesta no es que no oscile, ¡es que no existe!. La excitación es tan rápida que el sistema es incapaz de seguirla, por lo que se queda quieto.

  • A medida que nos acercamos a la frecuencia natural del sistema, el factor dinámico de carga crece. De hecho, cuando estamos exactamente en la frecuencia natural, el factor se va al infinito. Se dice entonces que estamos en resonancia. Si resolvemos la ecuación para este caso obtenemos la respuesta siguiente:

im6-resonancia

Sin embargo antes de llegar a la resonancia el sistema nos habría avisado con otro comportamiento llamativo: El batimiento, que tiene esta pinta:

im7-batimiento

Como vemos, el batimiento se caracteriza por una vibración a dos frecuencias: Una de ellas es la de excitación que está muy próxima a la frecuencia natural, y que por tanto es rápida. La otra es una frecuencia que es tanto más lenta cuanto más cerca estamos de la resonancia. Lo que sucede aquí es que la frecuencia de excitación se está combinando con la frecuencia natural, y la combinación es tanto más llamativa cuanto más cerca estamos de la resonancia.

FRENA UN POCO EL CARRO: VIBRACIONES FORZADAS AMORTIGUADAS

Vemos como en resonancia la amplitud de las oscilaciones crece, y crece y crece y así hasta el infinito. Bueno, en realidad hasta el infinito no. Pasan dos cosas: La primera es que bastante antes de llegar al infinito las oscilaciones dejan de ser pequeñas, con lo que la aproximación lineal en la que hemos sustentado todo esto deja de ser válida y, por tanto, no podemos usarla para hacer predicciones sobre lo que le va a pasar a nuestro sistema. La segunda es que, aunque la matemática nos diga que algo se va al infinito, la naturaleza tiene suele tener varios métodos para evitarlo. Uno de ellos es romper cosas o, como sucedería en este caso, doblar el hilo. Otro es introducir algún mecanismo disipativo, como nuestro amortiguamiento. Retomemos, entonces, la ecuación completa:

\ddot{\theta}+2\zeta\omega_n \dot{\theta}+\omega^2_n\theta=\dfrac{M_0}{m_{eq}exp(i \Omega t)}

Para resolver esto de hacemos como antes, proponemos una solución y sustituimos. Sin embargo, ahora es más sutil porque el amortiguamiento hace que tengamos tanto componente real como imaginaria; vamos a verlo.

Proponemos como solución:

\theta = A_\theta exp(i\Omega t)

Derivando una vez:  \dot{\theta}=i\Omega A_\theta exp(i\Omega t)=i\Omega \theta

Derivamos otra vez más:  \ddot{\theta}=i^2\Omega^2 A_\theta exp(i\Omega t)=-\Omega^2 \theta

Y sustituimos en la ecuación:

\left[-\Omega^2+2i\Omega\omega_n\zeta+\omega^2_n \right]A_\theta exp(i\Omega t)=\dfrac{M_0}{m_{eq} exp(i\Omega t)}

Ahora ya podemos cancelar los términos con exponenciales e introducir la deflexión estática y la relación de frecuencias. Si lo hacemos llegamos a:

\dfrac{A_\theta}{A_{\theta,est}=\dfrac{1}{(1-\tau^2)+2i\tau\zeta}=\beta}

Que, como adelantamos, ahora es una función compleja. Las funciones complejas a lo mejor nos asustan al principio, pero a poco que les echemos un poquito de cariño se convierten en amigas incondicionales y herramientas muy poderosas. Y para que lleguen a ser eso, lo primero es tener en cuenta que un número complejo puede representarse como un vector, y como tal, tiene dos propiedades básicas: Módulo (o magnitud) y ángulo (o fase).

im8-plano_z

Si tenemos eso en la cabeza es sencillo calcular el módulo del factor de amplificación y su desfase:

|\beta|=\dfrac{1}{\sqrt{(1-\tau^2)^2+(2\tau\zeta)^2}}

tan(\Phi)=\dfrac{-2\tau\zeta}{1-\tau^2}

Si representamos el factor de amplificación ahora tenemos lo siguiente:

im9-fdc-amortiguado

Mientras que la fase tiene esta pinta:

im10-desfase

Empecemos por el factor de amplificación. Ahí vemos que, para amortiguamientos no nulos, el factor está limitado, en el caso de resonancia al valor:

\beta(\tau=1)=\dfrac{1}{2\zeta}

Curiosamente, en el caso amortiguado el máximo del factor dinámico de carga no se produce en la resonancia, sino muy cerca. Concretamente en:

\tau(\beta_{max}=\sqrt{1-2\zeta^2})

Alcanzándose en ese caso:

\beta_{max}=\dfrac{1}{2\zeta\sqrt{1-\zeta^2}}

En ingeniería, los valores normales de la relación de amortiguamiento se encuentran entre el 2% y el 5%, por tanto, podemos, sin cogernos muchos los dedos, asumir que los sistemas amortiguados tienen el máximo cuando se excitan a la frecuencia natural.

En cuanto a la fase, vemos que toma valores (en la figura está representada en ángulos sexagesimales) negativos. Esto indica que la respuesta SIEMPRE va retrasada respecto a la excitación, lo cual es lógico: El sistema no va a abrir la puerta antes de enterarse de que le están llamando. También vemos que, a medida que sube el amortiguamiento, el retraso se va haciendo cada vez mayor para una misma frecuencia de excitación. Y por último, vemos que en el caso de la resonancia, el retraso siempre vale 90º.

Este es el momento de recordar la nota que dejé cuando presentamos la figura del factor de amplificación para el caso no amortiguado. Ahí, la curva siempre tiene valores positivos, pero para frecuencias superiores a la de resonancia el factor (1-τ2) es negativo. ¿Por qué la curva siempre es positiva? ¡Pues porque a partir de la resonancia el desfase vale 180º! De hecho, si os fijáis en esa figura, en el eje Y aparece que se representa el módulo de β, pero no hemos dicho nada de la fase. Si representamos, en el caso sin amortiguamiento, la solución del problema junto con la excitación antes, durante y después de la resonancia, vemos claramente este comportamiento:

im11-respuesta_tau0-2

im12-respuesta_tau1-0

im13-respuesta_tau3-0

Y eso es todo lo que quería contar hoy. Espero que esta introducción a la mecánica de vibraciones os haya gustado. Al final sólo hemos tratado con un sistema de un grado de libertad, y lineal, pero nos ha servido para aprender muchas cosas que es bueno tener en cuenta desde el principio. Alguien podrá decir: «Eh, es que un sistema de un grado de libertad no tiene aplicación práctica». Si alguien tiene esa inquietud yo le contestaría lo siguiente: Es verdad que, en la práctica, los sistemas de 1 gdl no existen, ya que tratamos con sistemas que en realidad tienen infinitos grados de libertad, y por tanto infinitas frecuencias naturales (e infinitos modos de vibración). Sin embargo, la primera frecuencia natural es la más baja de todas, y por tanto es la que más conviene conocer para asegurar que en el funcionamiento normal no vamos a excitar el sistema de forma dinámica. Y esa frecuencia natural se estima, de la forma más segura posible, (aquí los interesados pueden buscar una cosita que se llama el cociente de Rayleigh) analizando un sistema de un grado de libertad.

De todos modos, queda pendiente para próximas entradas revisar varias cosas: ¿Qué pasa si nos alejamos de la aproximación lineal? ¿Y qué pasa con los sistemas de varios grados de libertad, o que tienen infinitos de ellos? También nos hemos dejado fuera el tema de las vibraciones auto-excitadas, relacionado con la estabilidad de las soluciones… Como veis, aquí hay mucha tela que cortar… Y el patrón maestro nos lo ha vuelto a proporcionar, una vez más, el oscilador armónico.

3 Respuestas a “Ecos en la caverna – Forzando al oscilador

  1. Sí, han aparecido un par de erratas en las ecuaciones. Esa es una de ellas, y las otras están en los cálculos que definen el factor dinámico de carga (en los casos sin y con amortiguamiento). El editor de ecuaciones de vez en cuando gasta estas bromas.

  2. Gracias, muy interesante, aunque reporto un pequeño gazapo tipográfico: en la primera expresión que aparece de la energía cinética T, creo que el segundo sumando no ha de ser «x punto al cuadrado», sino que debería ser «y punto al cuadrado»
    Saludos.

Deja una respuesta

Introduce tus datos o haz clic en un icono para iniciar sesión:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Salir /  Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Salir /  Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Salir /  Cambiar )

Conectando a %s