Gravitational Radiation from Close Binaries with Time-Varying Masses

TL;DR

This paper extends the understanding of gravitational waves from close binaries by including effects of time-varying masses, showing these can significantly amplify gravitational radiation beyond orbital contributions.

Contribution

It introduces phase-dependent formulas for gravitational radiation from binaries with changing masses, expanding upon the classical Peters and Mathews model.

Findings

Mass variation effects can dominate gravitational wave emission.

Derived formulas reduce to classical results for constant masses.

Mass changes can increase gravitational radiation by orders of magnitude.

Abstract

In the quadrupole approximation of General Relativity in the weak-field limit, a time-varying quadrupole moment generates gravitational radiation. Binary orbits are one of the main mechanisms for producing gravitational waves and are the main sources and backgrounds for gravitational-wave detectors across the multi-band spectrum. In this Paper, we introduce additional contributions to the gravitational radiation from close binaries that arise from time-varying masses along with those produced by orbital motion. We derive phase-dependent formulae for these effects in the quadrupolar limit for binary point masses, which reduce to the formulae that Peters and Mathews (1963) derived when the mass of each component is taken to be constant. We show that gravitational radiation from mass variation can be orders of magnitude greater than that of orbital motion.

Equations

$$\ddot{\mathbf{r}}=-\frac{GM}{r^3}\mathbf{r},$$

$$\mathbf{f}=-\frac{1}{2}\frac{\dot{M}}{M}\dot{\mathbf{r}},$$

$$\frac{\mathrm{d}a}{\mathrm{d}t}$$

$$\frac{\mathrm{d}e}{\mathrm{d}t}$$

$$\frac{\mathrm{d}\omega }{\mathrm{d}t}$$

$$\frac{\mathrm{d}\nu }{\mathrm{d}t}$$

$$r=\frac{a\left(1-e^2\right)}{1+e\cos \nu }.$$

$$I_{xx}$$

$$I_{yy}$$

$$I_{xy}$$

$${\stackrel{...}{I}}_{xx}$$

$$+3\ddot{\mu }{r}^2\dot{\theta }\left(\frac{2e\sin \nu {\cos }^2\theta }{1+e\cos \nu }-\sin 2\theta \right)+\stackrel{...}{\mu }{r}^2{\cos }^2\theta ,$$

$${\stackrel{...}{I}}_{yy}$$

$$+3\ddot{\mu }{r}^2\dot{\theta }\left(\frac{2e\sin \nu {\sin }^2\theta }{1+e\cos \nu }-\sin 2\theta \right)+\stackrel{...}{\mu }{r}^2{\sin }^2\theta ,$$

$${\stackrel{...}{I}}_{xy}$$

$$+3\ddot{\mu }{r}^2\dot{\theta }\left(\frac{e\sin \nu \sin 2\theta }{1+e\cos \nu }+\cos 2\theta \right)+\frac{1}{2}\stackrel{...}{\mu }{r}^2\sin 2\theta .$$

$${\frac{\mathrm{d}E}{\mathrm{d}t}|}_{\mathrm{GW}}$$

$${\frac{\mathrm{d}{J}_i}{\mathrm{d}t}|}_{\mathrm{GW}}$$

$$\begin{gathered} L_{\mathrm{GW}}=\frac{32}{5}\frac{G^4}{c^5}\frac{{\mu }^2{M}^3}{a^5}+\frac{72}{5}\frac{G^3}{c^5}\frac{{\dot{\mu }}^2{M}^2}{a^2}-\frac{48}{5}\frac{G^3}{c^5}\frac{\ddot{\mu }\mu M^2}{a^2} \\ -\frac{12}{5}\frac{G^2}{c^5}\stackrel{...}{\mu }\dot{\mu }Ma+\frac{18}{5}\frac{G^2}{c^5}{\ddot{\mu }}^{2}Ma+\frac{1}{5}\frac{G}{c^5}{\stackrel{...}{\mu }}^2{a}^4, \end{gathered}$$

$$h_{jk}^{\mathrm{TT}}=\frac{2G}{c^4{D}_{\mathrm{L}}}{\ddot{\mathcal{I}}}_{jk},$$

$$h_0=\frac{4G^{5/3}}{c^4}\frac{\mu M^{2/3}{\Omega }^{2/3}}{D_{\mathrm{L}}}+\frac{2G^{5/3}}{c^4}\frac{\dot{\mu }{M}^{2/3}}{D_{\mathrm{L}}{\Omega }^{1/3}}+\frac{G^{5/3}}{c^4}\frac{\ddot{\mu }{M}^{2/3}}{D_{\mathrm{L}}{\Omega }^{4/3}},$$

