Ondas gravitacionales: del cuadrupolo a GW150914

1. Motivación: escuchar el espacio-tiempo

El 14 de septiembre de 2015, a las 09:50:45 UTC, los dos interferómetros de LIGO registraron, con 6.9 ms de diferencia, una señal que subía en frecuencia de 35 a 250 Hz en 0.2 s y alcanzaba una deformación de $10^{-21}$: un cambio de longitud de una milésima del diámetro de un protón en brazos de 4 km. Era la onda gravitacional emitida por dos agujeros negros de 36 y 29 masas solares al fundirse en uno de 62, a 410 Mpc; las tres masas solares restantes se habían convertido en energía gravitacional radiada, con una potencia máxima superior a la luz de todas las estrellas del universo observable. GW150914 confirmó, cien años después de que Einstein las predijera (y dudara de ellas), que las ondas gravitacionales existen y se pueden medir.

Pero la evidencia indirecta llevaba cuarenta años acumulándose: el púlsar binario PSR B1913+16, descubierto por Hulse y Taylor en 1974, pierde energía orbital exactamente al ritmo que dicta la fórmula del cuadrupolo. En esta lección derivamos esa fórmula y la confrontamos con ambos sistemas.

2. Gravedad linealizada y la ecuación de ondas

2.1 Perturbaciones de Minkowski

Escribimos $g_{\mu\nu} = \eta_{\mu\nu} + h_{\mu\nu}$ con $|h_{\mu\nu}|\ll1$ y conservamos términos lineales en $h$. Los símbolos de Christoffel son $\Gamma^\mu_{\nu\rho} = \tfrac12\eta^{\mu\sigma}(\partial_\nu h_{\sigma\rho} + \partial_\rho h_{\sigma\nu} - \partial_\sigma h_{\nu\rho})$ y el tensor de Einstein linealizado, en términos de la traza invertida $\bar h_{\mu\nu}\equiv h_{\mu\nu}-\tfrac12\eta_{\mu\nu}h$ (con $h = \eta^{\mu\nu}h_{\mu\nu}$), es

con $\Box = \eta^{\mu\nu}\partial_\mu\partial_\nu = -\partial_t^2/c^2 + \nabla^2$.

2.2 Libertad de gauge y gauge de Lorenz

Un cambio infinitesimal de coordenadas $x^\mu\to x^\mu+\xi^\mu$ transforma $h_{\mu\nu}\to h_{\mu\nu}-\partial_\mu\xi_\nu-\partial_\nu\xi_\mu$ sin alterar la física. Eligiendo $\xi^\mu$ tal que $\Box\xi_\nu = \partial^\mu\bar h_{\mu\nu}$ se impone el gauge de Lorenz $\partial^\mu\bar h_{\mu\nu} = 0$, que anula los tres últimos términos de $G_{\mu\nu}$. Las ecuaciones de Einstein $G_{\mu\nu} = 8\pi GT_{\mu\nu}/c^4$ se reducen a

una ecuación de ondas con fuente: la gravedad se propaga a la velocidad $c$. En el vacío, $\Box\bar h_{\mu\nu} = 0$ admite ondas planas $\bar h_{\mu\nu} = A_{\mu\nu}e^{ik_\rho x^\rho}$ con $k_\rho k^\rho = 0$.

2.3 Gauge transverso y sin traza (TT)

Con la libertad de gauge residual ($\Box\xi^\mu = 0$) se pueden imponer además $\bar h = 0$, $\bar h_{0\mu} = 0$ y $\partial_ih_{ij} = 0$: el gauge TT. Para una onda que viaja en $z$ quedan solo dos componentes independientes, las polarizaciones

El efecto físico sobre dos masas libres separadas $L$ en el plano transversal es $\delta L/L = \tfrac12h$: una onda $h_+$ estira en $x$ y comprime en $y$ (y al medio periodo, al revés); $h_\times$ hace lo mismo girado $45°$. Un interferómetro de Michelson con brazos perpendiculares mide directamente la diferencia $\delta L_x - \delta L_y = hL$: con $L = 4$ km y $h = 10^{-21}$, $4\times10^{-18}$ m.

3. La fórmula del cuadrupolo

3.1 Solución retardada

La solución de la ecuación de ondas con fuente es la de la electrodinámica:

En la zona lejana ($r\gg$ tamaño de la fuente $\gg$ y además $r\gg\lambda$) aproximamos $|\mathbf x-\mathbf x'|\approx r$ en el denominador y, si la fuente es lenta ($v\ll c$, de modo que $\lambda\gg$ tamaño), también en el argumento retardado:

3.2 De $T_{ij}$ al momento cuadrupolar

Usando la conservación $\partial_\mu T^{\mu\nu} = 0$ dos veces (identidad de Laplace, problema 1):

con $\rho = T^{00}/c^2$ la densidad de masa e $I^{ij}$ el segundo momento de masa. Así,

Es la fórmula del cuadrupolo (Einstein 1918). Proyectando al gauge TT solo sobrevive la parte sin traza, el momento cuadrupolar reducido $Q_{ij} = I_{ij}-\tfrac13\delta_{ij}I_{kk}$. No hay radiación dipolar: el dipolo de masa es $\int\rho\,\mathbf x$, cuya segunda derivada es la derivada del momento total, nula por conservación. Por eso la radiación gravitacional es intrínsecamente más débil que la electromagnética y por eso una esfera que pulsa radialmente o un cilindro que gira sobre su eje no radian.

3.3 Potencia radiada

El flujo de energía de la onda es $t^{00} = \frac{c^4}{32\pi G}\langle\dot h_{ij}^{\rm TT}\dot h_{ij}^{\rm TT}\rangle$ (pseudotensor de Isaacson, promediado sobre varios periodos). Integrando sobre la esfera y proyectando sobre todas las direcciones (problema 2) resulta la fórmula de luminosidad del cuadrupolo:

El prefactor $G/c^5 = 2.76\times10^{-53}$ W$^{-1}$ explica la debilidad: solo $\dddot Q\sim10^{26}$ en unidades SI (masas estelares a velocidades relativistas) produce potencias apreciables. Su inverso, $c^5/G = 3.63\times10^{52}$ W, es la luminosidad de Planck, cota natural de cualquier fuente; GW150914 alcanzó una milésima de ella.

4. Binarias: Peters y Mathews

4.1 Órbita circular

Dos masas $m_1$, $m_2$ ($M = m_1+m_2$, masa reducida $\mu = m_1m_2/M$) en órbita circular de radio $a$ y frecuencia orbital $\Omega = \sqrt{GM/a^3}$ equivalen a una masa $\mu$ en $\mathbf x = a(\cos\Omega t,\sin\Omega t,0)$. El cuadrupolo reducido tiene componentes $Q_{xx} = \mu a^2(\cos^2\Omega t-\tfrac13)$, $Q_{yy} = \mu a^2(\sin^2\Omega t-\tfrac13)$, $Q_{xy} = \mu a^2\sin\Omega t\cos\Omega t$, que oscilan a frecuencia $2\Omega$ (de ahí $f_{\rm GW} = 2f_{\rm orb}$). Derivando tres veces y promediando, $\langle\dddot Q_{ij}\dddot Q_{ij}\rangle = 32\mu^2a^4\Omega^6$ y

4.2 Órbita excéntrica

Peters y Mathews (1963) promediaron sobre una órbita kepleriana de excentricidad $e$ y obtuvieron el mismo resultado multiplicado por el factor de realce

que vale 1 en órbita circular y 11.86 para $e = 0.617$ (Hulse–Taylor): la emisión se concentra en el periastro. Peters (1964) dedujo de ahí la evolución secular de $a$ y $e$:

Con la tercera ley de Kepler, $\dot P_b/P_b = \tfrac32\dot a/a$, se obtiene la fórmula de Peters para el periodo orbital:

Todo esto está implementado en src/lib/math/ondas-gravitacionales.ts.

5. El púlsar de Hulse–Taylor

5.1 El sistema

PSR B1913+16, descubierto en Arecibo en julio de 1974 (Hulse & Taylor 1975), es un púlsar de 59 ms en órbita de 7.75 h con otra estrella de neutrones. El seguimiento de los tiempos de llegada de sus pulsos durante más de cuarenta años (9257 medidas; Weisberg & Huang 2016) da los parámetros keplerianos con precisión extraordinaria y tres efectos relativistas medibles: el avance del periastro $\dot\omega$, la dilatación temporal $\gamma$ y el decaimiento del periodo $\dot P_b$. Los dos primeros fijan las masas; el tercero es la prueba.

5.2 Nuestro cálculo

Con derivadaPeriodoOrbital(m1, m2, Pb, e) y las constantes IAU/CODATA de constants.ts: $a = 1.949\times10^9$ m $= 2.80\,R_\odot$, $f(e) = 11.857$ y

que reproduce el valor publicado con error relativo $2\times10^{-4}$ (diferencia atribuible al redondeo de las masas). El test ondas-gravitacionales.test.ts lo comprueba. En cuatro décadas, el periastro llega adelantado más de un minuto respecto a una órbita de periodo constante: es la famosa parábola de Weisberg y Taylor. La potencia radiada es $P = 7.8\times10^{24}$ W (2 % de la luminosidad solar), el semieje se contrae 3.5 m al año, y la integración de las ecuaciones de Peters da una coalescencia dentro de unos 300 millones de años.

5.3 Análisis de error

