Open Access
Issue
A&A
Volume 703, November 2025
Article Number A23
Number of page(s) 8
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202555169
Published online 31 October 2025

© The Authors 2025

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

Nearly all massive galaxies contain a super-massive black hole (SMBH) in their centre, as detailed in the review by Kormendy & Ho (2013). However, the growth channels of SMBHs and their initial seed masses are not well understood. Black holes can grow through merging processes between two SMBHs or through gas accretion when they are in an active galactic nucleus (AGN) state. In recent years, JWST has shed new light on SMBHs in the early universe, for instance, with the discovery of an active SMBH at a redshift of 10 (Maiolino et al. 2024); hence, the early stages of SMBH mass growth become more accessible to observations. From a theoretical modelling perspective, the merging of the SMBH in the galactic centres would be the bright source of gravitational wave emission (Reisswig et al. 2009; Holley-Bockelmann et al. 2020; Volonteri et al. 2003).

Despite the scepticism of the early 2000s about the SMBH binary (SMBHB) merging timescales (the so-called final-parsec problem), up-to-date realistic high-resolution cosmological and zoomed cosmological galactic-scale simulations point to a rather short timescale for the BHs merging in galactic centres. As leading examples of such works, we can refer to Illustris TNG50 (Nelson et al. 2019) and Romulus25 (Tremmel et al. 2017). In addition to the obvious great achievements of these simulations in the explanation of the formation of galactic-scale objects and long-term dynamical evolution, they also display some evident resolution issues (gravitational softening) with respect to detailed investigations of the possible SMBH mergers (see details in Khan et al. 2016). Also, based on cosmologically motivated zoomed high resolution simulations (Khan et al. 2016, 2018a; Koehn et al. 2023) confirmed the fast merging of the galaxies’ nuclei. In this process, the central galactic BHs are transferred to the common centre of the forming system on a dynamical friction timescale.

The SMBH binary merging process in galactic centers can dynamically be separated into a few stages (Begelman et al. 1980). During the initial stage of the galaxy merger, the SMBHs fall to the centre of the system due to dynamical friction, which can be expressed as

V ˙ df = 4 π G 2 ρ M bh V bh 2 χ ln Λ with χ = ρ ( < V bh ) ρ . $$ \begin{aligned} \dot{V}_{\rm df} = \frac{-4\pi G^2\rho M_{\rm bh}}{V_{\rm bh}^2}\chi \ln \Lambda \quad \mathrm{with} \quad \chi =\frac{\rho ( < V_{\rm bh})}{\rho }. \end{aligned} $$(1)

In general, the functions, χ, and the Coulomb logarithm, Λ, depend on the velocity of the massive object and the properties of the background system Binney & Tremaine (1987). For more details on this initial dynamical friction-dominated stage, we refer to Just et al. (2011). As a final product of this stage, we obtained the system of initially loosely bound BHs embedded in a large number of field stars.

The next dynamical stage is the tightening of the SMBH binary through repeated energy exchanges with the background stars (Merritt & Ferrarese 2001; Merritt 2001). At the end of this stage, we expect to see the hard binary SMBH already close to the mpc separation where the relativistic post-Newtonian effects have the capacity to further shrink the binary up to the GW merging phase (Brem et al. 2013; Sobolenko et al. 2017; Gualandris et al. 2022; Khan et al. 2015). In addition, more details on the available regularisation techniques to follow dynamical friction, scouring, and post-Newtonian evolution without gravitational softening in galaxy-scale and cosmological simulations can be found in literature (Rantala et al. 2017; Mannerkoski et al. 2023; Partmann et al. 2024; Zhou et al. 2025).

Tamfal et al. (2018), Liao et al. (2024), Gualandris et al. (2022), and Rawlings et al. (2023) explored the role of initial conditions and environmental factors in determining density profiles, nuclear star formation, and stochastic effects. During the final merging phase of the galactic SMBH, stochastic or random effects such as gravitational fluctuations are present. Such effects typically lead to the formation of a hard binary characterised by high eccentricity, whose evolution continues to become even more eccentric (for more details see Figs. 5 in Gualandris et al. 2022 and Sobolenko et al. 2022). A substantial body of numerical investigations has also been devoted to the interaction of the forming SMBH binary embedded in a mostly gaseous circum-nuclear or circum-binary disk (Fiacconi et al. 2013; Mayer 2013; Farris et al. 2014; Ryan & MacFadyen 2017).

NGC 7727 (Arp 222) is a peculiar spiral galaxy located about 27 Mpc away in the constellation Aquarius (Schweizer et al. 2018). It is notable for its strong tidal features, which are likely the result of a recent merger. Recently, it was discovered that the galaxy contains two separated nuclei and each of them hosts a SMBH (Voggel et al. 2022). One nucleus (hereafter Nuc 1) is in the photometric centre of the galaxy, whereas the second nucleus (Nuc 2 hereafter) is offset from the centre by only 500 pc. The second nucleus is probably the stripped nuclear star cluster of the merged galaxy. In Voggel et al. (2022), it was predicted that due to the small apparent separation of the two nuclei, they will likely merge within the next 250 million years. Therefore, this system could offer a rare, detailed insight into a double SMBH system that is close to merging its two SMBHs under the emission of gravitational waves (Moore et al. 2015).

