La ecuación del oscilador armónico amortiguado

Toda la física que necesitas para animar cabe en dos frases: un muelle tira hacia su reposo con una fuerza proporcional a la deformación, y un amortiguador se opone al movimiento con una fuerza proporcional a la velocidad. Súmalas, aplica que la fuerza es masa por aceleración y ya tienes una ecuación diferencial lineal de segundo orden con coeficientes constantes, que es de las pocas que se resuelven exactamente y a mano. Lo que sigue es esa derivación completa, y la reescritura que convierte tres parámetros con unidades en dos números que sí significan algo para quien diseña.

De dos fuerzas a una ecuación

Llamamos x(t) a la posición del objeto y x_obj a la posición de reposo, la que el muelle considera “sin deformar”. Lo cómodo es trabajar con la deformación, u = x - x_obj, porque así el reposo está en cero y las fórmulas quedan limpias. Si el objetivo no se mueve, u'' = x'' y no perdemos nada.

La ley de Hooke dice que la fuerza del muelle es F_k = -k·u. La constante k se llama rigidez y el signo negativo es lo único importante: la fuerza apunta siempre hacia el reposo, se haya estirado el muelle hacia donde se haya estirado. Sin ese signo no hay oscilación, hay una explosión exponencial.

El amortiguamiento viscoso dice que la fuerza del amortiguador es F_c = -c·u', proporcional a la velocidad y opuesta a ella. La constante c se llama coeficiente de amortiguamiento. Es el modelo del pistón en aceite, no el del rozamiento seco: el rozamiento seco tiene módulo constante y produce una ecuación no lineal bastante peor.

Segunda ley de Newton, F = m·a, con la aceleración a = u'':

m·u'' = -k·u - c·u'

Reordenando, la forma canónica:

m·u'' + c·u' + k·u = 0

Tres parámetros: masa, amortiguamiento y rigidez. Todos con unidades, todos acoplados entre sí, ninguno interpretable por separado. Duplicar k no significa “el doble de rápido” mientras no sepas cuánto vale m. Esa es la razón por la que ajustar muelles tocando estos tres números a ojo es tan frustrante, y la razón por la que el siguiente paso no es opcional.

Los dos números que sí significan algo

Divide toda la ecuación por m:

u'' + (c/m)·u' + (k/m)·u = 0

Los dos coeficientes que quedan tienen dimensiones concretas. k/m es el cuadrado de una frecuencia angular, y c/m es una frecuencia angular. Se definen entonces:

ω₀ = √(k/m)            frecuencia natural, en radianes por segundo
ζ  = c / (2·√(k·m))    factor de amortiguamiento, adimensional

Con eso, c/m = 2·ζ·ω₀ y la ecuación se convierte en la forma que aparece en todos los libros:

u'' + 2·ζ·ω₀·u' + ω₀²·u = 0

Tres parámetros con unidades se han convertido en dos, uno de los cuales no tiene unidades en absoluto. La ganancia no es estética.

ω₀ es la velocidad del sistema. Sin amortiguamiento, el objeto oscilaría con periodo T₀ = 2π/ω₀ segundos. Con k = 200 y m = 1 sale ω₀ = 14,14 rad/s, o sea 2,25 Hz, o sea un periodo de 0,444 s. Ese número lo puedes cronometrar y lo puedes comparar con lo que quiere el diseño.

ζ es cuánta energía se pierde por ciclo, y es adimensional, lo cual significa que su valor es comparable entre sistemas con masas y rigideces completamente distintas. ζ = 0,3 se siente igual de rebotón en un muelle rápido que en uno lento. Esa invariancia es la razón por la que ζ es el parámetro que se documenta en un sistema de diseño y c no.

El denominador de ζ, 2·√(k·m), tiene nombre propio: es el amortiguamiento crítico c_c, el valor exacto de c que separa los sistemas que oscilan de los que no. ζ es literalmente la fracción del amortiguamiento crítico que estás usando.

La ecuación característica

La ecuación es lineal, homogénea y de coeficientes constantes, así que se resuelve probando u(t) = e^(r·t). Sustituyendo y sacando factor común e^(r·t), que nunca es cero:

r² + 2·ζ·ω₀·r + ω₀² = 0

Una ecuación de segundo grado corriente. Sus raíces:

r = -ζ·ω₀ ± ω₀·√(ζ² - 1)

Todo el comportamiento del muelle está en el signo de ζ² - 1, es decir, en si ζ es menor, igual o mayor que uno. Esas son las tres regiones, y salen de la aritmética, no de una convención.

