Resonancias y caos

1. Motivación

El cinturón de asteroides tiene huecos en distancias muy concretas al Sol; Plutón cruza la órbita de Neptuno y sin embargo nunca se le acerca; Júpiter tiene miles de troyanos en sus puntos $L_4$ y $L_5$ pero casi ninguno en los demás. Los tres hechos tienen la misma causa: la resonancia entre periodos orbitales, que una perturbación pequeña convierte en el motor de cambios grandes. Y cuando varias resonancias se acercan, aparece el caos: la pérdida de predictibilidad a largo plazo que Poincaré intuyó en el problema de tres cuerpos y que hoy medimos numéricamente en el Sistema Solar entero.

2. Resonancias de movimiento medio

Dos cuerpos con movimientos medios $n$ y $n_P$ están en resonancia de movimiento medio $p{:}q$ (con $p, q$ enteros pequeños) si $n/n_P \approx p/q$, es decir, si el cuerpo da $p$ vueltas mientras el perturbador da $q$. Por la tercera ley de Kepler, con $a_P$ el semieje del perturbador,

Con $a_J = 5.20289$ UA (elementos medios J2000, JPL) las resonancias con Júpiter caen en

Los huecos de Kirkwood (1867) son las resonancias 3:1, 5:2, 7:3 y 2:1, prácticamente vacías; pero 3:2 (Hilda), 4:3 y 1:1 están pobladas. La resonancia por sí sola no vacía ni llena: lo que decide es la dinámica del ángulo resonante.

2.1 El ángulo resonante y el modelo del péndulo

Sean $\lambda = M + \varpi$ la longitud media (lección 2) y $\varpi$ la longitud del periapsis del cuerpo perturbado. Para una resonancia $(p+q){:}p$, con $\lambda_{\text{ext}}$ la longitud media del cuerpo exterior (el más lento) y $\lambda_{\text{int}}$ la del interior, el ángulo resonante en su forma más simple es

cuyos coeficientes suman cero (regla de d’Alembert) y cuya derivada $\dot\phi = (p+q)n_{\text{ext}} - p\,n_{\text{int}}$ se anula exactamente en la resonancia.

Lejos de la resonancia $\phi$ circula (crece sin límite, como $\lambda_P$ y $\lambda$); en la resonancia libra alrededor de un valor fijo. Al desarrollar la función perturbadora $\mathcal R$ (lección 4) en serie de Fourier de los ángulos, el término $\propto e^q\cos\phi$ es el único que no promedia a cero cerca de la resonancia, porque $\dot\phi \approx 0$. Conservando ese término, las ecuaciones de Lagrange se reducen a un grado de libertad (Murray & Dermott, §8.8): con $J$ una acción proporcional a la desviación de $a$ respecto de $a_{res}$,

que es el hamiltoniano del péndulo. Las curvas $H =$ const son libraciones (dentro de la separatriz) o circulaciones (fuera). La separatriz $H = \epsilon$ tiene semianchura en $J$ y la libración de pequeña amplitud tiene frecuencia

Como $\epsilon\propto e^q$, la anchura de la resonancia crece con la excentricidad: una resonancia de orden $q$ alto es estrecha para órbitas casi circulares. Y como $\omega_0 \ll n$, la libración es lenta: miles o decenas de miles de años en el Sistema Solar.

2.2 Plutón y Neptuno