The galaxy’s proximity to Earth makes it a great study case for galaxy merging and interaction as well as for the merger of two SMBHs. With the upcoming LISA mission (Moore et al. 2015), the understanding of double SMBH systems has become even more important, as their frequency is a direct information needed in the prediction of LISA event rates.

At the time of discovery, this was the closest known double SMBH system and, even more crucially, the closest in distance to us. This proximity enables us to study it in unprecedented detail, revealing a wealth of dynamical information that is not accessible with other, more distant candidates (e.g. Kollatschny et al. 2020).

The goal of this paper is to determine the shortest time required for the two SMBHs to merge. Earlier estimates based on a classic Chandrasekhar dynamical friction set the time of orbital migration by either of the SMBHs to be of the order of 250 Myr, for more details see Section 6.1 in Voggel et al. (2022). These estimates do not take into account the impact of anisotropic stellar velocities surrounding the SMBHs, which could be due, for example, to the residual rotation of the galactic system or streams of stars still observed on a kpc scale. The physical distance of ≈500 pc seen in the projection does not set strong constraints on the current orbit of the two sub-systems. In this contribution, we hypothesise that the projected distance is close to the actual three-dimensional relative positions, so that the orbit is in the plane of the sky.

The paper is organised as follows. In Section 2, we briefly describe the physical characteristics of NGC 7727, construct a simple physical model and corresponding numerical models. We also briefly discuss an N-body numerical code for the simulations. In Section 3 we present the dynamical timescales for the modelled SMBHB and depict their GW signal. We also discuss merger statistic from such types of objects. In Section 4, we summarise our results and conclusions.

2. Model construction for NGC 7727

2.1. Physical properties of NGC 7727

In a recent and comprehensive study of NGC 7727, Voggel et al. (2022) combined VLT kinematic data with HST photometry and Jeans anisotropic equilibrium galaxy models to confirm the presence of an SMBH in both nuclei. According to observational data, NGC 7727 has a total stellar mass of M* = 1.3 × 1011 M and resides in a dark matter (DM) halo with a mass of around 1013 M and a typical stellar-to-halo-mass relation (Behroozi et al. 2013).

Voggel et al. (2022) also determined the masses of the SMBHs using high spatial resolution spectroscopy from the MUSE instrument. They found that Nuc 1, located at the photometric centre of the galaxy, hosts a black hole with a mass of M1 = 1.54 × 108 M. The second nucleus Nuc 2 contains a smaller black hole with mass M2 = 6.33 × 106 M. A careful assessment of the background contamination by field stars led Voggel et al. (2022) to rule out a BH-free Nuc 2 at a 4.5σ confidence level. The two SMBHs are characterised by a mass ratio of ≈1 : 24, which makes the NGC 7727 a typical example of the minor merger SMBH system.

With an integrated bulge mass of 2.1 × 10^8 M, the Nuc 2’s SMBH makes up to 4.5% of the second nucleus, well above the ∼0.1% expected from the BH-bulge scaling relation (Kormendy & Ho 2013). For that reason, Nuc 2 has been identified as the surviving nuclear star cluster left over from a smaller galaxy that merged with NGC 7727. Additionally, optical emission lines indicate that Nuc 2 is a Seyfert-type low-luminosity AGN. The spread in colours and age templates reaching down to ≃1.5 Gyr makes it likely that the progenitor galaxy carried with it a gas-rich disc, which accounts for the younger stellar population forming throughout the merger. Thus Nuc 2 is in an advanced in-spiral phase, which will eventually lead to the merger of the two SMBHs: analytical estimates relying on Chandrasekhar dynamical friction suggest that the final binding and merger of the two nuclei will take place in less than 1 Gyr, based on a plausibly short 250 Myr timescale.

2.2. NGC 7727 Initial condition

To build the initial physical model, the estimated parameters of both nuclei and SMBHs, along with the stellar and gas components, were collected from existing observational data (Schweizer et al. 2018; Voggel et al. 2022).

The stellar mass of Nuc 1 is 1.14 × 109 M, while the Bulge is 51.6 × 109 M, giving us a total mass of around 5.7 × 109 M. We estimated the mass of Nuc 2 as 0.2 × 109 M. According to Voggel et al. (2022), the observed masses of the SMBHs are MBH = 1.54 × 108 M and MBH2 = 6.33 × 106 M. The fitted physical parameters of the central part of NGC 7727 are summarised in Table 1. The observed surface brightness distribution of the combined three Sérsic profiles is given in Fig. 1.

thumbnail Fig. 1.

Visualisation of the NGC 7727. Left image is from Voggel et al. (2022). Right image is from ESO’s VLT Survey Telescope.

Table 1.

NGC 7727 physical model parameters.

