Una manera alternativa de escribir la ecuación de movimiento del problema de los N cuerpos, es reconociendo que cada una las fuerzas gravitacionales que actúan sobre las partículas derivan de un potencial (son fuerzas conservativas). Así:
$$ \left\{m_i\ddot{\vec r}_i= -\frac{\partial}{\partial\vec{r}_i}\left(\sum_{j\neq i} m_i V_{j}\right)\right\}_{N}, $$donde $V_j=-\mu_j/r_{ij}$ es el potencial gravitacional de la partícula $j$ (ver potencial_campo_gravitacional), medido en la posición de la partícula $i$, y la derivada vectorial $\partial/\partial\vec{r}_i$ es, por definición en este libro, una representación del operador gradiente, $\vec{\nabla}_i$, calculado con respecto a las coordenadas de la partícula $i$ (ver la Sección Derivadas.)
Si introducimos lo que podemos llamar la energía potencial gravitacional mutua $U_{ij}\equiv m_i V_j$ (ver Ec. energia_fuerza_campo_gravitacional), la ecuación de movimiento se puede escribir como:
\begin{equation} \label{eq:ncuerpos_edm_potencial_individual} \left\{m_i\ddot{\vec r}_i= -\frac{\partial}{\partial\vec{r}_i} \left(\sum_{j\neq i} U_{ij}\right)\right\}_{N}, \end{equation}Si bien, la función de energía potencia $\sum_{j\neq i} U_{ij}$, que aparece en el lado derecho de las Ecs. (ncuerpos_edm_potencial_individual), es unica para cada partícual del sistema, existe una interesante propiedad del problema de los N cuerpos (y otros problemas con fuerzas de interacción centrales) que permite calcular las fuerzas en el sistema, como el gradiente de una única función de energía potencial. Aunque esta propiedad es presentada como un resultado trivial en casi todos los libros de mecánica celeste, esta lejos de serlo por lo que profundizaremos un momento en ella.
Consideremos el caso particular de un sistema de tres cuerpos. En términos explícitos, las e.d.m. de este sistema (Ecs. ncuerpos_formulacion_ecuaciones) tienen la forma :
\begin{eqnarray} \nonumber m_1 \ddot{r}_1 & = & -\frac{m_1\mu_2}{r_{12}^3}\vec{r}_{12}-\frac{m_1\mu_3}{r_{13}^3}\vec{r}_{13}\\ \nonumber m_2 \ddot{r}_2 & = & -\frac{m_2\mu_1}{r_{21}^3}\vec{r}_{21}-\frac{m_2\mu_3}{r_{23}^3}\vec{r}_{23}\\ \nonumber m_3 \ddot{r}_3 & = & -\frac{m_3\mu_1}{r_{31}^3}\vec{r}_{31}-\frac{m_3\mu_2}{r_{32}^3}\vec{r}_{32}\\ \end{eqnarray}Nótese, por ejemplo, que las fuerzas $\vec{F}_{12}=-m_1\mu_2 \vec{r}_{12}/r_{12}^3$ y $\vec{F}_{21}=-m_2\mu_1 \vec{r}_{21}/r_{21}^3$, en realidad pueden obtenerse de una sola energía potencial, a saber, $U_{12}=-m_1\mu_2/r_{12}$:
$$ \begin{array}{rcl} \vec{F}_{12} & = & \partial_{\vec r_1}U_{12} \\ \vec{F}_{21} & = & \partial_{\vec r_2}U_{12} \\ \end{array} $$La fuerza $\vec{F}_{12}$ se obtiene de $U_{12}$ si se deriva esta cantidad respecto a $\vec r_1$; por su lado, la fuerza $\vec{F}_{21}$ aparece deribando $U_{12}$ respecto a $\vec r_2$. Una situación similar ocurre con los pares de fuerzas $\vec{F}_{13}$, $\vec{F}_{31}$ y $\vec{F}_{23}$, $\vec{F}_{32}$.
Es decir, 6 fuerzas diferentes derivan, en realidad, de solo 3 términos de energía potencial: $U_{12}$, $U_{13}$ y $U_{23}$. Explícitamente:
\begin{eqnarray} \nonumber m_1 \ddot{r}_1 & = & -\frac{\partial U_{12}}{\partial\vec{r}_1}-\frac{\partial U_{13}}{\partial\vec{r}_1}\\ \nonumber m_2 \ddot{r}_2 & = & -\frac{\partial U_{12}}{\partial\vec{r}_2}-\frac{\partial U_{23}}{\partial\vec{r}_2}\\ \nonumber m_3 \ddot{r}_3 & = & -\frac{\partial U_{13}}{\partial\vec{r}_3}-\frac{\partial U_{23}}{\partial\vec{r}_3}\\ \end{eqnarray}Si reconocemos que $\partial U_{23}/\partial \vec{r}_1=\partial U_{13}/\partial \vec{r}_2=\partial U_{12}/\partial \vec{r}_3$, las e.d.m. se pueden escribir de forma más completa como:
\begin{eqnarray} \nonumber m_1 \ddot{r}_1 & = & -\frac{\partial U_{12}}{\partial\vec{r}_1}-\frac{\partial U_{13}}{\partial\vec{r}_1}-\frac{\partial U_{23}}{\partial\vec{r}_1}\\ \nonumber m_2 \ddot{r}_2 & = & -\frac{\partial U_{12}}{\partial\vec{r}_2}-\frac{\partial U_{13}}{\partial\vec{r}_2}-\frac{\partial U_{23}}{\partial\vec{r}_2}\\ \nonumber m_3 \ddot{r}_3 & = & -\frac{\partial U_{12}}{\partial\vec{r}_3}-\frac{\partial U_{13}}{\partial\vec{r}_3}-\frac{\partial U_{23}}{\partial\vec{r}_3}\\ \end{eqnarray}Allí identificamos, finalmente, una sola función de energía potencial para todo el sistema:
\begin{equation} \label{eq:potencial_tres_cuerpos} U=U_{12}+U_{13}+U_{23}, \end{equation}de la cual puede calcularse la fuerza sobre cualquier partícula del sistema usando la prescripción general:
\begin{equation} \label{eq:ncuerpos_edm_potencial} \left\{m_i\ddot{\vec{r}}_i=-\frac{\partial U}{\partial\vec{r}_i}\right\}_N \end{equation}Esta es la forma más compacta de escribir las ecuaciones de movimiento del problema de los N cuerpos.
Estudie cuidadosamente las diferencia entre esta forma de las ecuaciones, que usa una sola función de energía potencial y la menos económica forma original escrita en las Ecs. (ncuerpos_edm_potencial_individual).
NOTA: Motivación física de la energía potencial de N cuerpos. La deducción que hicimos aquí de la energía potencial de tres cuerpos, es de naturaleza eminentemente matemática. Sin embargo, una vez obtenida, es posible dar a ella una interpretación que nos permite entender las razones físicas profundas por las que la energía potencial del sistema tiene en realidad un número restringido de términos.
Es claro, físicamente, que al cuantificar la energía total de interacción de un sistema, solo debemos tener en cuenta una vez la contribución que cada pareja de partículas. Así, el sistema de partículas 1 y 2 contribuyen con $U_{12}$ a la energía total, es obvio entonces que sumar a esa energía el término $U_{21}$ sería redundante y físicamente incorrecto. Lo mismo pasa con las energías de interacción $U_{13}$ que descartan la necesidad de incluir la $U_{31}$ y la $U_{23}$ que una vez sumada hace innecesario agregar la energía $U_{32}$.
Este razonamiento, es la razón física por la cual la funcuón de energía potencial en la Ec. (potencial_tres_cuerpos) total solo contiene tres términos (en lugar de 6 como podría intuirse de la forma de la Ec. ncuerpos_edm_potencial_individual).
¿Podría el lector generalizar este resultado a cuatro cuerpos? O más interesante aún ¿podría ofrecer una fórmula general para el número de términos que tiene la función de energía potencial de un sistema de N cuerpos? (intente encontrar la respuesta antes de mirar el pie de nota) $^1$$. Así, por ejemplo, con $N=3$ y $k=2$, el número de combinaciones es $C_{3,2}=3!/2!=3$ que coincide con el ejemplo desarrollado aquí (ver Ec. potencial_tres_cuerpos.)]
La energía potencial de N cuerpos $U$, que escribimos en la Ec. (potencial_tres_cuerpos) para un sistema de tres partículas, puede ahora escribirse, en general, para un sistema de N cuerpos, de dos maneras diferentes:
Restringiendo las sumatorias:
\begin{equation} \label{eq:U_restringido} U=-\sum_{i=1}^{i<j}\sum_{j=1}^{j=N} \frac{m_i\mu_j}{r_{ij}} \end{equation}
Una fórmula que se escribe en casi todos los textos y de forma simplificada como:
$$ U=-\sum_{i<j} \frac{m_i\mu_j}{r_{ij}} $$
Esta expresión es relativamente simple, pero expandir las sumatorias, en la práctica, puede ser muy complicado (el lector puede intentar, por ejemplo, escribir la energía potencial de un sistema de cuatro cuerpos, usando la fórmula anterior y no el razonamiento físico que esbozamos en la última Nota.)
Sin restringir las sumatorias:
\begin{equation} \label{eq:U_no_restringido} U=-\frac{1}{2}\sum_i\sum_{j\neq i} \frac{m_i\mu_j}{r_{ij}} \end{equation}
El factor $1/2$ en el lado derecho, aunque extraño desde el punto de vista de la física, proviene del hecho de que al no restringir la sumatoria sobre la $i$, apareceran los términos duplicados $m_i\mu_j/r_{ij}$ y $m_j\mu_i/r_{ji}$ que son exactamente iguales.
En lo que queda de este texto usaremos la parametrización de la Ec. (U_no_restringido).
Si multiplicamos las Ecs. (ncuerpos_edm_potencial) por el factor integrante $m_i\vec{r}_i$ (usando el producto punto o producto escalar) y sumamos sobre todas las partículas (como lo hemos hecho para obtener las otras constantes de movimiento) obtenemos:
$$ \sum_i m_i \dot{\vec{r}}_i\cdot\ddot{\vec{r}}_i = -\sum_i\dot{\vec{r}}_i\cdot\frac{\partial U}{\partial\vec{r}_i} $$El lado derecho de esta ecuación corresponde a la regla de la cadena (ver Ec. regla_cadena_general):
$$ \frac{\mathrm{d}U}{\mathrm{d}t}=\frac{\partial U}{\partial \vec{r}_1}\cdot\frac{\mathrm{d}\vec{r}_1}{\mathrm{d}t}+\frac{\partial U}{\partial \vec{r}_2}\cdot\frac{\mathrm{d}\vec{r}_2}{\mathrm{d}t}+\ldots+\frac{\partial U}{\partial \vec{r}_n}\cdot\frac{\mathrm{d}\vec{r}_n}{\mathrm{d}t} $$lo que permite conseguir el objetivo del factor integrante: escribir la ecuación diferencial en cuadraturas:
$$ \frac{d}{dt}\left(\frac{1}{2}\sum_i m_i \dot{\vec r_i}^2\right)=-\frac{dU}{dt} $$De esta ecuación podemos identificar, finalmente, una nueva constante de movimiento en el problema de los N cuerpos:
\begin{equation} \label{eq:ncuerpos_energia} \frac{1}{2}\sum_i m_i \dot{\vec r_i}^2+U(\{\vec{r}_i\})=E \end{equation}Con $U(\{\vec{r}_i\})=(1/2)\sum_i\sum_{j\neq i} m_i\mu_j/r_{ij}$.
Físicamente el valor de $E$ corresponde a la energía mecánica total del sistema de N cuerpos y esta integral no es otra cosa que una expresión del teorema de convervación de la energúa (ver Teo. regla_cadena_general) en la ausencia de fuerzas no conservativas.
En un terreno más práctico, sería interesante, usar este resultado, para preguntarse, por ejemplo, cuánta energía mecánica existe en el sistema Tierra-Luna.
Para ello podemos apoyarnos otra vez de SPICE, como lo hicimos antes para calcular el momentum angular del sistema Sol-Jupiter (ver Algoritmo masas_estado_planetas):
import spiceypy as spy
#Constante de gravitación universal
G=6.67e-20 # km^3 / kg s^2
#Asumimos un tiempo cualquiera, en este caso t=J2000.0
tef=0
#Carga kernels con posiciones (bsp) y masas (tpc)
spy.furnsh('pymcel/data/de430.bsp')
spy.furnsh('pymcel/data/de430.tpc')
#Parámetro gravitacional, posiciones y velocidades
mutierra=spy.bodvrd("EARTH","GM",1)[1][0]
tierra,tluz=spy.spkezr("EARTH",tef,
"ECLIPJ2000","None","EARTH_BARYCENTER")
rtierra=tierra[:3]
vtierra=tierra[3:]
muluna=spy.bodvrd("MOON","GM",1)[1][0]
luna,tluz=spy.spkezr("MOON",tef,
"ECLIPJ2000","None","EARTH_BARYCENTER")
rluna=luna[:3]
vluna=luna[3:]
Una diferencia importante entre este algoritmo y el que vimos para el cálculo del momentum angular, es que las posiciones y velocidades de la Tierra y la Luna, son calculadas respecto a la posición del centro de masa del sistema (EARTH_BARYCENTER) en lugar de hacerlo respecto del baricentro del sistema solar (SSB). Si bien, nada en la Ec. (ncuerpos_energia) nos dice que la energía mecánica total solo se conserva en el sistema de referencia del centro de masa, nuestra interpretación del valor de esta cantidad será más fácil allí. El lector puede intentar hacer el cálculo relativo al baricentro del sistema solar y comparar con el resultado obtenido aquí.
La energía mecánica total del sistema, será entonces:
#Masas
mtierra=mutierra/G
mluna=muluna/G
#Energía potencial
from numpy.linalg import norm
U=-G*mluna*mtierra/norm(rtierra-rluna)
#Energía cinética
K=0.5*mutierra*norm(vtierra)**2+0.5*mluna*norm(vluna)**2
#Energía mecánica total
E=K+U
print(f"Energía potencial U: {U} kg km^2/s^2")
print(f"Energía cinética K: {K} kg km^2/s^2")
print(f"Energía mecánica total E: {E} kg km^2/s^2")
¿Qué significa que la energía mecánica total del sistema Tierra-Luna sea negativa?
Como lo comentamos en la Sección Energía potencial gravitacional el signo negativo de la energía potencial gravitacional (responsable aquí del signo negativo de $E$) se debe a que escogimos el valor de la constante $C$ en la Ec. (energia_potencial_gravitacional) como cero.
Siguiendo un razonamiento similar al que usamos en esa sección para explicar físicamente el significado de esta elección arbitraria (y por tanto el signo de $U$), podemos decir aquí que el valor de $|E|$ en el sistema Tierra-Luna y en general en cualquier sistema de dos cuerpos, es el trabajo mínimo requerido para separar a la Tierra y a la Luna a una distancia inmensa. O, en términos escatológicos, $|E|$ es la energía mínima necesaria para destruir la unión entre nuestro planeta y su eterna compañera.
No es difícil mostrar (ver problemas al final del capítulo) que esta energía es igual a la energía cinética de un cuerpo con la masa de la Luna moviéndose a $~$1 km/s relativa al centro de masa del sistema. Es decir, una colisión elástica con un cuerpo gemelo de la Luna podría hacer el trabajo de separarnos de ella$^2$.
Para hacernos ahora a una idea más cercana del valor de esta energía, las unidades de energía utilizadas aquí, kg km$^2$/s$^2$, que siguen la convención en SPICE de usar km como unidad de longitud, no son las más convenientes. Para convertirlas a las unidades del sistema internacional, los joules, J = kg m$^2$/s$^2$, basta multiplicar las energía obtenidas anteriormente por un factor de conversión de $10^6$.
En conclusión podemos decir que (el valor absoluto de) la energía mecánica total del sistema Tierra-Luna es $~4\times 10^{28}$ j. Esta energía es aproximadamente igual a la energía consumida por toda la humanidad en casi 2 millones de años.
NOTAS AL PIE: