Massive Binaries Research Articles

The second data release from the European Pulsar Timing Array

Pulsar timing arrays offer a probe of the low-frequency gravitational wave spectrum (1–100 nHz), which is intimately connected to a number of markers that can uniquely trace the formation and evolution of the Universe. We present the dataset and the results of the timing analysis from the second data release of the European Pulsar Timing Array (EPTA). The dataset contains high-precision pulsar timing data from 25 millisecond pulsars collected with the five largest radio telescopes in Europe, as well as the Large European Array for Pulsars. The dataset forms the foundation for the search for gravitational waves by the EPTA, presented in associated papers. We describe the dataset and present the results of the frequentist and Bayesian pulsar timing analysis for individual millisecond pulsars that have been observed over the last ~25 yr. We discuss the improvements to the individual pulsar parameter estimates, as well as new measurements of the physical properties of these pulsars and their companions. This data release extends the dataset from EPTA Data Release 1 up to the beginning of 2021, with individual pulsar datasets with timespans ranging from 14 to 25 yr. These lead to improved constraints on annual parallaxes, secular variation of the orbital period, and Shapiro delay for a number of sources. Based on these results, we derived astrophysical parameters that include distances, transverse velocities, binary pulsar masses, and annual orbital parallaxes.

Binary masses and luminosities with Gaia DR3

Context. The recent third data release (DR3) of Gaia has brought some new exciting data about stellar binaries. It provides new opportunities to fully characterize more stellar systems and contributes to enriching our global knowledge of stellar behaviour. Aims. By combining the new Gaia non-single stars catalogue with double-lined spectroscopic binaries (SB2), we can determine the individual masses and luminosities of the components. To fit an empirical mass-luminosity relation in the Gaia G band, lower-mass stars must be added. These masses can be derived using Gaia-resolved wide binaries combined with literature data. Methods. Using the BINARYS tool, we combined the astrometric non-single star solutions in the Gaia DR3 with SB2 data from two other catalogues: the 9th Catalogue of Spectroscopic Binary orbits (SB9), and APOGEE. We also searched for low-mass stars that are resolved in Gaia with direct imaging and HIPPARCOS data or with a literature mass fraction. Results. The combination of Gaia astrometric non-single star solutions with double-lined spectroscopic data enabled us to characterize 43 binary systems with SB9 and 13 systems with APOGEE. We furthermore derived the masses of 6 low-mass binaries that are resolved with Gaia. We then derived an empirical mass-luminosity relation in the Gaia G band down to 0.12 ℳ⊙.

Born to Be Wide: The Distribution of Wide Binaries in the Field and Soft Binaries in Clusters

Most stars, binaries, and higher-multiplicity systems are thought to form in stellar clusters and associations that later dissociate. Very wide binaries can be easily disrupted in clusters due to dynamical evaporation (soft binaries) and/or tidal disruption by the gravitational potential of the cluster. Nevertheless, wide binaries are quite frequent in the field, where they can sometimes play a key role in the formation of compact binaries and serve as tools to study key physical processes. Here we use analytic tools to study the dynamical formation of soft binaries in clusters and their survival as field binaries following cluster dispersion. We derive the expected properties of very wide binaries both in clusters and in the field. We analytically derive their detailed distributions, including the wide binary fraction as a function of mass in different cluster environments, binary mass functions and mass ratios, and the distribution of their orbital properties. We show that our calculations agree well in most aspects with the results of N-body simulations but show some different binary fraction dependence on the cluster mass. We find that the overall fraction of wide binaries scales as , where N ⋆ is the size of the cluster, even for non-equal-mass stars. More massive stars are more likely to capture wide companions, with most stars above 5 M ⊙ likely to capture at least one stellar companion, and formation of triples is found to be frequent.

GRMHD study of accreting massive black hole binaries in astrophysical environment: A review

We present an overview of recent numerical advances in the theoretical characterization of massive binary black hole (MBBH) mergers in astrophysical environments. These systems are among the loudest sources of gravitational waves (GWs) in the universe and particularly promising candidates for multimessenger astronomy. Coincident detection of GWs and electromagnetic (EM) signals from merging MBBHs is at the frontier of contemporary astrophysics. One major challenge in observational efforts searching for these systems is the scarcity of strong predictions for EM signals arising before, during, and after merger. Therefore, a great effort in theoretical work to-date has been to characterize EM counterparts emerging from MBBHs concurrently to the GW signal, aiming to determine distinctive observational features that will guide and assist EM observations. To produce sharp EM predictions of MBBH mergers it is key to model the binary inspiral down to coalescence in a full general relativistic fashion by solving Einstein’s field equations coupled with the magnetohydrodynamics equations that govern the evolution of the accreting plasma in strong-gravity. We review the general relativistic numerical investigations that have explored the astrophysical manifestations of MBBH mergers in different environments and focused on predicting potentially observable smoking-gun EM signatures that accompany the gravitational signal.

Properties and merger signatures of galaxies hosting LISA coalescing massive black hole binaries

The gravitational wave (GW) antenna LISA will detect the signal from coalescing massive black hole binaries (MBHBs) of 104 − 107 M⊙, providing clues as to their formation and growth throughout cosmic history. Some of these events will be localised with a precision of several to less than a deg2, enabling the possible identification of their host galaxy. This work explores the properties of the host galaxies of LISA MBHBs below z ≲ 3. We generate a simulated lightcone using the semi-analytical model L-Galaxies applied to the merger trees of the high-resolution N-body cosmological simulation Millennium-II. The model shows that LISA MBHBs are expected to be found in optically dim (r > 20), star-forming (sSFR > 10−10 yr−1), gas-rich (fgas > 0.6), and disc-dominated (B/T < 0.7) low-mass galaxies of stellar masses 108 − 109 M⊙. However, these properties are indistinguishable from those of galaxies harbouring single massive black holes of comparable mass, making it difficult to select LISA hosts among the whole population of low-mass galaxies. Motivated by this, we explore the possibility of using merger signatures to select LISA hosts. We find that 40%−80% of the galaxies housing LISA MBHBs display merger features related to the interaction that brought the secondary MBH to the galaxy. Despite this, around 60% of dwarf galaxies placed in the surroundings of the LISA hosts will show these kinds of features as well, challenging the unequivocal detection of LISA hosts through the search for merger signatures. Consequently, the detection of an electromagnetic transient associated with the MBHB merger will be vital in order to pinpoint the star-forming dwarf galaxy where these binary systems evolve and coalesce.

The B-type Binaries Characterisation Programme – II. VFTS 291: a stripped star from a recent mass transfer phase

ABSTRACT Recent studies of massive binaries with putative black hole companions have uncovered a phase of binary evolution that has not been observed before, featuring a bloated stripped star that very recently ceased transferring mass to a main-sequence companion. In this study, we focus on the candidate system VFTS 291, a binary with an orbital period of 108 d and a high semi-amplitude velocity (K1 = 93.7 ± 0.2 km s−1). Through our analysis of the disentangled spectra of the two components, together with dynamical and evolutionary arguments, we identify a narrow-lined star of ∼1.5–2.5 $\, \mathrm{M}_\odot$ dominating the spectrum, and an early B-type main-sequence companion of 13.2 ± 1.5 $\, \mathrm{M}_\odot$. The low mass of the narrow-lined star, and the high mass ratio, suggest that VFTS 291 is a post-mass-transfer system, with the narrow-lined star being bloated and stripped of its hydrogen-rich envelope, sharing many similarities with other recently discovered stripped stars. Our finding is supported by our detailed binary evolution models, which indicate that the system can be well explained by an initial configuration consisting of an 8.1 $\, \mathrm{M}_\odot$ primary with an 8 $\, \mathrm{M}_\odot$ companion in a 7 d orbital period. While some open questions remain, particularly concerning the surface helium enrichment of the stripped star and the rotational velocity of the companion, we expect that high-resolution spectroscopy may help reconcile our estimates with theory. Our study highlights the importance of multi-epoch spectroscopic surveys to identify and characterize binary interaction products, and provides important insights into the evolution of massive binary stars.

Large-scale Evolution of Seconds-long Relativistic Jets from Black Hole–Neutron Star Mergers

We present the first numerical simulations that track the evolution of a black hole–neutron star (BH–NS) merger from premerger to r ≳ 1011 cm. The disk that forms after a merger of mass ratio q = 2 ejects massive disk winds (3–5 × 10−2 M ⊙). We introduce various postmerger magnetic configurations and find that initial poloidal fields lead to jet launching shortly after the merger. The jet maintains a constant power due to the constancy of the large-scale BH magnetic flux until the disk becomes magnetically arrested (MAD), where the jet power falls off as L j ∼ t −2. All jets inevitably exhibit either excessive luminosity due to rapid MAD activation when the accretion rate is high or excessive duration due to delayed MAD activation compared to typical short gamma-ray bursts (sGRBs). This provides a natural explanation for long sGRBs such as GRB 211211A but also raises a fundamental challenge to our understanding of jet formation in binary mergers. One possible implication is the necessity of higher binary mass ratios or moderate BH spins to launch typical sGRB jets. For postmerger disks with a toroidal magnetic field, dynamo processes delay jet launching such that the jets break out of the disk winds after several seconds. We show for the first time that sGRB jets with initial magnetization σ 0 > 100 retain significant magnetization (σ ≫ 1) at r > 1010 cm, emphasizing the importance of magnetic processes in the prompt emission. The jet–wind interaction leads to a power-law angular energy distribution by inflating an energetic cocoon whose emission is studied in a companion paper.

Iteratively comparing gravitational-wave observations to the evolution of massive stellar binaries

Iteratively comparing gravitational-wave observations to the evolution of massive stellar binaries

In search of the maximum lost momentum

We performed a series of 1381 full numerical simulations of high energy collision of two black holes to search for the maximum recoil velocity after their merger. We studied equal mass binaries with opposite spins pointing along the orbital plane to maximize asymmetric gravitational radiation and performed a search of spin orientations in the plane, impact parameters, and initial linear momenta to find a maximum recoil velocity extrapolated to the extreme spinning case of [Formula: see text][Formula: see text]km/s, thus tightly bounding recoil by [Formula: see text] the speed of light.

Ultimate Black Hole Recoil: What is the Maximum High-Energy Collision Kick?

We performed a series of 1381 full numerical simulations of high energy collision of black holes to search for the maximum recoil velocity after their merger. We consider equal mass binaries with opposite spins pointing along their orbital plane and perform a search of spin orientations, impact parameters, and initial linear momenta to find the maximum recoil for a given spin magnitude s. This spin sequence for s=0.4, 0.7, 0.8, 0.85, 0.9 is then extrapolated to the extreme case, s=1, to obtain an estimated maximum recoil velocity of 28,562±342 km/s, thus approximately bounded by 10% of the speed of light.

Implications of pulsar timing array observations for LISA detections of massive black hole binaries

ABSTRACT Pulsar timing arrays (PTAs) and the Laser Interferometer Space Antenna (LISA) will open complementary observational windows on massive black hole binaries (MBHBs), i.e. with masses in the range ${\sim} 10^6\!-\!10^{10}\, \rm M_{\odot }$. While PTAs may detect a stochastic gravitational wave background from a population of MBHBs, during operation LISA will detect individual merging MBHBs. To demonstrate the profound interplay between LISA and PTAs, we estimate the number of MBHB mergers that one can expect to observe with LISA by extrapolating direct observational constraints on the MBHB merger rate inferred from PTA data. For this, we postulate that the common signal observed by PTAs (and consistent with the increased evidence recently reported) is an astrophysical background sourced by a single MBHB population. We then constrain the LISA detection rate, $\mathcal {R}$, in the mass–redshift space by combining our Bayesian-inferred merger rate with LISA’s sensitivity to spin-aligned, inspiral–merger–ringdown waveforms. Using an astrophysically informed formation model, we predict a $95{{\, \rm per\, cent}}$ upper limit on the detection rate of $\mathcal {R} \lt 134\, {\rm yr}^{-1}$ for binaries with total masses in the range $10^7\!-\!10^8\, \rm M_{\odot }$. For higher masses, i.e. ${\gt} 10^8\, \rm M_{\odot }$, we find $\mathcal {R} \lt 2\, (1)\, \mathrm{yr}^{-1}$ using an astrophysically informed (agnostic) formation model, rising to $11 \, (6)\, \mathrm{yr}^{-1}$ if the LISA sensitivity bandwidth extends down to 10−5 Hz. Forecasts of LISA science potential with PTA background measurements should improve as PTAs continue their search.

Extreme mass ratios and fast rotation in three massive binaries

ABSTRACT The origin of rapid rotation in massive stars remains debated, although binary interactions are now often advocated as a cause. However, the broad and shallow lines in the spectra of fast rotators make direct detection of binarity difficult. In this paper, we report on the discovery and analysis of multiplicity for three fast-rotating massive stars: HD 25631 (B3V), HD 191495 (B0V), and HD 46485 (O7V). They display strikingly similar TESS light curves, with two narrow eclipses superimposed on a sinusoidal variation due to reflection effects. We complement these photometric data by spectroscopy from various instruments (X-Shooter, Espadons, FUSE...), to further constrain the nature of these systems. The detailed analyses of these data demonstrates that the companions of the massive OB stars have low masses (∼1 M⊙) with rather large radii (2–4 R⊙) and low temperatures (<15 kK). These companions display no UV signature, which would exclude a hot subdwarf nature, but disentangling of the large set of X-Shooter spectra of HD 25631 revealed the typical signature of chromospheric activity in the companion’s spectrum. In addition, despite the short orbital periods (P = 3−7 d), the fast-rotating OB-stars still display non-synchronized rotation and all systems appear young (<20 Myr). This suggests that, as in a few other cases, these massive stars are paired in those systems with non-degenerate, low-mass PMS companions, implying that fast rotation would not be a consequence of a past binary interactions in their case.

Internal kinematics of Gaia DR3 wide binaries: anomalous behaviour in the low acceleration regime

ABSTRACT The Gaia eDR3 catalogue has recently been used to study statistically the internal kinematics of wide binary populations using relative velocities of the two component stars, ΔV, total binary masses, mB, and separations, s. For s ≳ 0.01 pc, these binaries probe the low-acceleration a ≲ 2a0 regime where gravitational anomalies usually attributed to dark matter are observed in the flat rotation curves of spiral galaxies, where a0 ≈ 1.2 × 10−10 m s−2 is the acceleration scale of MOND. Such experiments test the degree of generality of these anomalies, by exploring the same acceleration regime using independent astronomical systems of vastly smaller mass and size. A signal above Newtonian expectations has been observed when a ≲ 2a0, alternatively interpreted as evidence of a modification of gravity, or as due to kinematic contaminants; undetected stellar components, unbound encounters, or spurious projection effects. Here I take advantage of the enhanced DR3 Gaia catalogue to perform a more rigorous study of the internal kinematics of wide binaries than what has previously been possible. Internally determined Gaia stellar masses and estimates of binary probabilities for each star using spectroscopic information, together with a larger sample of radial velocities, allow for a significant improvement in the analysis and careful exclusion of possible kinematic contaminants. Resulting ΔV scalings accurately tracing Newtonian expectations for the high acceleration regime, but markedly inconsistent with these expectations in the low acceleration one, are obtained. A non-Newtonian low acceleration phenomenology is thus confirmed.

Study of a sample of faint Be stars in the exofield of CoRoT

Context. The search and interpretation of non-radial pulsations from Be star light curves observed with the CoRoT satellite requires high-quality stellar astrophysical parameters. Aims. The present work is devoted to the spectroscopic study of a sample of faint Be stars observed by CoRoT in the fourth long run (LRA02). Methods. The astrophysical parameters were determined from the spectra in the λλ4000–4500 Å wavelength domain observed with the VLT/FLAMES instruments at ESO. Spectra were fitted with models of stellar atmospheres using our GIRFIT package. Spectra obtained in the λλ6400–7200 Å wavelength domain enabled the confirmation or, otherwise, a first identification of Be star candidates. Results. The apparent parameters (Teff, log g, Vsin i) for a set of 19 B and Be stars were corrected for the effects induced by the rapid rotation. These allowed us to determine: (1) stellar masses that are in agreement with those measured for detached binary systems; (2) spectroscopic distances that agree with the Gaia parallaxes; and (3) centrifugal/gravity equatorial force ratios of ~0.6–0.7, which indicate that our Be stars are subcritical rotators. A study of the Balmer Hα, Hγ and Hδ emission lines produced: (1) extents of the circumstellar disk (CD) emitting regions that agree with the interferometric inferences in other Be stars; (2) R– dependent exponents n(R) = ln[ρ(R)/ρo]/ln(Ro/R) of the CD radial density distributions; and (3) CD base densities ρo similar to those inferred in other recent works. Conclusions. The Hγ and Hδ emission lines are formed in CD layers close to the central star. These lines produced a different value of the exponent n(R) than assumed for Hα. Further detailed studies of Hγ and Hδ emission lines could reveal the physical properties of regions where the viscous transport of angular momentum to the remaining CD regions is likely to originate from. The subcritical rotation of Be stars suggests that their huge discrete mass-ejections and concomitant non-radial pulsations might have a common origin in stellar envelope regions that become unstable to convection due to rotation. If it is proven that the studied Be stars are products of binary mass transfer phases, the errors induced on the estimated Teff by the presence of stripped sub-dwarf O/B companions are not likely to exceed their present uncertainties.