We initiated the galaxy merger from the dynamical system of two unbound central SMBHs located inside Nuc 1 and Nuc 2 with an initial separation of ΔR = 480 pc. As the type of orbit for the SMBH2, we chose an eccentric one (assuming that we are now at the pericentre) with orbital eccentricity e = 0.5. To test the assumptions made about the physical configuration of the two nuclei’s positions and orbital eccentricity, along with their impact on the outcome of our models, we performed calculations (not included in this contribution) with a wider initial physical separation between the nuclei obtained for a projection angle i varying from 30 to 60 degrees, increasing the separation from 554 pc up to 960 pc. The main effect is to delay the formation of the binary by up to at most 35 Myr. Hence, our choice of i = 0o clearly leads to the shortest merger timescale and provides a limiting case at the same time.

2.3. Numerical realizations of the NGC 7727 initial condition

In our physical model, each of the SMBHs is surrounded by its own initially bound stellar systems Nuc 1 and Nuc 2, which are located in the common stellar system known as the Bulge. All three density distributions are described by the generalised Sérsic model (Sérsic 1963). For the numerical realisation of the Sérsic distribution of each of the three components, we used the popular Python-based AGAMA library (Vasiliev 2019). The main parameters of the distributions we present in Table 1. We centred the Nuc 1 and Bulge systems in the origin of the coordinate system. We placed the Nuc 2 system at an initial distance of ΔR = 480 pc along the Z axis. We note that we generated the combined system with three different Sérsic index models and with three slightly different axis flattening modes as well.

To check the robustness of the numerical modelling results for the dynamical SMBHB orbital hardening, we constructed three different numerical realisations of the same physical model with three different particle numbers: N = 325k, 650k, and 1300k (see Table 2 with the numerical indexes A, B, C). The particle numbers for each component were chosen in such a way that the individual masses of field particles would end up with nearly the same masses. In this way, we were able to suppress the mass segregation effects (due to the different masses of particles from different components) in our dynamical N-body system during evolution.

Table 2.

NGC 7727 numerical models.

For each of these numerical realisations, we created two particle distributions with two randomised initial positions and velocities. Thus, in total, we obtained six independent numerical realisations of the same physical model, as detailed in Table 2. In Fig. 2, we present the density distribution at the initial time for our physical NGC 7727 configuration for the numerical realization B = 650k, rnd2.

thumbnail Fig. 2.

Initial density distributions at projection planes (X, Y) and (Y, Z) for numerical model B = 650k, rnd2. Black and dark blue dots represent the BHs in Nuc 1 and Nuc 2, respectively.

In Fig. 3, we present the evolution of the cumulative mass distribution profiles for the three different components Nuc 1, Nuc 2, and Tot = Nuc 1 + Nuc 2 + Bulge. The cumulative mass distribution at the start of the simulation is presented as black solid and dashed lines. The two SMBHs’ individual positions and masses are also represented at the same time.

thumbnail Fig. 3.

Evolution of the total cumulative mass distribution (black lines for the initial masses), for the bound stage (blue shaded lines), for the merged stage (red shaded lines), and post-Newtonian phase (green shaded lines) for model Brnd2. The labels on the shades besides black represent the passed time in the simulation in Myr.

During the N-body modelling, we used the NB scaling units (as with Hénon units1, here we set the gravitational constant G equal to 1). These numerical simulation units are based on the main observed physical parameters of the system; that is, the distance between the components Nuc 1 and Nuc 2 and the total mass of the system. In this way, we obtained the next NB scaling units for mass, distance, velocity, and time, respectively expressed as

MNB = 1.0 × 109 M,RNB = 100 pc,

VNB  = 207.4 km s−1,TNB = 0.471 Myr.

In these units, the speed of light velocity is cNB = 1445.6.

2.4. N-body code

For the global dynamical integration, we used the high-order parallel N-body code φ−GPU2, which is based on the fourth-order Hermite integration scheme with hierarchical individual block time steps (Berczik et al. 2011, 2013). One more feature in our integration scheme is the minimisation and forced synchronisation of the SMBH particle time steps. In this way, we could ensure that the SMBH particle times were synchronised with each other and this also set a minimum time step between all particles. We already thoroughly tested and used our modified φ−GPU code for a wide range of models that involve the evolution of multiple SMBHs in our previous works (Arca Sedda et al. 2019; Koehn et al. 2023; Berczik et al. 2024).