Horizons da en J2000 $a_N = 30.0939$ UA, $a_P = 39.5713$ UA, $T_N = 60\,298$ d y $T_P = 90\,922$ d: $T_P/T_N = 1.5079$. La resonancia 3:2 exacta estaría en $a_{res} = a_N(3/2)^{2/3} = 39.434$ UA: Plutón está dentro de su anchura, y el cociente osculador oscila alrededor de 1.5 exactamente. El ángulo resonante $\phi = 3\lambda_P - 2\lambda_N - \varpi_P$ (Plutón es el cuerpo exterior: $\dot\phi = 3n_P - 2n_N = 0$) libra alrededor de $180°$ con amplitud $\approx 80°$ y periodo $\approx 20\,000$ años (Cohen & Hubbard 1965, que lo descubrieron integrando numéricamente). Consecuencia geométrica: en una conjunción Plutón–Neptuno ($\lambda_P = \lambda_N$) se tiene $\phi = \lambda_P - \varpi_P$, de modo que $\phi \approx 180°$ obliga a que las conjunciones ocurran cuando Plutón está en su afelio. Aunque el perihelio de Plutón (29.70 UA) queda dentro de la órbita de Neptuno (30.09 UA), las conjunciones se producen siempre cerca del afelio, a más de 17 UA de distancia. Malhotra (1993) explicó cómo Plutón fue capturado en esta resonancia durante la migración hacia fuera de Neptuno, un mecanismo que también da cuenta de los cientos de «plutinos» del cinturón de Kuiper.

2.3 Los huecos de Kirkwood

¿Por qué 3:1 está vacía y 3:2 llena? Wisdom (1983) construyó una aplicación simpléctica para la resonancia 3:1 y descubrió que la zona resonante es caótica: un asteroide con $e \approx 0.1$ puede saltar en $\sim10^4$–$10^5$ años a $e > 0.3$, cruzar la órbita de Marte y ser eliminado por un encuentro cercano. En 3:2, en cambio, la libración mantiene las conjunciones lejos de Júpiter (como con Plutón) y el ángulo resonante protege. La diferencia entre proteger y expulsar depende de la fase de las conjunciones y de la presencia de resonancias seculares superpuestas ($\nu_5$, $\nu_6$), que es precisamente la situación en que el péndulo deja de bastar y hace falta el mapa estándar.

3. El mapa estándar de Chirikov

3.1 Del rotor pateado al mapa

Consideremos un rotor libre (acción $p$, ángulo $\theta$) que recibe un impulso periódico de intensidad $K\sin\theta$:

Entre patadas, $p$ es constante y $\theta$ avanza $p$; en cada patada $p$ salta $K\sin\theta$. Tomando el estado justo después de cada patada se obtiene el mapa estándar

(implementado en src/lib/math/resonancias.ts). Es el modelo mínimo de una resonancia perturbada: para $K = 0$ es integrable ($p$ constante, cada valor de $p$ es un toro invariante sobre el que $\theta$ gira con número de rotación $p/2\pi$); el término $K\sin\theta$ es la resonancia $1{:}1$ centrada en $p = 0$, y sus réplicas en $p = 2\pi m$.

3.2 Propiedades

Simpléctico. El jacobiano es $J = \begin{pmatrix}1 + K\cos\theta & 1\\ K\cos\theta & 1\end{pmatrix}$ (filas $\partial(\theta', p')/\partial\theta$, $\partial(\theta',p')/\partial p$), con $\det J = 1$: conserva el área del espacio de fases, como todo flujo hamiltoniano (Problema 2). También es reversible: conjugado a su inversa por la reflexión $(\theta,p)\mapsto(-\theta,p)$.

Puntos fijos. $(\theta, p) = (0, 0)$ y $(\pi, 0)$ (y sus trasladados en $p = 2\pi m$). La estabilidad lineal viene de la traza de $J$: un punto fijo de una aplicación de área 1 es elíptico (estable) si $|\operatorname{tr}J| < 2$ e hiperbólico si $|\operatorname{tr}J| > 2$. En $(0,0)$, $\operatorname{tr}J = 2 + K > 2$: hiperbólico para todo $K > 0$. En $(\pi, 0)$, $\operatorname{tr}J = 2 - K$: elíptico para $0 < K < 4$, con frecuencia de libración $\omega$ dada por $2\cos\omega = 2 - K$, es decir, $\cos\omega = 1 - K/2$ (Problema 3); para $K \ll 1$, $\omega \approx \sqrt K$, que es la $\omega_0 = \sqrt{\epsilon\beta}$ del péndulo con $\epsilon = K$, $\beta = 1$.

