35. Integración numérica de órbitas y elección del paso temporal

35.1 De la ecuación continua a pasos discretos

La dinámica orbital satisface dr⃗/dt = v⃗ y dv⃗/dt = a⃗(r⃗). Un integrador construye estados aproximados en t, t + Δt, t + 2Δt y así sucesivamente.

Reducir Δt suele mejorar el resultado, pero aumenta la cantidad de evaluaciones y no corrige por sí solo un algoritmo mal implementado.

35.2 Qué representa el paso temporal

Durante Δt suponemos que una aproximación local describe suficientemente bien la evolución. El paso debe resolver la escala más rápida relevante, no solamente la duración total de la misión.

En una órbita excéntrica, el periapsis exige pasos menores porque la velocidad y la aceleración cambian con mayor rapidez.

35.3 Euler explícito

El método más directo usa las derivadas al comienzo del paso:

Es de primer orden: el error global disminuye aproximadamente en proporción a Δt.

35.4 Por qué Euler falla en órbitas

Euler desplaza la posición en la tangente y luego cambia la velocidad, sin curvar el recorrido dentro del paso. En una órbita circular suele agregar energía y generar una espiral hacia afuera.

Es excelente para comprender la discretización y construir pruebas, pero rara vez es una buena elección orbital de producción.

35.5 Euler-Cromer o semimplícito

Una modificación mínima actualiza primero la velocidad:

Sigue siendo de primer orden, pero en muchos sistemas ligados mantiene la energía acotada mejor que Euler explícito.

35.6 Verlet de velocidad

Verlet incorpora la curvatura de la posición y promedia aceleraciones:

Es de segundo orden y requiere evaluar la aceleración en la posición nueva.

35.7 Leapfrog

Leapfrog puede expresarse como medio impulso, deriva completa y otro medio impulso. Posiciones y velocidades “saltan” unas sobre otras en el tiempo.

Verlet de velocidad y leapfrog están estrechamente relacionados. Son simplécticos para Hamiltonianos separables y muy utilizados en dinámica gravitatoria de largo plazo.

35.8 Runge–Kutta de cuarto orden

RK4 combina cuatro estimaciones de la derivada dentro de cada paso. Su error global es de orden Δt⁴ para problemas suaves.

Ofrece gran precisión a corto y mediano plazo, pero requiere cuatro evaluaciones del campo y no es simpléctico. A tiempos muy largos puede mostrar deriva secular de invariantes.

35.9 Orden no significa estabilidad ilimitada

Un método de orden alto con Δt excesivo puede ser peor que uno sencillo con un paso razonable. El orden describe cómo disminuye el error al refinar el paso dentro del régimen de convergencia.

También importan estabilidad, costo por paso, conservación geométrica y comportamiento durante encuentros cercanos.

35.10 Integradores simplécticos

Un integrador simpléctico preserva la estructura geométrica del flujo hamiltoniano. No conserva exactamente la energía en cada paso, pero suele hacerla oscilar alrededor del valor correcto.

Esta propiedad es especialmente valiosa al estudiar miles o millones de períodos orbitales.

35.11 Elegir Δt a partir del período

Para una primera prueba puede usarse una fracción del período mínimo:

Valores como 100 a 1.000 pasos por período son puntos de partida, no garantías. La excentricidad, perturbaciones y precisión exigida pueden requerir mucho más.

35.12 Escala dinámica local

Una escala temporal útil cerca de una masa central es:

Es el período local dividido por 2π. Un criterio adaptativo puede limitar Δt a una fracción pequeña de τ y reducirlo automáticamente al disminuir r.

35.13 Pasos adaptativos

Los métodos adaptativos estiman el error local comparando aproximaciones de distinto orden o un paso completo contra dos medios pasos. Ajustan Δt para satisfacer tolerancias absolutas y relativas.

Son útiles en sobrevuelos y órbitas muy excéntricas, aunque combinar adaptatividad con propiedades simplécticas requiere técnicas específicas.

35.14 Error local y error global

El error local se introduce en un paso suponiendo exacto su estado inicial. El error global acumula y propaga los errores de todos los pasos anteriores.

Para un método de orden p, el error local suele ser O(Δt^p+1) y el global O(Δt^p).

35.15 Convergencia mediante refinamiento

Ejecutá la misma simulación con Δt, Δt/2 y Δt/4. Si las soluciones se acercan de forma compatible con el orden del método, existe evidencia de convergencia.

Comparar solamente una imagen de la órbita puede ocultar diferencias de fase; conviene medir el estado en tiempos comunes e invariantes físicos.

35.16 Euler y Euler-Cromer con JavaScript

function aceleracion(estado, mu) {
  const r3 = Math.hypot(estado.x, estado.y) ** 3;
  return { x: -mu * estado.x / r3, y: -mu * estado.y / r3 };
}

function pasoEuler(estado, dt, mu) {
  const a = aceleracion(estado, mu);
  return { x: estado.x + estado.vx * dt,
    y: estado.y + estado.vy * dt,
    vx: estado.vx + a.x * dt, vy: estado.vy + a.y * dt };
}