Measuring spin in coalescing binaries of neutron stars that show double precursors

Gamma-ray bursts resulting from binary neutron-star mergers are sometimes preceded by precursor flares. These harbingers can be ignited by quasi-normal modes, excited by orbital resonances, shattering the stellar crust of one of the inspiralling stars up to ≳10 s before coalescence. In the rare case when a system displays two precursors, successive overtones of either interface modes or g modes can be responsible for the overstrainings. Since the free-mode frequencies of these overtones have an almost constant ratio, and the inertial-frame frequencies for rotating stars are shifted relative to static ones, the spin frequency of the flaring component can be constrained as a function of the equation of state, the binary mass ratio, the mode quantum numbers, and the spin-orbit misalignment angle. As a demonstration of the method, we find that the precursors of GRB090510 hint at a spin frequency range of 2 ≲ ν⋆/Hz ≲ 20 for the shattering star if we allow for an arbitrary misalignment angle, assuming ℓ = 2 g modes account for the events.

Modelling stellar evolution in mass-transferring binaries and gravitational-wave progenitors with metisse

ABSTRACT Massive binaries are vital sources of various transient processes, including gravitational-wave mergers. However, large uncertainties in the evolution of massive stars, both physical and numerical, present a major challenge to the understanding of their binary evolution. In this paper, we upgrade our interpolation-based stellar evolution code metisse to include the effects of mass changes, such as binary mass transfer or wind-driven mass loss, not already included within the input stellar tracks. metisse’s implementation of mass loss (applied to tracks without mass loss) shows excellent agreement with the sse fitting formulae and with detailed mesa tracks, except in cases where the mass transfer is too rapid for the star to maintain equilibrium. We use this updated version of metisse within the binary population synthesis code bse to demonstrate the impact of varying stellar evolution parameters, particularly core overshooting, on the evolution of a massive (25 and 15 M⊙) binary system with an orbital period of 1800 d. Depending on the input tracks, we find that the binary system can form a binary black hole or a black hole–neutron star system, with primary (secondary) remnant masses ranging between 4.47 (1.36) and 12.30 (10.89) M⊙, and orbital periods ranging from 6 d to the binary becoming unbound. Extending this analysis to a population of isolated binaries uniformly distributed in mass and orbital period, we show that the input stellar models play an important role in determining which regions of the binary parameter space can produce compact binary mergers, paving the way for predictions for current and future gravitational-wave observatories.

Failed supernovae as a natural explanation for the binary black hole mass distribution

Context. As the number of detected gravitational wave sources increases, the better we can understand the mass distribution of binary black holes (BBHs). This “stellar graveyard” shows several features, including an apparent mass gap that makes the distribution bimodal. In turn, the observed chirp mass distribution appears to be trimodal. Aims. We aim to investigate the extent to which we can explain the observed mass distribution based on stellar evolution, specifically with the hypothesis that the mass gap is caused by the difference between successful and failed supernovae (SNe). Methods. We posed a hypothetical remnant function, based on the literature of stellar evolution simulations, which relates initial mass to remnant mass, while including a “black hole island” and producing a bimodal remnant distribution. Moreover, we looked at observed type II SN rates in an attempt to detect the effect of failed SNe. Finally, using a simplified estimation of binary evolution, we determined the remnant distribution resulting from our remnant function and compared it with observations. Results. We find that failed SNe lower type II SN rates by approximately 25%, but the inferred rate from SN surveys is not accurate enough to confirm this. Furthermore, our estimation based on the remnant function produces a mass distribution that matches the general shape of the observed distributions of individual as well as chirp masses. Conclusions. Based on our research, we conclude that the failed SN mechanism and the presence of the black hole island are a natural hypothesis for explaining the individual BBH mass distribution and chirp mass distribution. However, to obtain a firmer conclusion, more detailed simulations are needed.

Ton S 415: a close binary containing an intermediate helium subdwarf discovered with SALT and TESS