Anchura de la resonancia principal. Para $K$ pequeño, la separatriz del péndulo asociado tiene semianchura $\Delta p = 2\sqrt K$.

3.3 Solapamiento de resonancias y transición al caos

Las resonancias en $p = 0$ y $p = 2\pi$ tienen anchura total $4\sqrt K$ cada una y distan $2\pi$. El criterio de Chirikov dice que cuando sus separatrices se tocan, las órbitas pueden pasar de una a otra y el caos se extiende por todo $p$: $4\sqrt K = 2\pi$ da $K \approx (\pi/2)^2 = 2.47$. Esta estimación ignora las resonancias secundarias ($p = \pi$, de anchura $\propto K$) y la capa caótica que rodea cada separatriz; incluyéndolas, Chirikov (1979) obtuvo $K_c \approx 1$. El valor preciso lo dio Greene (1979): el último toro invariante que separa $p < 0$ de $p > 0$ es el de número de rotación «más irracional», el número áureo $(\sqrt5-1)/2$, y se destruye en

Para $K < K_g$ el caos existe pero está confinado entre toros invariantes (en $p$ no puede crecer más que $2\pi$); para $K > K_g$ hay difusión ilimitada en $p$. Para $K \gg 1$ casi todo el espacio es caótico y el exponente de Lyapunov tiende a $\lambda \approx \ln(K/2)$ (Chirikov 1979): con $K = 5$, $\lambda \approx 0.92$ por iteración, comprobado numéricamente en resonancias.test.ts.

Reproduce la tabla con el simulador: usa los presets de $K$, pulsa en el retrato para lanzar órbitas y lee $\lambda$ de la última. El panel de difusión compara $\langle(\Delta p)^2\rangle$ con la estimación cuasilineal $K^2 n/2$: por debajo de $K_g = 0.971635$ queda acotada; por encima crece sin límite.

4. El teorema KAM

El mapa estándar ilustra el contenido del teorema de Kolmogorov–Arnold–Moser (1954–1963). Sea $H = H_0(\mathbf J) + \epsilon H_1(\mathbf J, \boldsymbol\theta)$ un sistema casi integrable en variables acción-ángulo, con $H_0$ no degenerado ($\det\partial^2H_0/\partial\mathbf J^2 \ne 0$). Entonces, para $\epsilon$ suficientemente pequeño, sobreviven (ligeramente deformados) los toros invariantes de $H_0$ cuyas frecuencias $\boldsymbol\omega = \partial H_0/\partial\mathbf J$ son diofánticas:

es decir, «suficientemente alejadas» de toda resonancia $\mathbf k\cdot\boldsymbol\omega = 0$. La medida del conjunto de toros destruidos tiende a cero con $\epsilon$. Los toros resonantes y sus vecinos se rompen en cadenas de islas y capas caóticas (teorema de Poincaré–Birkhoff), y el último en romperse es el de frecuencia «más irracional»: exactamente lo que Greene encontró. Dos matices importantes en astronomía:

• En sistemas con dos grados de libertad (como el mapa estándar o el problema restringido plano) los toros son barreras: el caos queda confinado hasta que se rompe el último toro. Con tres o más, los toros no separan el espacio de fases y hay difusión de Arnold a través de la red de resonancias, aunque extraordinariamente lenta.
• El problema planetario es degenerado (el hamiltoniano kepleriano depende de una sola acción, $L = \sqrt{\mu a}$) y hay que recurrir a versiones del teorema para sistemas degenerados (Arnold 1963): la estabilidad KAM del Sistema Solar solo se ha probado para masas planetarias muchísimo menores que las reales.

En la práctica, el Sistema Solar es caótico: Laskar (1989) y Sussman & Wisdom (1992) midieron un tiempo de Lyapunov de unos 5 millones de años para los planetas interiores. No significa que Mercurio vaya a escapar mañana (las acciones $a$ cambian poco: el caos es sobre todo en las fases), sino que una incertidumbre de 15 m en la posición actual de la Tierra crece hasta 150 millones de km en 100 millones de años: las efemérides no pueden extenderse más allá de ese horizonte.

5. Análisis numérico

Exponente de Lyapunov de un mapa. El método de Benettin et al. (1980) propaga un vector tangente $\mathbf v$ con el jacobiano $J_k$ en cada iteración, mide su crecimiento y lo renormaliza: $\lambda = \lim_{N\to\infty}\frac1N\sum_{k=1}^N\ln\lVert J_k\mathbf v_{k-1}\rVert$ con $\lVert\mathbf v_k\rVert = 1$. La convergencia es lenta, $\propto N^{-1/2}$ (los sumandos son como variables aleatorias), y para órbitas regulares $\lambda_N \sim \ln N/N \to 0$: distinguir $\lambda = 0$ de $\lambda = 10^{-4}$ exige $N \gtrsim 10^5$ iteraciones; por eso en problemas reales se usan indicadores más rápidos (MEGNO, FLI, análisis de frecuencias de Laskar).

Redondeo y reversibilidad. El mapa en coma flotante no es exactamente simpléctico ni reversible (el módulo $2\pi$ y el redondeo rompen la estructura en $\sim10^{-16}$ por paso), pero, a diferencia de un integrador de EDOs, no introduce un error de truncamiento: cada iteración es «exacta» salvo redondeo. Para pruebas de reversibilidad se puede formular el mapa en enteros (aritmética modular en $\mathbb Z_N$), que es exactamente reversible. En órbitas caóticas el redondeo se amplifica en $e^{\lambda N}$ y la órbita calculada se separa de la verdadera en $\sim\ln(10^{16})/\lambda \approx 37/\lambda$ iteraciones; la justificación de que el resultado estadístico sigue siendo válido es el lema de sombreado (shadowing): existe una órbita verdadera cercana a la calculada, aunque con otra condición inicial.

Resonancias en integradores. Un integrador simpléctico de paso $h$ introduce su propia resonancia cuando $h$ es conmensurable con un periodo del sistema (resonancias de paso): el leapfrog con $h = T/4$ sobre una órbita kepleriana produce artefactos que desaparecen con $h = T/4.1$. Es un caso real de «mapa estándar» inducido por el método numérico.

6. Simulador

El simulador integra $N$ cuerpos con leapfrog y estima el exponente de Lyapunov a partir de la divergencia de dos condiciones iniciales vecinas. Preguntas:

1. En el ocho de Chenciner–Montgomery (órbita periódica estable), comprueba que la distancia entre trayectorias vecinas crece linealmente, no exponencialmente: $\lambda\to0$.
2. En el problema pitagórico, estima $\lambda$ y el horizonte de predictibilidad $t_h \approx \ln(\delta_{\max}/\delta_0)/\lambda$ para $\delta_0 = 10^{-10}$.
3. Reduce el paso a la mitad y comprueba que $\lambda$ no cambia apreciablemente (es una propiedad del sistema, no del integrador), mientras que la trayectoria concreta sí.

Problemas

1. (Demostración) Deduce $a_{res} = a_P(q/p)^{2/3}$ y calcula las resonancias 2:1, 3:1 y 5:2 con Júpiter. Un asteroide con $a = 2.50$ UA está en 3:1: ¿cuántas vueltas da entre dos conjunciones consecutivas con Júpiter y dónde ocurren éstas si $\phi$ libra en $0$?

2. (Demostración) Sea $T$ el mapa estándar. Demuestra que conserva el área y que es reversible: con la involución $R(\theta, p) = (-\theta,\ p + K\sin\theta)$ se cumple $R\circ T\circ R = T^{-1}$.