function pasoEulerCromer(estado, dt, mu) {
  const a = aceleracion(estado, mu);
  const vx = estado.vx + a.x * dt, vy = estado.vy + a.y * dt;
  return { x: estado.x + vx * dt, y: estado.y + vy * dt, vx, vy };
}

const inicial = { x: 1, y: 0, vx: 0, vy: 1 };
console.log("Euler:", pasoEuler(inicial, 0.1, 1));
console.log("Euler-Cromer:", pasoEulerCromer(inicial, 0.1, 1));

35.17 Verlet de velocidad con JavaScript

function pasoVerlet(estado, dt, mu) {
  const campo = ({ x, y }) => {
    const r3 = Math.hypot(x, y) ** 3;
    return { x: -mu * x / r3, y: -mu * y / r3 };
  };
  const a0 = campo(estado);
  const x = estado.x + estado.vx * dt + 0.5 * a0.x * dt ** 2;
  const y = estado.y + estado.vy * dt + 0.5 * a0.y * dt ** 2;
  const a1 = campo({ x, y });
  return { x, y,
    vx: estado.vx + 0.5 * (a0.x + a1.x) * dt,
    vy: estado.vy + 0.5 * (a0.y + a1.y) * dt };
}

let estado = { x: 1, y: 0, vx: 0, vy: 1 };
const dt = 0.01;
for (let i = 0; i < 1000; i++) estado = pasoVerlet(estado, dt, 1);
console.log(estado);

35.18 Laboratorio: método, paso y error

Todos los métodos parten del mismo estado circular normalizado: μ = 1, r = 1 y v = 1. Aumentá Δt para hacer visibles sus diferencias. Reiniciá al cambiar de método para una comparación justa.

35.19 Qué observar en el laboratorio

• Euler aumenta el radio y la energía de forma sistemática;
• Euler-Cromer deforma la órbita, pero suele mantenerla ligada;
• Verlet produce un error de energía pequeño y oscilatorio;
• RK4 resulta muy preciso con pasos moderados, aunque no es simpléctico;
• todos los métodos fallan si Δt deja de resolver la dinámica.

35.20 Tiempo simulado y tiempo de animación

La velocidad visual no debe cambiar las ecuaciones. Podemos ejecutar varios pasos físicos por cuadro, pero cada uno conserva el mismo Δt.

Usar directamente el intervalo irregular entre cuadros como paso físico dificulta la reproducibilidad. Es preferible un acumulador de tiempo y pasos físicos fijos.

35.21 Encuentros, colisiones y eventos

Un integrador puede saltar desde un lado de un planeta al otro sin detectar su superficie. Los eventos necesitan localizar cruces dentro del paso, reducir Δt o interpolar la trayectoria.

Una colisión bien detectada no implica que el modelo posterior sea correcto: todavía hay que definir rebote, fusión, fragmentación o impacto.

35.22 Ejercicio propuesto

En unidades μ = 1, una partícula parte de r⃗ = (1, 0), v⃗ = (0, 1). Aplicá un paso de Euler explícito con Δt = 0,1.

1. Calculá la aceleración inicial.
2. Obtené el nuevo estado.
3. Compará la energía específica inicial y final.

const dt = 0.1;
const estado = { x: 1, y: 0, vx: 0, vy: 1 };
const aceleracion = { x: -1, y: 0 };
const nuevo = {
  x: estado.x + estado.vx * dt,
  y: estado.y + estado.vy * dt,
  vx: estado.vx + aceleracion.x * dt,
  vy: estado.vy + aceleracion.y * dt
};
const energia0 = 0.5 * (estado.vx ** 2 + estado.vy ** 2) - 1;
const energia1 = 0.5 * (nuevo.vx ** 2 + nuevo.vy ** 2)
  - 1 / Math.hypot(nuevo.x, nuevo.y);

console.log(nuevo);
console.log("energía inicial:", energia0.toFixed(5));
console.log("energía final:", energia1.toFixed(5));
console.log("error relativo:", ((energia1 - energia0) / Math.abs(energia0) * 100).toFixed(2), "%");

35.23 Errores frecuentes

• elegir Δt por la velocidad deseada de la animación;
• actualizar x antes de terminar una fórmula que requiere el estado antiguo;
• comparar métodos con pasos o estados iniciales diferentes;
• creer que una órbita visualmente cerrada garantiza precisión;
• usar un paso fijo grande durante el periapsis o un encuentro;
• suponer que orden alto equivale a conservación perfecta a largo plazo;
• no repetir la simulación con pasos menores.

35.24 Ideas para recordar

• Integrar es aproximar la evolución continua mediante pasos finitos.
• Euler es didáctico, pero introduce una fuerte deriva orbital.
• Euler-Cromer mejora mucho el comportamiento ligado con un cambio mínimo.
• Verlet y leapfrog son opciones simplécticas de segundo orden.
• RK4 ofrece alta precisión por paso, con mayor costo y sin estructura simpléctica.
• El paso debe resolver la escala dinámica más rápida.
• La convergencia se comprueba repitiendo con Δt cada vez menor.

En el próximo tema aprenderemos a verificar la simulación mediante energía, momento angular y leyes de Kepler.