In the current version of the φ−GPU code, we used individual softening for different types of particles. For the mixed interaction between different types of particles, we used the equation for effective softening: εij = εji = 0.5 × (εi + εj). For the Nuc 1 + Nuc 2 particles, we set the gravitational softening εNuc = 0.1 pc and for the Bulge we set εBulge = 0.1 pc. As a result, for the mixed softening between Bulge and Nuc particles, we get εNuc,Bulge = 0.055 pc. For the SMBH interactions, we applied exactly zero softening for the SMBH particles: εBH = 0 pc. We used the same mixed (effective) softening for any interaction between field particles (star) and SMBH gravitational interactions, but with a special extra reduction factor: εBH,i = 0.01 × 0.5 × (εBH + εi). The effective softening for these types of interactions is given as εBH,Nuc = 5 · 10−5 pc and εBH,Bulge = 5 · 10−4 pc, respectively. Initially, we ran the simulation in a pure Newtonian regime and turned on the post-Newtonian (PN) terms during the progressive merging stage and stopped the simulation when the separation of the SMBHs fell below four Schwarzschild radii. This final time was assumed to be the merging time for the SMBHB. In the current implementation of the code, we added the PN formalism up to the 3.5PN term (including the spin-orbit and spin-spin interactions as well) for the relativistic orbit calculation of the SMBHs (Sobolenko et al. 2017, 2022). Among the PN terms, 2.5PN is the leading energy dissipation term, which ‘carries away’ most of the binary SMBHB energy from the system due to the emission of gravitational waves (Khan et al. 2018b,a; Berczik et al. 2024). In our current set-up, we used the non-Spin version of the PN formalism and limited the PN calculation up to the 2.5PN term.

3. Dynamical evolution and gravitational waves

3.1. Bound stage

We show the evolution of the SMBHB orbit parameters, such as the separation, ΔR, inverse semi-major axis, 1/a, and eccentricity, e, in Fig. 4 as part of our analysis of the evolution of the binary system. At T = 0.0 Myr, the SMBHs separated by ΔR are not gravitationally bound initially.

thumbnail Fig. 4.

Evolution of the black hole binary parameters at the stage of the pre-bound and bound systems. Left panel: Separations, ΔR, between SMBHB for the numerical models A, B, and C as solid lines. The dotted colour lines represent the second initial randomisation seeds for the particle distributions. Horizontal black dashed lines represent the level of the Schwarzschild radii for SMBHB with summarised mass: MSMBHB = MBH1 + MBH2. Middle pane: Evolution of the inverse semimajor axis for the same models. Black dashed lines show the simple linear fit for the ‘hardening’ of the orbits for A, B, and C models, respectively. Right panel: Evolution of eccentricity, e, for the SMBHB system for all numerical models.

In all of our models, the SMBHBs become bound after ≈60–65 Myr of dynamical evolution (Table 3), showing a quick inspiral time in the pure Newtonian regime (middle panel of Fig. 4). After this time, the inverse semi-major axis is constantly increasing almost as a linear function. The six numerical models presented here (in addition to significantly different particle numbers and different randomisation seeds for particle distributions) have consistent hardening rates ≈0.15 pc−1 Myr−1. We then see the formation of the hard binary within a timescale of ≈30 Myr afterwards, with a high eccentricity of e > 0.9, except in one model (Arnd1 ) where it remained smaller. In our case, this high eccentricity also strongly influences the final PN merging time of the system (Peters & Mathews 1963).

We conclude that all numerical models lead to a quite similar SMBHB evolution, independent of the variations in particle number or the range of random seed realizations. In all models, the pre-bound phase takes roughly the same time (≈60 − 65 Myr). After this phase, the Newtonian hardening phase (due to the three-body scattering of stars around SMBHB) has a similar hardening rate (close to ≈0.15 pc−1 Myr−1). Therefore, after ≈65 − 70 Myr of evolution, we arrive at similar conditions for our PN simulations in the final merging stage.

In Fig. 3, we present the detailed evolution of the cumulative mass distribution of our system. The initial distribution is presented in black. The prebound epoch is presented in a dark blue colour. The distributions during the bound stages are plotted in light blue colour. In the hard-binary phase, we present the distributions as a red palette, while green shows the results for the final PN phase of our runs. During the complex dynamical evolution of our NGC 7727 model, we see a slight change in the slope of the combined total cumulative mass distribution Tot from ≈2 to ≈3. This slope change can be interpreted as the formation of a dense core inside the central ≈1 − 2 pc of our simulation.

3.2. Merge stage

After the bound stage formed, we activated the PN description in the code (1PN, 2PN, and 2.5PN) for the three models A, B, C, with the first set of randomization rnd1. All three runs, with different particle numbers, produce the final merging of SMBHB around the same time. The PN terms were switched on between ≈65 − 70 Myr, then the SMBHB orbital pericentre separation was smaller than a few thousand Schwarzschild radii. In all systems, the final merging occurs at ≈120 − 140 Myr (see also Table 3).

Table 3.

NGC 7727 bound and merger times.

The detailed evolution of the SMBHB parameters for our largest N model simulation Crnd1 , is shown in Fig. 5. The general evolution is very similar to the plots we present in Fig. 4. Here, we also show the PN phase of the dynamical evolution of SMBHB. In this case, the dissipative PN terms become significant after ≈105 Myr, which we can clearly see in the e evolution, which corresponds to the SMBHB pericenter separation of ≈1000 Rsw.

thumbnail Fig. 5.