La razón $0.9983\pm0.0016$ está a $1\sigma$ de la unidad. Su incertidumbre no la fija el púlsar (el error de $\dot P_b^{\rm obs}$ es $4\times10^{-4}$ relativo), sino la corrección $\dot P_b^{\rm gal}$ por las aceleraciones relativas del púlsar y del Sol en la Galaxia (efecto Shklovskii y gradiente del potencial), que depende de la distancia al púlsar, mal conocida ($\sim$ 20 %), y de $R_0$, $\Theta_0$ (Weisberg & Huang usan $8.34$ kpc y $240$ km/s). La lección de curvas de rotación mostró cuánto varían esos parámetros entre autores: aquí cambian el resultado en la tercera cifra. Un paralaje VLBA del púlsar reduciría la incertidumbre.

6. Chirp: masa de chirp y evolución de la frecuencia

6.1 La masa de chirp

Para una binaria circular, la energía orbital es $E = -G\mu M/2a$ y $f_{\rm GW} = \Omega/\pi$. Igualando $dE/dt = -P$ y expresando $a$ en función de $f$:

La masa de chirp $\mathcal M$ es la única combinación de masas que determina, a este orden, cómo sube la frecuencia; medir $f$ y $\dot f$ en la señal la da directamente: $\mathcal M = \frac{c^3}{G}\left[\frac{5}{96}\pi^{-8/3}f^{-11/3}\dot f\right]^{3/5}$ (masaChirpDesdeChirp). Integrando $\dot f$, el tiempo hasta la coalescencia desde la frecuencia $f$ es

y la amplitud en orientación óptima, de la fórmula del cuadrupolo,

6.2 GW150914

La señal (Abbott et al. 2016, PRL 116, 061102) fue detectada con razón señal/ruido 24, llegó primero a Livingston y $6.9^{+0.5}_{-0.4}$ ms después a Hanford, barrió de 35 a 250 Hz en ~0.2 s (unos 8 ciclos) con deformación máxima $1.0\times10^{-21}$, y coincide con la forma de onda de relatividad numérica para la inspiración, fusión y relajación («ringdown») de dos agujeros negros. Parámetros (intervalos de credibilidad del 90 %, marco de la fuente):

Del par $(f,\dot f)$ en la fase de inspiración los autores obtienen una masa de chirp en el marco del detector $\mathcal M\simeq30\,M_\odot$, lo que implica $M\gtrsim70\,M_\odot$ en ese marco. Con las masas de la tabla, masaChirp(36, 29) da $28.1\,M_\odot$ en la fuente, es decir $(1+z)\times28.1 = 30.6\,M_\odot$ en el detector: consistente. (Las masas observadas aparecen multiplicadas por $1+z$ porque la frecuencia se dilata cosmológicamente.)

Por qué agujeros negros. A $f_{\rm GW} = 150$ Hz, la tercera ley de Kepler con $M = 70\,M_\odot$ da una separación $a = [GM/(\pi f)^2]^{1/3} = 347$ km, mientras que el radio de Schwarzschild de $70\,M_\odot$ es 207 km y la velocidad orbital, $0.43c$. Solo dos objetos más compactos que estrellas de neutrones (que a esas masas no existen) pueden orbitar tan cerca sin tocarse: agujeros negros. Es el argumento del artículo original.

Órdenes de magnitud. Con $\mathcal M = 30\,M_\odot$ y $D = 410$ Mpc, amplitudDeformacion da $h = 2.4\times10^{-21}$ a 150 Hz en orientación óptima; promediando sobre orientaciones y la respuesta del detector se reduce a $\sim10^{-21}$, lo observado. El tiempo desde 35 Hz hasta la fusión es $\tau = 0.17$ s y $\dot f(35\,\text{Hz}) = 77$ Hz/s. La energía radiada, $3\,M_\odot c^2 = 5.4\times10^{47}$ J, es la que el Sol emite en $4\times10^{13}$ años; la luminosidad máxima, $3.6\times10^{49}$ W, es $10^{-3}c^5/G$ y supera a toda la luz del universo observable ($\sim10^{49}$ W).

El simulador integra la fase newtoniana del chirp hasta $f_{\rm ISCO}$. Con el preset GW150914 comprueba $\tau(35\,\text{Hz})$ y el orden de $h$; el pie usa la masa de chirp en el marco del detector, $(1+z)\mathcal M$, por eso su $\dot f$ sale algo mayor que el valor de arriba. Cambia a GW170817 y explica por qué la señal dura minutos en lugar de una fracción de segundo (ec. de $\tau \propto \mathcal M^{-5/3} f^{-8/3}$).

7. Problemas

1. Identidad del cuadrupolo (demostración). A partir de $\partial_\mu T^{\mu\nu} = 0$, demuestra que $\int T^{ij}d^3x = \frac{1}{2c^2}\frac{d^2}{dt^2}\int T^{00}x^ix^j\,d^3x$ para una fuente de soporte compacto.

2. El factor 1/5 (demostración). Partiendo del flujo $\frac{c^3}{32\pi G}\langle\dot h^{\rm TT}_{ij}\dot h^{\rm TT}_{ij}\rangle r^2d\Omega$, usa el proyector TT $\Lambda_{ij,kl}(\hat n)$ y las integrales angulares $\int\hat n_i\hat n_j\,d\Omega = \frac{4\pi}{3}\delta_{ij}$, $\int\hat n_i\hat n_j\hat n_k\hat n_l\,d\Omega = \frac{4\pi}{15}(\delta_{ij}\delta_{kl}+\delta_{ik}\delta_{jl}+\delta_{il}\delta_{jk})$ para obtener $P = \frac{G}{5c^5}\langle\dddot Q_{ij}\dddot Q_{ij}\rangle$.

3. Potencia de la binaria circular. Calcula $\dddot Q_{ij}$ para la binaria de §4.1 y verifica $\langle\dddot Q_{ij}\dddot Q_{ij}\rangle = 32\mu^2a^4\Omega^6$.

4. Masa de chirp (demostración). Deduce $\dot f$ de §6.1 a partir de $E = -G\mu M/2a$, $P$ de §4.1 y la tercera ley de Kepler.

5. Hulse–Taylor en Python (computacional). Implementa la fórmula de Peters y las ecuaciones de evolución de $a$ y $e$. Reproduce $\dot P_b^{\rm RG}$ e integra hasta la coalescencia; dibuja $e(t)$ y comprueba que la órbita se circulariza antes de fundirse.

import numpy as np
from scipy.integrate import solve_ivp
G=6.6743e-11; c=299792458.0; Ms=1.98841e30; yr=365.25*86400
m1,m2=1.438*Ms,1.390*Ms; M=m1+m2; Pb=0.322997448918*86400; e0=0.6171340
f=lambda e:(1+73/24*e**2+37/96*e**4)/(1-e**2)**3.5
Pbdot=-192*np.pi/5*(2*np.pi*G/(c**3*Pb))**(5/3)*m1*m2/M**(1/3)*f(e0)
print(Pbdot)   # -2.402e-12
a0=(G*M*Pb**2/(4*np.pi**2))**(1/3)
def rhs(t,y):
    a,e=y; k=G**3*m1*m2*M/c**5
    return [-64/5*k/a**3*f(e), -304/15*k*e/a**4*(1+121/304*e**2)/(1-e**2)**2.5]
sol=solve_ivp(rhs,[0,1e9*yr],[a0,e0],events=lambda t,y:y[0]-1e6,rtol=1e-8)
print(sol.t[-1]/yr/1e6, "Myr", sol.y[1,-1])

6. GW150914 a mano. De los datos de la PRL: (a) calcula $\mathcal M$ para $(36, 29)\,M_\odot$ y pásala al marco del detector; (b) con $\mathcal M_{\rm det} = 30\,M_\odot$ estima $\dot f$ a 35 Hz y el tiempo hasta la fusión; (c) la separación y la velocidad orbital a 150 Hz para $M = 70\,M_\odot$; (d) la amplitud $h$ a 410 Mpc y 150 Hz.

7. Una binaria que no se detecta. Dos estrellas de $1\,M_\odot$ en órbita circular de 1 día a 100 pc. Calcula $f_{\rm GW}$, $P$, $h$ y el tiempo de coalescencia. ¿Por qué LIGO no la ve y sí la vería LISA (banda $10^{-4}$–$1$ Hz)?

8. Autoevaluación

9. Lecturas

• Schutz (2009), cap. 9: la derivación más accesible de la gravedad linealizada y el cuadrupolo. Carroll (2004), cap. 7. Hartle (2003), caps. 16 (ondas) y 23 (fuentes). Weinberg (1972), cap. 10.
• Maggiore (2008), vol. 1, cap. 3 (fórmula del cuadrupolo con todos los detalles) y cap. 4 (binarias, Peters–Mathews, chirp).
• Peters & Mathews (1963), Phys. Rev. 131, 435; Peters (1964), Phys. Rev. 136, B1224.
• Hulse & Taylor (1975), ApJ 195, L51; Taylor & Weisberg (1982), ApJ 253, 908; Weisberg & Huang (2016), ApJ 829, 55, arXiv:1606.02744.
• Abbott et al. (2016), PRL 116, 061102 (acceso abierto), y PRL 116, 241102 para los parámetros. Datos públicos en el Gravitational Wave Open Science Center (gwosc.org): descarga GW150914 y reproduce la figura 1.