ABSTRACT Hot subdwarfs are early-type low-mass stars lying between the main sequence and the white dwarf sequence. The majority have helium-burning cores and helium-poor surfaces. A minority have helium-rich surfaces. A few, including Ton S 415, have surfaces with between 10 per cent and 90 per cent helium (by number). Their properties are diverse and their origins mostly unknown. Ton S 415 was identified as a suspected binary in a Southern African Large Telescope (SALT) survey of hot subdwarfs from its large radial velocity. Using SALT spectroscopy and Transiting Exoplanet Survey Satellite light curves, we confirm that Ton S 415 is a close binary with an orbital period $p=84.6460\pm 0.0004\, {\rm min}$ and velocity semi-amplitude $K=175.5\pm 1.0\, {\rm km\, s^{-1}}$. Analysis of the SALT spectrum and broad-band spectral energy distribution shows the visible star to have an effective temperature $T_{\rm eff}=43\, 300\pm 1000\, {\rm K}$, surface gravity $\log g/{\rm cm\, s^{-2}}=5.89\pm 0.10$, surface helium-to-hydrogen ratio log y = −0.62 ± 0.10, radius $R=0.1074\pm 0.0025\, R_{\odot }$, and mass $M=0.33\pm 0.09\, {\rm {\rm M}_{\odot }}$. By modelling the light curve, we estimate the binary mass ratio to be q = 0.7 ± 0.3 and infer an unseen white dwarf companion with a mass of $M_{\rm WD}=0.47\pm 0.24\, {\rm {\rm M}_{\odot }}$. The results are consistent with a post-common-envelope binary that lost its hydrogen envelope before reaching the peak of the red giant branch, becoming a hot subdwarf with a non-canonical mass. We predict that the system will evolve into a double white dwarf binary before eventually either forming a stably accreting AM CVn system, or merging to form an R CrB star.

Evidence for mass-dependent peculiar velocities in compact object binaries: towards better constraints on natal kicks

ABSTRACT We compile a catalogue of low-mass and high-mass X-ray binaries, some recently reported binaries that likely host a neutron star (NS) or a black hole (BH), and binary pulsars (a pulsar and a non-degenerated companion) that have measured systemic radial velocities (γ). Using Gaia and radio proper motions together with γ, we integrate their Galactic orbits and infer their post-supernova (post-SN) 3D peculiar velocities ($v_\mathrm{pec }^{z=0}$ at Galactic plane crossing); these velocities bear imprints of the natal kicks that compact objects received at birth. With the sample totalling 85 objects, we model the overall distribution of $v_\mathrm{pec }^{z=0}$ and find a two-component Maxwellian distribution with a low- ($\sigma _v \approx 21\, \mathrm{km s^{-1}}$) and a high-velocity ($\sigma _v \approx 107\, \mathrm{km s^{-1}}$) component. A further comparison between distributions of binary subgroups suggests that binaries hosting high-mass donors/luminous companions mostly have $v_\mathrm{pec }^{z=0}\lesssim 100\, \mathrm{km s^{-1}}$, while binaries with low-mass companions exhibit a broader distribution that extends up to $\sim 400\, \mathrm{km s^{-1}}$. We also find significant anticorrelations of $v_\mathrm{pec }^{z=0}$ with binary total mass (Mtot) and orbital period (Porb), at over 99 per cent confidence. Specifically, our fit suggests $v_\mathrm{pec }^{z=0}\propto M_\mathrm{tot}^{-0.5}$ and $v_\mathrm{pec }^{z=0}\propto P_\mathrm{orb}^{-0.2}$. Discussions are presented on possible interpretation of the correlations in the context of kinematics and possible biases. The sample should enable a range of follow-up studies on compact object binary kinematics and evolution.

LitePIG: a lite parameter inference system for the gravitational wave in the millihertz band

We present a Python based parameter inference system for the gravitational wave (GW) measured in the millihertz band. This system includes the following features: the GW waveform originated from the massive black hole binaries (MBHB), the stationary instrumental Gaussian noise, the higher-order harmonic modes, the full response function from the time delay interferometry and the Gaussian likelihood function with the dynamic nested parameter sampler. In particular, we highlight the role of higher-order modes. By including these modes, the luminosity distance estimation precision can be improved roughly by a factor of 50, compared with the case with only the leading order (ℓ = 2, ∣m∣ = 2) mode. This is due to the response functions of different harmonic modes on the inclination angle are different. Hence, it can help to break the distance-inclination degeneracy. Furthermore, we show the robustness of testing general relativity (GR) by using higher-order harmonics. Our results show that the GW from MBHB can simultaneously constrain four of the higher harmonic amplitudes (deviation from GR) with 95% confidence level of , , and , respectively.