Combined evolution of the black holes’ parameters for the pre-merging and merging stages (including the PN phases) of the NGC 7727 SMBHB. Left upper panel: Separation ΔR between SMBHB for the numerical model Crnd1 . Horizontal black dashed lines represent 10, 100, and 1000 times the Schwarzschild radii for the SMBHBs’ mass. The solid black line shows the size of the semi-major axis at the time when the binary starts to harden. Right upper panel: Evolution of the inverse semi-major axis for the same model. Left and right bottom panels: Evolution of eccentricity, e, and pericentre for the SMBHB system. Different colours represent the different phases of the SMBHB evolution, with dark blue: unbound phase, blue: pre-merging, and light blue: merging, including the PN terms.

During the final phase of the merger, the pericentre distance of SMBHB falls below ≈10 Rsw (see the light blue part of the plots in Fig. 5). We stopped our runs when the separation of the SMBHB fell below 4 × (Rsw1 + Rsw2). Before the very final merging stage after ≈105 Myr, the SMBHB completely detaches from the surrounding stellar system. Thus, for showing detailed GW frequency plots, we switched to the two-body version of the same φ−GPU code with the very detailed particle output resolution, with a frequency of ≈100 points per orbit. The similarly detailed evolution of the other two Arnd1 and Brnd1 dynamical models are presented in Fig. 6. The three models have slightly different final SMBHB merging times, but all are within the range of 130 ± 10 Myr.

thumbnail Fig. 6.

Same as on Fig. 5, but for the numerical models Arnd1 and Brnd1 .

During the last phase of NGC 7727 modelling, we create a more frequent output of the SMBHB parameters to enable the calculation of the system’s GW signal. For the waveform calculation, we used simple GW quadrupole term expressions from Kidder (1995) and also Brem et al. (2013), Sobolenko et al. (2017):

h ij = 2 G μ D L c 4 [ Q ij + P 0.5 Q ij + P Q ij + P 1.5 Q ij + . . . ] . $$ \begin{aligned} h^{ij} = \frac{2G\mu }{D_{\rm L}c^4}~\left[Q^{ij}+P^{0.5}Q^{ij}+PQ^{ij}+P^{1.5}Q^{ij}+~...\right]. \end{aligned} $$(2)

Here, P is a correction term for the corresponding 𝒫𝒩 order, μ is the reduced mass, DL is the luminosity distance to the galaxy, and Qij is the quadrupole term. The last one can be written as

Q ij = 2 [ v i v j G M BH 12 r n i n j ] , $$ \begin{aligned} Q^{ij}=2\left[v^{i}v^{j}-\frac{GM_{\rm BH12}}{r}n^{i}n^{j}\right], \end{aligned} $$(3)

where vi and ni are the relative velocity and normalised position vectors in this reference frame, respectively.

For illustrative purposes, we neglected the higher order terms and, in this case, we can calculate the GW tensor in the source frame simply as

h ij 4 G μ D L c 4 [ v i v j G M BH 12 r n i n j ] . $$ \begin{aligned} h^{ij}\approx \frac{4G\mu }{D_{\rm L}c^{4}}\left[v^{i}v^{j}-\frac{GM_{\rm BH12}}{r}n^{i}n^{j}\right]. \end{aligned} $$(4)

We chose the virtual detector to be oriented in such a way that the coordinate axes coincide with the source frame, which allows us to not make any further coordinate transformations. We computed h+ and h× from hij, which give the relevant measurable strains in the ‘+’ and ‘×’ polarisations (Brem et al. 2013; Sobolenko et al. 2017).

In Fig. 7, we show the waveforms and amplitude-frequency maps for the last 100, 50, 10, and 5 years for the SMBHB system before its merger. We calculated the waveforms for h+ and h× polarisations and created the amplitude-frequency maps for the final phase of our model Crnd1 . It is worth noting that the PN approximation works reasonably well only for describing the early in-spiral SMBHB. Finally, numerical relativity and perturbation theory should be used for the full waveform map of the final SMBHB merging event.

thumbnail Fig. 7.

Time-frequency representations of the strain data for predicted gravitational waveforms of h+ polarisation from SMBHB merging of NGC 7727 for the last 100, 50, 10, and 5 yr in panels a, b, c, and d, respectively. Data refer to the Crnd1 model.

The frequencies obtained from our simulations for the SMBHB merging events (from ≈2 μHz up to 10 μHz, see Fig. 7) from SMBHBs with masses 1.5 × 108 M and 6.3 × 106 M at distances of DL = 27 Mpc are located just in between the sensitive curves of the pulsar timing array (PTA) and Laser Interferometer Space Antenna (LISA) consortiums. In other words, in principle, only during the very final phase of the merging would the GW signal from the NGC 7727 system fall within the detection range of LISA. The amplitude of the possible GW dimensionless signal strain (≈10−15 − 3 × 10−15) is well above the LISA detection limit for the final 10 μHz frequency range.

3.3. The merger statistics

We have shown that the SMBH merger would be around or just above the detection limit of the LISA interferometer. The event rates of such mergers that will be detected are clearly related to the galaxy merger rate. Early estimates of ≈0.008 per Gyr and per galaxy in the Local Universe (Toomre 1977; Tremaine 1981) are supported by more recent studies based on the ΛCDM paradigm (e.g. Conselice 2006; Huško et al. 2022). The merger rate per galaxy increases rapidly with cosmological redshift. For example, Huško et al. (2022) gives values of ≈10−2 per Gyr at z = 0.1, rising to ≈10−1 at z = 0.5; whereas Conselice (2006) gives an estimate of ≈1 per Gyr for a redshift range of z ≈ 1.5 to ≤4. These numbers are summed over galaxy DM halos in the brackets [108, 1012] M. Based on such galactic halo statistics and JWST detection of quasar activity already at z ≈ 7 and above, Liu & Inayoshi (2025) offers an estimate of 10−1 SMBH binary mergers per year when integrating out to z ≃ 5. For systems such as NGC 7727, the timescale for the SMBH to bind and then merge ( ≈ 130 Myr) is a factor of 2 shorter than previous estimates based on dynamical friction for that system. In future works, we plan to use other galaxy configurations to explore the question of whether the merging timescale is systematically shorter than estimates based on Chandrasekhar’s dynamical friction theory. We note that identifying which late-stage mergers also host a double BH system similar to NGC 7727 is difficult observationally, as it requires resolving the sphere of influence of the black holes (see Voggel et al. 2022 for a discussion of this issue). Generally, one may only find recent merging systems in wide-field images via their tidal features without knowing whether the merging galaxies each hosted an SMBH. Therefore, it is not clear to what extent a significant population of candidate merging SMBH binaries might be missing in terms of the Local Universe. The new data from the DESI survey for galaxy clustering (Yang et al. 2023; Adame et al. 2025) should provide a solid test bed from which to further constrain the number of sources detectable by LISA.

4. Conclusions

NGC 7727 is a particularly notable source because it offers a wealth of information on post-merger galaxy structures, showcasing the remnants of a smaller galaxy merging with a larger host. The presence of a secondary nucleus with an SMBH suggests complex dynamics and the survival of a stripped nuclear star cluster, a process that is often challenging to observe.

We carried out a dynamical modelling of the evolution of the SMBH system in the dense stellar environment in the nucleus of the NGC 7727 galaxy. Our physical model was constrained based on the latest available observational data (Voggel et al. 2022) for this object. We ran a set of numerical simulations using the well-tested dynamical N-body φ−GPU code. Our code, which also includes post-Newtonian corrections, allows us to follow the dynamical evolution of the SMBHB up to a few Schwarzschild radii separations. Based on the general physical model, we constructed three numerical models with different particle number realisations, namely, 325k, 650k, and 1300k, where each model also had one more additional seed of randomisations.

At the start of the simulation, the two black holes are not gravitationally bound. The system spends more than half of its final merging time of ≈130 Myr in the dynamical friction phase. Then, a quickly bound binary formsm with a high eccentricity of ≈0.98, and the two-body scattering phase proceeds for ≈60 − 120 Myr. Over the last ≈10 Myr, the BHs separation is seen to be rapidly shrinking due to the gravitational wave emission. The final merging time in our model is around 130 ± 10 Myr when starting from the actual separation today. Our study shows the principal possibility of the NGC 7727 binary black holes merging on a relatively short timescale. Such sources of gravitational waves are potential targets for the future LISA campaign during the final merging in the 10 μHz frequency range.


2

N-body code φ−GPU: https://github.com/berczik/phi-GPU-mole

Acknowledgments

The authors thank the anonymous referee for a very constructive report and suggestions that helped significantly improve the quality of the manuscript. The authors gratefully acknowledge computing resources used for this modelling on the Gauss Center for Supercomputing e.V. (GCS), through the John von Neumann Institute for Computing (NIC), with the GCS supercomputers JUWELS-booster at Julich Supercomputing Centre (JSC), Germany. A significant part of the simulations was also performed on the Luxembourg national supercomputer MeluXina (project p200574). The authors appreciate the LuxProvide teams for their expert support. PB, MI, OV, and MS thank the support from the special program of the Polish Academy of Sciences and the U.S. National Academy of Sciences under the Long-term program to support Ukrainian research teams, grant No. PAN.BFB.S.BWZ.329.022.2023. PB acknowledges the support of the National Science Foundation of China (NSFC) under grant No. 12473017. We thank the “Development of the concept for the first Kazakhstani orbital cislunar telescope - Phase I”, financed by the Ministry of Science and Higher Education of the Republic of Kazakhstan, No. BR24992759. This collaboration was made possible thanks to the support from the French-Ukrainian “DNIPRO” through grant number Ts 17/2–1 awarded in 2022-23 jointly to PB, MI, and CMB: we are grateful to all these agencies.