3. (Demostración) Linealiza el mapa en $(\pi, 0)$ y demuestra que los multiplicadores son $e^{\pm i\omega}$ con $\cos\omega = 1 - K/2$ si $0 < K < 4$, y reales si $K > 4$. ¿Qué ocurre en $K = 4$? Calcula $\omega$ para $K = 0.5$ y $K = 2$.

4. (Péndulo) Para $H = \tfrac12\beta J^2 - \epsilon\cos\phi$, demuestra que la separatriz es $J = \pm2\sqrt{\epsilon/\beta}\cos(\phi/2)$ y que el periodo de libración de pequeña amplitud es $2\pi/\sqrt{\epsilon\beta}$. Si para Plutón el periodo de libración es $\approx 20\,000$ años y $T_P = 248$ años, ¿cuál es el cociente $\omega_0/n_P$?

5. (Chirikov) Reproduce la estimación $K \approx (\pi/2)^2$ del solapamiento de las resonancias principales y mejórala incluyendo la resonancia secundaria en $p = \pi$ (orden 2, anchura $\propto K$): demuestra que la condición de contacto se convierte en $2\sqrt K + K/2 \cdot c = \pi$ para cierta constante $c$ y estima $K$ con $c = 1$.

6. (Computacional) Programa el mapa estándar y: (a) dibuja el retrato de fases ($10^3$ iteraciones de 50 condiciones iniciales) para $K = 0.5, 0.97, 2, 5$; (b) calcula $\lambda$ por el método de Benettin con $N = 10^5$ para una órbita caótica en cada $K$ y compara con $\ln(K/2)$; (c) para $K = 0.9$ y $K = 1.05$, lanza una órbita en $(\theta,p) = (0.1, 0.1)$ y mide $\max|p|$ tras $10^6$ iteraciones: ¿qué cambia al cruzar $K_g$?

import numpy as np
def mapa(th, p, K):
    p = p + K*np.sin(th); return (th + p) % (2*np.pi), p
def lyapunov(th, p, K, N):
    v = np.array([1.0, 0.0]); s = 0.0
    for _ in range(N):
        c = K*np.cos(th); J = np.array([[1+c, 1.0], [c, 1.0]])
        v = J @ v; nv = np.linalg.norm(v); s += np.log(nv); v /= nv
        th, p = mapa(th, p, K)
    return s/N
for K in [0.5, 0.97, 2, 5]:
    print(K, lyapunov(0.1, 0.1, K, 100_000), np.log(K/2))

7. (Kirkwood y Hilda) Usa las fórmulas de la sección 2 para explicar cualitativamente por qué en 3:2 las conjunciones con Júpiter ocurren en el perihelio del asteroide (protección) mientras que en 3:1 alternan perihelio y afelio. Busca en Horizons (o en el MPC) los elementos de (153) Hilda y comprueba $T/T_J$.

Autoevaluación

Lecturas

• Murray & Dermott (1999), cap. 8: resonancias de movimiento medio, el modelo del péndulo (§8.8) y los casos de Plutón e Hilda; cap. 9: solapamiento de resonancias, mapa estándar y caos en el Sistema Solar (§9.3–9.6).
• Lichtenberg & Lieberman (1992), caps. 3–4: el mapa estándar, el criterio de Chirikov, el método de Greene y el teorema KAM con demostración esquemática.
• Chirikov (1979), DOI 10.1016/0370-1573(79)90023-1, y Greene (1979), DOI 10.1063/1.524170: los dos artículos fundacionales.
• Wisdom (1983), DOI 10.1016/0019-1035(83)90127-6: la explicación caótica del hueco 3:1.
• Cohen & Hubbard (1965), DOI 10.1086/109674, y Malhotra (1993), DOI 10.1038/365819a0: la resonancia de Plutón y su origen.
• Laskar (1989), DOI 10.1038/338237a0, y Sussman & Wisdom (1992), DOI 10.1126/science.257.5066.56: el tiempo de Lyapunov del Sistema Solar.