Self consistent model for the evolution of eccentric massive black hole binaries in stellar environments: implications for gravitational wave observations

We construct evolutionary tracks for massive black hole binaries (MBHBs) embedded in a surrounding distribution of stars. The dynamics of the binary is evolved by taking into account the erosion of the central stellar cusp bound to the massive black holes, the scattering of unbound stars feeding the binary loss cone, and the emission of gravitational waves (GWs). Stellar dynamics is treated in a hybrid fashion by coupling the results of numerical 3-body scattering experiments of bound and unbound stars to an analytical framework for the evolution of the stellar density distribution and for the efficiency of the binary loss cone refilling. Our main focus is on the behaviour of the binary eccentricity, in the attempt of addressing its importance in the merger process and its possible impact for GW detection with the planned Laser Interferometer Space Antenna ({\it LISA}), and ongoing and forthcoming pulsar timing array (PTA) campaigns. We produce a family of evolutionary tracks extensively sampling the relevant parameters of the system which are the binary mass, mass ratio and initial eccentricity, the slope of the stellar density distribution, its normalization and the efficiency of loss cone refilling. We find that, in general, stellar dynamics causes a dramatic increase of the MBHB eccentricity, especially for initially already mildly eccentric and/or unequal mass binaries. When applied to standard MBHB population models, our results predict eccentricities in the ranges $10^{-3}-0.2$ and $0.03-0.3$ for sources detectable by {\it LISA} and PTA respectively. Such figures may have a significant impact on the signal modelling, on source detection, and on the development of parameter estimation algorithms.