Cuando ζ < 1 el radicando es negativo y las raíces son un par conjugado complejo. Escribiendo ω_d = ω₀·√(1 - ζ²), la frecuencia amortiguada, las raíces son -ζ·ω₀ ± i·ω_d y la solución real es:

u(t) = e^(-ζ·ω₀·t) · [ A·cos(ω_d·t) + B·sin(ω_d·t) ]

Una oscilación de frecuencia ω_d dentro de una envolvente exponencial decreciente. Observa que ω_d < ω₀ siempre: amortiguar no solo apaga la oscilación, también la ralentiza. Con ζ = 0,5 la frecuencia baja un 13%; con ζ = 0,9, un 56%.

Cuando ζ = 1 el radicando es cero y hay una raíz doble r = -ω₀. Una raíz doble en una ecuación de segundo orden exige un segundo modo linealmente independiente, que es t·e^(r·t):

u(t) = (A + B·t) · e^(-ω₀·t)

Cuando ζ > 1 el radicando es positivo y hay dos raíces reales distintas, ambas negativas porque √(ζ² - 1) < ζ:

u(t) = A·e^(r₁·t) + B·e^(r₂·t)

Dos exponenciales decrecientes superpuestas, ninguna oscilación.

Fijar las constantes con las condiciones iniciales

A y B no son propiedades del muelle: son las condiciones de partida. Y aquí está lo que hace que un muelle sea interrumpible, porque las condiciones son dos: dónde está el objeto y a qué velocidad va.

Con u(0) = u₀ y u'(0) = v₀, para el caso subamortiguado:

A = u₀
B = (v₀ + ζ·ω₀·u₀) / ω_d

Para el caso crítico:

A = u₀
B = v₀ + ω₀·u₀

Y para el sobreamortiguado, con r₁ y r₂ las dos raíces reales:

A = (v₀ - r₂·u₀) / (r₁ - r₂)
B = (r₁·u₀ - v₀) / (r₁ - r₂)

Con eso ya tienes un muelle exacto, sin integración numérica y sin acumulación de error, evaluable en cualquier instante en tiempo constante:

// Solucion analitica del oscilador amortiguado.
// Devuelve una funcion t -> {x, v} valida para cualquier t >= 0.
function resolverMuelle({ objetivo = 0, x0 = 1, v0 = 0, omega0 = 12, zeta = 0.4 }) {
  const u0 = x0 - objetivo;

  if (zeta < 1) {
    const wd = omega0 * Math.sqrt(1 - zeta * zeta);
    const A = u0;
    const B = (v0 + zeta * omega0 * u0) / wd;
    return (t) => {
      const env = Math.exp(-zeta * omega0 * t);
      const c = Math.cos(wd * t), s = Math.sin(wd * t);
      const u = env * (A * c + B * s);
      const du = env * ((B * wd - zeta * omega0 * A) * c
                      - (A * wd + zeta * omega0 * B) * s);
      return { x: objetivo + u, v: du };
    };
  }

  if (zeta === 1) {
    const A = u0;
    const B = v0 + omega0 * u0;
    return (t) => {
      const env = Math.exp(-omega0 * t);
      return { x: objetivo + (A + B * t) * env, v: (B - omega0 * (A + B * t)) * env };
    };
  }

  const raiz = omega0 * Math.sqrt(zeta * zeta - 1);
  const r1 = -zeta * omega0 + raiz;
  const r2 = -zeta * omega0 - raiz;
  const A = (v0 - r2 * u0) / (r1 - r2);
  const B = (r1 * u0 - v0) / (r1 - r2);
  return (t) => ({
    x: objetivo + A * Math.exp(r1 * t) + B * Math.exp(r2 * t),
    v: A * r1 * Math.exp(r1 * t) + B * r2 * Math.exp(r2 * t),
  });
}

// Un muelle que arranca a 200px del objetivo, en reposo.
const muelle = resolverMuelle({ x0: 200, v0: 0, omega0: 14.14, zeta: 0.5 });
console.log(muelle(0).x.toFixed(2));    // 200.00
console.log(muelle(0.25).x.toFixed(2)); // -32.46  (ya ha pasado de largo)
console.log(muelle(1).x.toFixed(2));    // 0.13    (practicamente en reposo)

Ese código funciona pegado tal cual en cualquier consola. La derivada del caso subamortiguado sale de derivar el producto de la exponencial por la combinación de senos y cosenos; conviene comprobarla numéricamente la primera vez, porque un signo mal puesto ahí produce un muelle que parece correcto en posición y está mal en velocidad, que es el peor error posible: solo se manifiesta al interrumpir.