$${\zeta }_1\equiv \frac{\dot{\mu }/\mu }{2\Omega },$$

$${\zeta }_2\equiv \frac{\ddot{\mu }/\mu }{{\left(2\Omega \right)}^2},$$

$$I_{jk}=\left(\begin{array}{ccc}\mu x^2& \mu xy& \mu xz\\ \mu xy& \mu y^2& \mu yz\\ \mu xz& \mu yz& \mu z^2\end{array}\right),$$

$${\stackrel{...}{I}}_{xx}$$

$${\stackrel{...}{I}}_{yy}$$

$${\stackrel{...}{I}}_{zz}$$

$${\stackrel{...}{I}}_{xy}$$

$$+6\dot{\mu }\dot{x}\dot{y}+3\mu \ddot{x}\dot{y}+3\dot{\mu }x\ddot{y}+3\mu \dot{x}\ddot{y}+\mu x\stackrel{...}{y},$$

$${\stackrel{...}{I}}_{xz}$$

$$+6\dot{\mu }\dot{x}\dot{z}+3\mu \ddot{x}\dot{z}+3\dot{\mu }x\ddot{z}+3\mu \dot{x}\ddot{z}+\mu x\stackrel{...}{z},$$

$${\stackrel{...}{I}}_{yz}$$

$$+6\dot{\mu }\dot{y}\dot{z}+3\mu \ddot{y}\dot{z}+3\dot{\mu }y\ddot{z}+3\mu \dot{y}\ddot{z}+\mu y\stackrel{...}{z}.$$

$$m_2(t)=m_{2,\mathrm{f}}+\frac{1}{2}\left(m_{2,\mathrm{i}}-m_{2,\mathrm{f}}\right)\left(1+\frac{2}{\pi }{\tan }^{-1}\left[(t-t_0)/\tau \right]\right),$$

Full text

gravitational waves — close binaries

††software: Matplotlib (Hunter, 2007), Numpy (Walt et al., 2011)

1 Introduction

Einstein’s theory of General Relativity (GR) predicted the existence of gravitational waves (GWs) (Einstein, 1916, 1918). The existence of GWs was later indirectly confirmed from observations of the orbital decay of a pulsar binary (Hulse & Taylor, 1975). Recently, the Laser Interferometer Gravitational Wave Observatory (LIGO) has ushered in the era of GW astronomy with detections of GWs from compact-object mergers (LIGO-Virgo Collaboration, 2016a, b, 2017a, 2017b, 2017c, 2017d). The Laser Interferometer Space Antenna (LISA) will then open up the millihertz frequency band (Danzmann et al., 1996; Amaro-Seoane et al., 2017). The combination of ground-based and space-based detectors holds promise for advances in multi-band GW astronomy, standard-siren cosmology, and tests of GR (Sesana, 2016). At nanohertz frequencies, pulsar timing arrays are currently operating and will eventually detect the stochastic GW background from the cosmic population of inspiraling supermassive black-hole binaries (Arzoumanian et al., 2018). Even though a detection has yet to be made, the upper limits are providing astrophysically relevant constraints on the cosmic supermassive black-hole binary population (Sesana et al., 2018; Holgado et al., 2018b). The decihertz band may be opened with future space-based detectors. The main astrophysical GW sources and GW backgrounds that these detectors across the multi-band GW spectrum aim to detect are from the orbits of compact binaries. In this Paper, we consider an additional effect: time-varying masses in a binary due to wind mass loss or accretion from mass transfer, for example, generate time-varying quadrupole moments that generate gravitational radiation in addition to those produced by orbital motion.

The analytic expressions for gravitational radiation from point-mass binaries in the quadrupole approximation of GR in the weak-field limit (Peters & Mathews, 1963; Peters, 1964) have been used widely in the astrophysics community for modeling the evolution of close stellar binaries, generating compact-binary distributions from binary population synthesis codes (e.g., Belczynski et al., 2002, 2008), and estimating the strain from sub-parsec supermassive black-hole binary candidates. In this Paper, we consider the more general treatment of the two-body problem in which the mass of the binary is time-dependent. This is motivated by the variety of binary systems that interact, i.e., low-mass X-ray binaries, cataclysmic variables, and supernovae in close binaries. Binary interactions include Roche-lobe mass transfer (e.g., Paczyński, 1971; Eggleton, 1983), mass loss, and common-envelope evolution (e.g., Holgado et al., 2018a). A number of previous works in the astrophysics and fundamental-physics communities have considered close binaries with time-varying masses in different contexts (e.g., Macedo et al., 2013; Barausse et al., 2014; Tauris, 2018). Cheng et al. (2004) considered the effect of the time-varying masses on the quadrupole moment in the context of neutron-star spindown, but did not consider the phase-dependent effects that we present here. We derive new expressions for gravitational radiation from binaries with variable masses and we quantify the differences between treatments with and without the additional contributions for an explosive mass-loss case. The expressions that arise from this new analysis reduce to the Peters and Mathews (PM herein) formulae (Peters & Mathews, 1963; Peters, 1964) when the mass of each component is taken to be constant.