References

  1. Adame, A. G., Aguilar, J., Ahlen, S., et al. 2025, JCAP, 2025, 021 [CrossRef] [Google Scholar]
  2. Arca Sedda, M., Berczik, P., Capuzzo-Dolcetta, R., et al. 2019, MNRAS, 484, 520 [Google Scholar]
  3. Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980, Nature, 287, 307 [Google Scholar]
  4. Behroozi, P. S., Wechsler, R. H., & Conroy, C. 2013, ApJ, 770, 57 [NASA ADS] [CrossRef] [Google Scholar]
  5. Berczik, P., Nitadori, K., Zhong, S., et al. 2011, International conference on High Performance Computing, 8 [Google Scholar]
  6. Berczik, P., Spurzem, R., Wang, L., Zhong, S., & Huang, S. 2013, in Third International Conference "High Performance Computing", HPC-UA, 2013, 52 [Google Scholar]
  7. Berczik, P., Sobolenko, M., & Ishchenko, M. 2024, A&A, 687, L18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Binney, J., & Tremaine, S. 1987, Galactic dynamics (N.J, Princeton: Princeton University Press) [Google Scholar]
  9. Brem, P., Amaro-Seoane, P., & Spurzem, R. 2013, MNRAS, 434, 2999 [NASA ADS] [CrossRef] [Google Scholar]
  10. Conselice, C. J. 2006, ApJ, 638, 686 [NASA ADS] [CrossRef] [Google Scholar]
  11. Farris, B. D., Duffell, P., MacFadyen, A. I., & Haiman, Z. 2014, ApJ, 783, 134 [NASA ADS] [CrossRef] [Google Scholar]
  12. Fiacconi, D., Mayer, L., RoÅ¡kar, R., & Colpi, M. 2013, ApJ, 777, L14 [NASA ADS] [CrossRef] [Google Scholar]
  13. Gualandris, A., Khan, F. M., Bortolas, E., et al. 2022, MNRAS, 511, 4753 [NASA ADS] [CrossRef] [Google Scholar]
  14. Holley-Bockelmann, K., Bellovary, J., Bender, P., et al. 2020, ArXiv e-prints [arXiv:2012.02650] [Google Scholar]
  15. HuÅ¡ko, F., Lacey, C. G., & Baugh, C. M. 2022, MNRAS, 509, 5918 [Google Scholar]
  16. Just, A., Khan, F. M., Berczik, P., Ernst, A., & Spurzem, R. 2011, MNRAS, 411, 653 [Google Scholar]
  17. Khan, F. M., Holley-Bockelmann, K., & Berczik, P. 2015, ApJ, 798, 103 [Google Scholar]
  18. Khan, F. M., Fiacconi, D., Mayer, L., Berczik, P., & Just, A. 2016, ApJ, 828, 73 [Google Scholar]
  19. Khan, F. M., Capelo, P. R., Mayer, L., & Berczik, P. 2018a, ApJ, 868, 97 [Google Scholar]
  20. Khan, F. M., Berczik, P., & Just, A. 2018b, A&A, 615, A71 [EDP Sciences] [Google Scholar]
  21. Kidder, L. E. 1995, Phys. Rev. D, 52, 821 [Google Scholar]
  22. Koehn, H., Just, A., Berczik, P., & Tremmel, M. 2023, A&A, 678, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Kollatschny, W., Weilbacher, P. M., Ochmann, M. W., et al. 2020, A&A, 633, A79 [EDP Sciences] [Google Scholar]
  24. Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511 [Google Scholar]
  25. Liao, S., Irodotou, D., Johansson, P. H., et al. 2024, MNRAS, 528, 5080 [Google Scholar]
  26. Liu, H., & Inayoshi, K. 2025, Phys. Rev. D, 111, 043012 [Google Scholar]
  27. Maiolino, R., Scholtz, J., Witstok, J., et al. 2024, Nature, 627, 59 [NASA ADS] [CrossRef] [Google Scholar]
  28. Mannerkoski, M., Rawlings, A., Johansson, P. H., et al. 2023, MNRAS, 524, 4062 [Google Scholar]
  29. Mayer, L. 2013, Classical Quantum Gravity, 30, 244008 [Google Scholar]
  30. Merritt, D. 2001, ApJ, 556, 245 [Google Scholar]
  31. Merritt, D., & Ferrarese, L. 2001, MNRAS, 320, L30 [NASA ADS] [CrossRef] [Google Scholar]
  32. Moore, C. J., Cole, R. H., & Berry, C. P. L. 2015, Classical Quantum Gravity, 32, 015014 [Google Scholar]
  33. Nelson, D., Pillepich, A., Springel, V., et al. 2019, MNRAS, 490, 3234 [Google Scholar]
  34. Partmann, C., Naab, T., Rantala, A., et al. 2024, MNRAS, 532, 4681 [NASA ADS] [CrossRef] [Google Scholar]
  35. Peters, P. C., & Mathews, J. 1963, Phys. Rev., 131, 435 [Google Scholar]
  36. Rantala, A., Pihajoki, P., Johansson, P. H., et al. 2017, ApJ, 840, 53 [Google Scholar]
  37. Rawlings, A., Mannerkoski, M., Johansson, P. H., & Naab, T. 2023, MNRAS, 526, 2688 [NASA ADS] [CrossRef] [Google Scholar]
  38. Reisswig, C., Husa, S., Rezzolla, L., et al. 2009, Phys. Rev. D, 80, 124026 [Google Scholar]
  39. Ryan, G., & MacFadyen, A. 2017, ApJ, 835, 199 [Google Scholar]
  40. Schweizer, F., Seitzer, P., Whitmore, B. C., Kelson, D. D., & Villanueva, E. V. 2018, ApJ, 853, 54 [NASA ADS] [CrossRef] [Google Scholar]
  41. Sérsic, J.L. 1963, Bol. Asociacion Argentina Astronomia La Plata Argentina, 6, 41 [Google Scholar]
  42. Sobolenko, M., Berczik, P., Spurzem, R., & Kupi, G. 2017, Kin. Phys. Celestial Bodies, 33, 21 [Google Scholar]
  43. Sobolenko, M., Kompaniiets, O., Berczik, P., et al. 2022, MNRAS, 517, 1791 [NASA ADS] [CrossRef] [Google Scholar]
  44. Tamfal, T., Capelo, P. R., Kazantzidis, S., et al. 2018, ApJ, 864, L19 [NASA ADS] [CrossRef] [Google Scholar]
  45. Toomre, A. 1977, in Evolution of Galaxies and Stellar Populations eds. B. M. Tinsley, R. B. G. Larson, & D. Campbell, 401 [Google Scholar]
  46. Tremaine, S. 1981, in Structure and Evolution of Normal Galaxies eds. S. M. Fall, & D. Lynden-Bell, 84 [Google Scholar]
  47. Tremmel, M., Karcher, M., Governato, F., et al. 2017, MNRAS, 470, 1121 [NASA ADS] [CrossRef] [Google Scholar]
  48. Vasiliev, E. 2019, MNRAS, 482, 1525 [Google Scholar]
  49. Voggel, K. T., Seth, A. C., Baumgardt, H., et al. 2022, A&A, 658, A152 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Volonteri, M., Haardt, F., & Madau, P. 2003, ApJ, 582, 559 [Google Scholar]
  51. Yang, J., Fan, X., Gupta, A., et al. 2023, ApJS, 269, 27 [NASA ADS] [CrossRef] [Google Scholar]
  52. Zhou, Y., Mukherjee, D., Chen, N., et al. 2025, ApJ, 980, 79 [Google Scholar]