2 The Perturbed Two-Body Problem

If two constant point masses are in an orbit bound by their mutual gravitational interaction, their motion can be described as a Keplerian orbit in the absence of any internal and external perturbations. The six orbital elements $(a,e,i,\varpi,\omega,\nu)$ fully define the Keplerian orbit, where $a$ is the semi-major axis, $e$ is the eccentricity, $i$ is the inclination, $\varpi$ is the longitude of the ascending node, $\omega$ is the argument of periapsis, and $\nu$ is the true anomaly. We consider a gravitationally bound point-mass binary with time-dependent component masses $m_{1}(t)$ and $m_{2}(t)$ with a reduced mass $\mu=m_{1}m_{2}/M(t)\equiv\mu(t)$, where $M(t)=m_{1}+m_{2}$ is the total mass. Allowing for time-dependent masses introduces perturbative forces on the two bodies.

We consider the case where the primary’s mass remains constant and the secondary undergoes isotropic mass loss: $\dot{m_{1}}=0$ and $\dot{m}_{2}<0$. The resulting total mass-loss rate of the binary is $\dot{M}=\dot{m}_{2}$. For the isotropic case, only the self-gravity between the two components governs the equation of motion for the relative separation vector ${\bf r}$

[TABLE]

where $G$ is the gravitational constant, $r$ is the magnitude of ${\bf r}$, and where we take the same coordinate system as Peters & Mathews (1963) with the origin located at the barycenter. There is an implicit perturbing force that induces the following acceleration

[TABLE]

which is implicit in the isotropic variation of the mass-losing body. This has been derived in several previous works (Hadjidemetriou, 1963, 1966a, 1966b; Verhulst, 1975). If the mass variation is anisotropic, i.e., either from non-conservative mass exchange and asymmetric ejecta, then the equation of motion will have an explicit perturbing force (e.g., Dosopoulou & Kalogera, 2016a; Dosopoulou & Kalogera, 2016b). The orbital evolution can be described as an osculating orbit, which is the Keplerian orbit the binary would follow if the perturbation ceased immediately. The full set of equations for the time evolution of the osculating orbital elements from the isotropically varying mass are (e.g., Hadjidemetriou, 1963)

[TABLE]

where the specific angular momentum is conserved, i.e., the quantity $r^{2}\dot{\theta}=\left[GMa\left(1-e^{2}\right)\right]^{1/2}$ is a constant of the motion. Since there is no net force in the direction perpendicular to the orbital plane, $\dot{i}=\dot{\varpi}=0$, and the orbital motion thus remains coplanar. The mass variation causes the argument of periapsis to precess at a rate $\dot{\omega}$ (see Eq. 3c).

3 Quadrupole Gravitational Radiation

For eccentric Keplerian orbits, the orbital separation $r$ between the binary masses obeys

[TABLE]

The components of the mass quadrupole tensor $I_{jk}$ for an eccentric binary are

[TABLE]

where $\theta=\nu+\omega$, which simplifies to $\theta=\nu$ if the precession of the argument of periapsis is taken to be negligible. We can choose the initial $\theta$ to begin where $\omega=0$. For isotropic wind mass loss, the periapsis distance $r_{\rm p}=a(1-e)$ is constant, so ${\dot{r}}_{\rm p}=0$. When including gravitational radiation reaction, the periapsis distance will indeed decrease, though one can take ${\dot{r}}_{\rm p}\approx 0$ on timescales less than the mass-loss timescale $M/|\dot{M}|$ when the GW coalescence timescale is much longer than the mass-loss timescale. The third time derivatives of the mass-quadrupole tensor components are {widetext}

[TABLE]

If the periapsis precession is taken to be negligible, $\dot{\omega}=0$, then the first two terms of Eqs. (6) simplify to the canonical PM expressions. The additional terms in Eqs. (6) contain time derivatives of the masses of the binary, which may also be thought of as time derivatives of the binary’s chirp mass ${\cal M}=(m_{1}m_{2})^{3/5}/M^{1/5}$.