All Tables

Table 1.

NGC 7727 physical model parameters.

Table 2.

NGC 7727 numerical models.

Table 3.

NGC 7727 bound and merger times.

All Figures

thumbnail Fig. 1.

Visualisation of the NGC 7727. Left image is from Voggel et al. (2022). Right image is from ESO’s VLT Survey Telescope.

In the text
thumbnail Fig. 2.

Initial density distributions at projection planes (X, Y) and (Y, Z) for numerical model B = 650k, rnd2. Black and dark blue dots represent the BHs in Nuc 1 and Nuc 2, respectively.

In the text
thumbnail Fig. 3.

Evolution of the total cumulative mass distribution (black lines for the initial masses), for the bound stage (blue shaded lines), for the merged stage (red shaded lines), and post-Newtonian phase (green shaded lines) for model Brnd2. The labels on the shades besides black represent the passed time in the simulation in Myr.

In the text
thumbnail Fig. 4.

Evolution of the black hole binary parameters at the stage of the pre-bound and bound systems. Left panel: Separations, ΔR, between SMBHB for the numerical models A, B, and C as solid lines. The dotted colour lines represent the second initial randomisation seeds for the particle distributions. Horizontal black dashed lines represent the level of the Schwarzschild radii for SMBHB with summarised mass: MSMBHB = MBH1 + MBH2. Middle pane: Evolution of the inverse semimajor axis for the same models. Black dashed lines show the simple linear fit for the ‘hardening’ of the orbits for A, B, and C models, respectively. Right panel: Evolution of eccentricity, e, for the SMBHB system for all numerical models.

In the text
thumbnail Fig. 5.

Combined evolution of the black holes’ parameters for the pre-merging and merging stages (including the PN phases) of the NGC 7727 SMBHB. Left upper panel: Separation ΔR between SMBHB for the numerical model Crnd1 . Horizontal black dashed lines represent 10, 100, and 1000 times the Schwarzschild radii for the SMBHBs’ mass. The solid black line shows the size of the semi-major axis at the time when the binary starts to harden. Right upper panel: Evolution of the inverse semi-major axis for the same model. Left and right bottom panels: Evolution of eccentricity, e, and pericentre for the SMBHB system. Different colours represent the different phases of the SMBHB evolution, with dark blue: unbound phase, blue: pre-merging, and light blue: merging, including the PN terms.

In the text
thumbnail Fig. 6.

Same as on Fig. 5, but for the numerical models Arnd1 and Brnd1 .

In the text
thumbnail Fig. 7.

Time-frequency representations of the strain data for predicted gravitational waveforms of h+ polarisation from SMBHB merging of NGC 7727 for the last 100, 50, 10, and 5 yr in panels a, b, c, and d, respectively. Data refer to the Crnd1 model.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.