Gravitational radiation carries away the orbital energy and angular momentum of the binary at the following rates

[TABLE]

where $c$ is the speed of light and ${\cal I}_{jk}$ is the reduced (trace-free) quadrupole moment tensor ${\cal I}_{jk}=I_{jk}-\frac{1}{3}\delta_{jk}\delta_{\ell m}I_{\ell m}$. If we take the orbit to be circular, then the GW luminosity simplifies to

[TABLE]

where the positive terms in Eq. (8) correspond to GW energy loss from the orbit and the negative terms correspond to GW energy gain in the orbit. From the strain tensor in the transverse-traceless gauge,

[TABLE]

where $D_{\rm L}$ is the luminosity distance to the source, the strain amplitude for a circular binary is

[TABLE]

where $\Omega$ is the angular frequency of the binary, $\Omega=2\pi f_{\rm orb}=2\pi/P_{\rm orb}$. In the transverse-traceless gauge, the plus and cross polarizations of the strain are $h_{+}=h_{xx}^{\rm TT}=-h_{yy}^{\rm TT}$ and $h_{\times}=h_{xy}^{\rm TT}=h_{yx}^{\rm TT}$. From Eq. (10), we can define the following parameters

[TABLE]

which quantify the relative importance of mass-variation effects on the GW amplitude relative to the PM term of a given source. For $\zeta_{1}=1$, the mass-transfer timescale needs to be comparable to the orbital period. The true trajectory may differ significantly from that of an osculating orbit since the perturbing force in this regime is comparable to the self-gravity.

We emphasize that even though the above analysis is solely applicable to a binary whose mass is varying isotropically, it is illustrative of how one includes time-varying masses in the PM analysis and eventually extensions to the post-Newtonian analysis. In the limit that the time-derivatives of the masses are zero, then the canonical PM expressions remain.

4 Explosive Mass Loss in a Close Binary

Supernova explosions in close binaries result in extreme mass loss. If enough mass is ejected during the explosion, then the binary can be disrupted, i.e., the final energy of the binary is $\geq 0$ and the orbit becomes unbound. We again consider the case of isotropic mass loss such that the kick imparted to the remnant of the exploding body is solely due to the acceleration from the ejected mass (i.e., a “Blaauw kick” (Blaauw, 1961)), such that the equation of motion remains the same as Eq. (1). In this regime, we choose to work in Cartesian coordinates such that ${\bf r}=(x,y,z)$, $\dot{\bf r}=(\dot{x},\dot{y},\dot{z})$, and $\ddot{\bf r}=(\ddot{x},\ddot{y},\ddot{z})$. We also ignore any interaction between the companion and the ejecta. When the mass-loss timescale is shorter than the orbital period, the perturbing force can become greater than the gravitational attractive force and the motion can no longer be accurately described as an osculating Keplerian orbit. One then needs to integrate the equation of motion directly. The mass quadrupole tensor of the binary in these coordinates is

[TABLE]

and the third derivatives of each component are

[TABLE]

The GW luminosity and strain waveform can then be calculated with these expressions using Eqs. (7a) and (10), respectively. Eqs. (13) are general for any point-mass binary with a given equation of motion and time-dependent masses. These equations can also be translated into any choice of coordinate system. Since there are no forces acting in the $z$ direction, then the $z$ components will be zero and the orbital motion will remain coplanar in the $x$-$y$ plane.

In binary population synthesis codes, supernova mass loss in close binaries is often treated as being effectively instantaneous and conservative of the total orbital energy and orbital angular momentum. When treating the gravitational radiation, the mass-loss timescale should be bounded by the light-crossing timescale of the binary. The total energy and total angular momentum can be treated as conserved as long as the amount of energy and angular momentum radiated away in GWs is included in the balance equations. Supernovae in close binaries are the most extreme cases of mass loss and are thus expected to be the loudest sources of GWs from mass variation.

4.1 Numerical Comparison

We consider a close binary consisting of a neutron star and a helium star in an initially circular orbit with an initial semi-major axis of $a_{0}=10^{-2}\ {\rm au}$. We fix the neutron-star mass to be $m_{\rm NS,1}=1.4M_{\odot}$ and the helium-star mass to be $m_{\rm He}=2.4M_{\odot}$. The helium star may be large enough to be tidally deformed, but we ignore this effect and are solely interested in the GWs sourced by the supernova mass loss since this process operates on dynamical timescales, while the tidal effects on the orbit operate on secular timescales. The helium star will undergo a supernova explosion, which we take to be isotropic for simplicity and will form a neutron star with the same mass as the primary $m_{\rm NS,2}=m_{\rm NS,1}=1.4M_{\odot}$. The supernova explosion itself will be a source of GWs in the LIGO band, but here, we solely focus on the GWs sourced by the binary mass loss, whose characteristic frequency corresponds to the mass-loss timescale. We choose the following smooth function as a model for the mass-losing body such that the time derivatives of the mass are continuous:

[TABLE]

where $m_{2,\rm i}$ is the initial mass, $m_{2,\rm f}$ is the final mass, $\tau$ is a given timescale for the mass loss and $t_{0}$ is the the time at which the mass loss peaks. Given that the only constraint for $\tau$ is the light-crossing time $\tau_{\rm c}$ of the binary, we consider two cases: $\tau=[\tau_{\rm c},100\tau_{\rm c}]$. The light-crossing time relative to the orbital period for this particular system is $\Omega\tau_{\rm c}\approx 0.002$, so the mass-loss timescale for both cases are less than the orbital period. Other smooth functions may be used to model the sharp decrease in mass, but such choices will not affect our overall results.

We plot the orbital trajectory for the two cases in the top panel Figure 1. For each case, we integrate the equation of motion (see Eq. 1) with the mass of the secondary obeying Eq. (14). The peak of the mass loss is taken to occur at $t_{0}=P_{0}$, where $P_{0}$ is the initial orbital period of the binary (obtained from Kepler’s Third Law). The GW luminosity evaluated from Eqs. (7a) and (13) and the GW luminosity predicted from the PM formalism [i.e., excluding terms in Eqs. (13) with ${\cal O}(\geq\dot{\mu})$] are plotted in the bottom panel of Figure 1.

A spike in the GW luminosity occurs during the explosive mass loss and the amplitude of the spike increases with a decreasing mass-loss timescale. The dashed lines in the bottom panel of Figure 1 correspond to the GW luminosity computed with the PM formalism, i.e., ignoring the time-derivatives of the reduced mass in Eqs. (13). The maximum difference between our predicted GW luminosity and that of PM is about 4.18 and 1.37 for the $\tau_{\rm c}$ and $100\tau_{\rm c}$, respectively, in units of $\log_{10}\left(L_{\rm GW}/L_{\rm GW,0}\right)$. Even for the longer mass-loss timescale case, our predicted GW luminosity is about an order of magnitude larger than the PM prediction. The blue dashed line, corresponding to the PM case with $\tau=\tau_{\rm c}$ exhibits a peak which comes from the acceleration on the binary from the mass loss, i.e., the Blauuw kick. Since the mass loss is isotropic, the final orbital configuration is in the same orbital plane as the initial one. A more systematic and comprehensive survey of the parameter space is left for future study, especially since additional physics such as natal kicks from asymmetric supernovae will affect the GW luminosity and final state of the binary (e.g., Hills, 1983; Kalogera, 1996).

5 Discussion and Conclusions

In this Paper, we have shown that time-varying masses in close binaries contribute additional gravitational radiation along with the radiation produced by orbital motion. We have also shown that gravitational radiation from time-varying masses can be orders of magnitude greater than that of the orbital motion. This analysis is performed in the quadrupolar approximation of GR in the weak-field limit and the derived formulae reduce to the PM formulae when the mass of each binary component is taken to be constant.

Our analysis may be useful for implementations in binary population synthesis and generating approximate waveforms. It may also provide new testable predictions for the evolution of compact binaries. Applying the analysis of this Paper to higher-order post-Newtonian orders for higher-compactness binary systems and to the gravitational self-force for high to extreme mass-ratio binaries is the next extension of this work. These extensions may reveal richer physical behavior than the quadrupole approximation that is presented here. Future work will also include a more comprehensive survey of the parameter space including natal kicks from asymmetric supernovae.

Additional future studies will investigate how these new contributions to the gravitational radiation will affect the formation of compact binaries and the rates and distributions of compact-object mergers. The analysis we present can also be extended to hierarchical triple systems. A detection of GWs from mass transfer would not only validate this effect, but may also provide new and complementary tests of GR at different regimes of curvature and field strengths.

AMH is supported by the DOE NNSA Stewardship Science Graduate Fellowship under grant number DE-NA0003864. AMH thanks Charalampos Markakis, Eliu Huerta, Gabrielle Allen, Xin Liu, Charles Gammie, Nico Yunes, Fred Rasio, and Vicky Kalogera for helpful comments and fruitful discussions. PMR acknowledges support from NSF under AST 14-13367. .