EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 687 (July 2024) A&A, 687 (2024) L18 Full HTML
Open Access
Issue
A&A
Volume 687, July 2024
Article Number L18
Number of page(s) 13
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/202450141
Published online 16 July 2024

© The Authors 2024

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

Our understanding today of the co-evolution of galaxies and their supermassive black holes (SMBHs; Kormendy & Ho 2013) has been well established. In nearly all nearby spatially resolved galactic centres, an SMBH has been detected. Such a tight connection between SMBHs and their host galaxies sheds light on the mutual cosmological and dynamical existence of these components. The observed (and numerically modelled) galactic mergers inevitably lead us to the presence of SMBH mergers (Begelman et al. 1980; Kauffmann & Haehnelt 2000; Haehnelt & Kauffmann 2000). In previously proposed models of SMBH mergers, there has been a well known “final parsec problem” reported, as described in Berczik et al. (2005) and Milosavljević & Merritt (2003a,b). A possible solution to this problem was originally found in the non-axisymmetric (or even triaxial) dynamical loss cone repopulation mechanisms (Berczik et al. 2006, 2022; Nixon et al. 2011; Preto et al. 2011; Khan et al. 2013; Gualandris et al. 2022; Zhu & Chen 2024). The typical timescales of such mergers cover billions of years, whereby new galactic mergers may occur concurrently to an unfinished prior merger. Such a combination of events may result in the presence of three (or more) SMBHs (Valtonen et al. 1994; Hoffman & Loeb 2007; Ryu et al. 2018).

The study of the triple SMBH phase of merging galaxies is a complex and time-consuming dynamic problem even today. We know of several very significant studies devoted to this topic, namely, those of Makino & Ebisuzaki (1994), Hoffman & Loeb (2007), Amaro-Seoane et al. (2010) and Koehn et al. (2023). As a general conclusion from these Letters, we can assume that the triple SMBH stage is quite common over galactic merger history. Generally, as a result of such a triple dynamical evolution, we have two outcomes: ejection or coalescence (Bonetti et al. 2016, 2019; Mannerkoski et al. 2021). Although the triple interaction could potentially lead to the ejection of one or even all of the SMBHs (Valtonen et al. 1994), most systems are long-lived (≈1 Gyr; Hoffman & Loeb 2007); thus, a final coalescence is more common than ejection, confirming the analytical results of Makino & Ebisuzaki (1996). Based on the reports cited above, we can conclude that the coalescence scenario is more favourable compared with the ejection scenario.

The NGC 6240 galaxy (z = 0.0243, DL = 111.2 Mpc1) is one of the firstly directly detected dual active galactic nuclei (AGNs), which shows the signature of the galactic merging (Komossa et al. 2003; Puccetti et al. 2016). The galaxy centre shows intensive activity across different spectral ranges: the presence of H2O maser (Hagiwara et al. 2002, 2003), ultra-luminous infrared emission (Gerssen et al. 2004; Sanders et al. 2003; Iono et al. 2007), and X-ray emission (Pasquali et al. 2003; Fabbiano et al. 2020). The dual-AGN nature of the NGC 6240 nuclei has been seen in the extensive Multi-Element Radio Linked Interferometer Network observations at 1.4 Ghz and 5 Ghz by Beswick et al. (2001). The detected two main radio sources: north (N) and south (S), which are also well matched with the compact X-ray sources. The N nucleus can be clearly classified as an AGN according to the characteristics in the radio band observations. The S nucleus spectrum contains emission from the AGN and circumnuclear starburst and supernova remnants (Gallimore & Beswick 2004; Hagiwara et al. 2011). Recently, the results from Kollatschny et al. (2020) and Fabbiano et al. (2020) on the double structure of the S nucleus have further contributed to the discussion.

The SMBH mass in the S nuclei falls within the range of (0.87 − 2.0)×109 M⊙, based on high-resolution stellar kinematic results from Medling et al. (2011). Using K-band data from the Very Large Telescope and the classical SMBH-σ relation proposed by Tremaine et al. (2002), the SMBH masses for the N and S nuclei were estimated as (1.4 ± 0.4)×108 M⊙ and (2.0 ± 0.4)×108 M⊙, respectively (Engel et al. 2010). Engel et al. (2010) studied the motion of molecular gas through CO emission and associated it with circular movements. However, Treister et al. (2020) later linked it to turbulence motion. More recently, the velocity dispersion obtained with the MUSE instrument was found to correspond to an N nucleus black hole mass of (3.6 ± 0.8)×108 M⊙ and a combined S = S1+S2 nucleus black hole mass of (8.0 ± 0.8)×108 M⊙ (Kollatschny et al. 2020).

This object has already been investigated by our team in a few Letters (Sobolenko et al. 2021, 2022). In both Letters, we restricted our numerical modelling to the binary SMBH case. Based on our extensive post-Newtonian study, we subsequently estimated maximum merging timescale as ≈70 Myr for N and S SMBH components (see Fig. 5 in Sobolenko et al. 2022).

The Letter is organised as follows. In Sect. 2 we introduce the initial conditions for our physical and numerical model of the NGC 6240 triple SMBH system and shortly describe the dynamical integration procedure. In Sect. 3, we present the dynamical evolution of our SMBH triple system with von Zeipel–Lidov–Kozai (ZLK) oscillations. In Sect. 4, we present the results of our numerical model integration and summarise our findings.

2. Triple SMBH initial conditions

2.1. Physical model

To build the initial physical model, the estimated parameters of the central black holes, along with the stellar and gas components, were collected based on existing observational data (Kollatschny et al. 2020; Medling et al. 2011; Engel et al. 2010; Tecza et al. 2000; Tacconi et al. 1999). The selected physical parameters which were adopted in our numerical models are presented in Table 1. The values for the galaxy components stellar mass M* (S1, S2, N) were selected from the (see In Sect. 4.2 on page 10, Kollatschny et al. 2020). These masses were derived from the components dynamical Jeans modelling using the corresponding stellar velocity dispersion maps (see Fig. 10, Kollatschny et al. 2020).

Table 1.

NGC 6240 parameters of the physical model.

In Fig. 1 we present the main physical components of our model. The masses and Plummer radii of physical components are presented in Table 1. The observed separation between S1 and S2 is ΔR12 = 198 pc. The observed separation between N and the dynamical centre of binary S1+S2 is ΔR23 = 856 pc, according to Kollatschny et al. (2020). Because these distances are the projected separations between the components, it is obvious, that these values are the lower limits of the real 3D separations between SMBHs. To represent the large amount of gas which is placed between the S and N components, we added the so-called “stellar bulge” (B – between, in Table 1 and yellow colour in Fig. 1).

thumbnail Fig. 1.

NGC 6240 initial configuration for our triple SMBH physical model. The colour spheres represent the stellar bulges S1, S2, and N with central SMBHs (circles with dash edges) and stellar bulge B without a central black hole. SMBH initial velocities are noted VS1, S2, N. Initial separation for pair S1+S2 is noted as ΔR12 and for pair S+N is noted as ΔR23. The initial velocities are calculated in the assumption of initial orbital eccentricity e0 = 0.5 for both binary pairs. The radii of the spheres are scaled according to their bulge sizes (cutting the radius of the Plummer spheres).

Each component is surrounded by its bound stellar systems with simple Plummer density distribution (Plummer 1911) as follows:

ρ ( r ) = 3 M 0 4 π a 3 ( 1 + r 2 a 2 ) − 5 2 , $$ \begin{aligned} \rho (r) = \frac{3 M_0}{4 \pi a^3} \left( 1+\frac{r^2}{a^2} \right)^{-\frac{5}{2}}, \end{aligned} $$(1)

which produces the cumulative mass distribution:

M ( < r ) = M 0 r 3 ( r 2 + a 2 ) 3 / 2 , $$ \begin{aligned} M({ < }r) = M_0 \frac{r^3}{(r^2 + a^2)^{3/2}}, \end{aligned} $$(2)

where M0 is the total mass of each galactic bulge (total stellar mass M* in Table 1) and a is a scale factor that characterises the size of each nucleus (Plummer radius in Table 1). The a scale factors for each component were calculated according with the M0 masses, taking in account the central stellar densities reported in Kollatschny et al. (2020). Due to the flat central distribution of the Plummer profile, our estimated triple SMBH hardening will be smaller compared to the more centrally peaked core distribution profiles (Jaffe 1983; Hernquist 1990; Dehnen 1993). Thus, using the Plummer distribution we estimated the minimum numerical hardening for our S1 versus S2 and S versus N components.

2.2. Numerical realisations

Based on our physical model, we constructed our numerical realisations with different total particle numbers: N = 135 k, 270 k and 540 k. For each of these numerical realisations, we created the particle distributions with initial position and velocities for five different random realisations. Thus, we obtained a total of 15 independent numerical realisations of the same physical model (see Table 1 and Fig. 1). The dynamical evolution of one of the 540 k systems is shown in the video published on YouTube2.

To scale our dynamical modelling N-body or Hénon units3 (NB), we used the observed physical parameters of the system; namely the distance between the southern components S1 and S2 and the total mass of the S system. In this way, the NB units for mass, distance, velocity, and time can be respectively presented as

M NB = 1.3 × 10 10 M ⊙ , R NB = 198 pc , V NB = 531.4 km s − 1 , T NB = 0.364 Myr . $$ \begin{aligned}&M_{\rm NB} = 1.3\times 10^{10}\,M_{\odot }, \quad \quad R_{\rm NB} = 198\,\mathrm{pc},\\&V_{\rm NB} = 531.4\,\mathrm{km\,s}^{-1},\quad \quad T_{\rm NB} = 0.364\,\mathrm{Myr}. \end{aligned} $$

For the global dynamical integration, we used the high-order parallel N-body code φ-GPU4 which is based on the fourth-order Hermite integration scheme with hierarchical individual block time steps (Berczik et al. 2011, 2013). The current version of the φ-GPU code uses native GPU support and direct code access to the GPU using the NVIDIA native CUDA library. The current version of the code allows to have for each particle the individual softening. For the SMBH particles, we used exactly zero softening. For the SMBH versus field particle (star) interactions, we used some small effective softening: 10−6 × RNB ≈ 0.2 mpc. For the field particles self-interaction we also used a reasonably small softening parameter 10−4 × RNB ≈ 20 mpc.

The correct dynamical integration of the particles can be sensitive to the global integration time step parameter, η (Makino & Aarseth 1992). In the particular case of the fourth-order Hermite integration scheme, a time step can be written as

Δ t i = η · A i ( 1 ) A i ( 2 ) , $$ \begin{aligned}&\Delta t_{i} = \eta \cdot \frac{A_i^{(1)}}{A_i^{(2)}},\end{aligned} $$(3)

A i ( k ) = | a i ( k − 1 ) | | a i ( k + 1 ) | + | a i ( k ) | 2 , $$ \begin{aligned}&A_i^{(k)} = \sqrt{|\boldsymbol{a}_{i}^{(k-1)}||\boldsymbol{a}_{i}^{(k+1)}|+| \boldsymbol{a}_{i}^{(k)}|^{2}}, \end{aligned} $$(4)

where a(k) is the kth derivative of the i-th particle acceleration. Thus, the time steps are directly proportional to the η parameter, which is ultimately responsible for the total integration accuracy itself. For the highest-order integration Hermite scheme, the generalised Aarseth criterion can be found in Nitadori & Makino (2008). Based on our earlier experience (Sobolenko et al. 2022, 2021), for our set of models, we used η = 0.1 as a general parameter.

We also modified the standard Hermite hierarchical time step scheme, adding the extra ηSMBH = 0.05 for the SMBH particle criteria. One more feature in our integration scheme was a minimisation and forcing synchronisation of the SMBH particle time steps. This way, we always had the SMBH particle times synchronised between each other and with 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 (Arca Sedda et al. 2019; Mirza et al. 2017; Khan et al. 2018; Koehn et al. 2023).

Since in this case the post-Newtonian forces were not used in the simulations, we cannot talk about the physical merging of the black holes, only about the formation/not formation of the stable bound or unbound hierarchical triple systems (HTSs).

3. Dynamical evolution of the triple SMBH system with the ZLK oscillation

The results of the formation/non-formation of SMBH HTSs are presented for three numerical models with particle numbers 135 k, 270 k, and 540 k. Each of them had five different randomisation seeds for the initial positions and velocities of particles. As an illustration of the dynamical evolution of the system, we present a video of the density distribution for the highest particle number, 540 k, on YouTube5.

A summary of the results of modelling the SMBH dynamics in the NGC 6240 galaxy are presented in Table 2. As we see, we have two separate physical phases of system evolution. During the first phase of evolution, the SMBHs S1 and S2 become bound. This event, independently of particle numbers, happens around Tb, S1+S2 ≈ 3.6 Myr from the beginning of the integration. The second phase started when the N component became bound with the S1+S2 system, hereafter S. This happens around Tb, S+N ≈ 18 Myr from the beginning of integration, again independently of the particle numbers.

Table 2.

Simulation results of the formation or non-formation of hierarchical triple SMBH systems and follow-up ZLK oscillations.

In Fig. A.5, as an illustration, we present the evolution of separation between SMBH pairs for five realisations of the numerical model with 540 k particles: for S1+S2 and S+N. As a guiding line, we plot a level of ten times the sum of the three SMBH Schwarzschild radii. We can see that after the third SMBH (N component) arrives at the S component and becomes bound with S, there is a significant oscillation of ΔRS1+S2. One of the simplest explanations of this behaviour is a possible manifestation of the ZLK mechanism (von Zeipel 1910; Lidov 1961, 1962; Kozai 1962, for a historical overview and nomenclature, we refer to Ito & Ohtsuka 2019).

The system dynamical evolution for the other two field particles, numbers 135 k and 270 k, is also presented in Appendix A. We see the HTS formation and the subsequent ZLK oscillations only in four cases out of the five numerical randomisations in the case of the 135 k and 270 k particle numbers. For the highest particle number, 540 k, for all initial randomisation seeds, the HTS and ZLK mechanisms are observed. From this simple fact, we can conclude that in the case of a large enough particle number, the HTS and ZLK formation is inevitable in our initial model configuration.

In Fig. A.5, at the late phase of the dynamical evolution for the randomisations 1 (blue) and 2 (green) we see a classical ejection effect in these HTSs. The lowest-mass SMBH (S2) during the end phase of the simulation was ejected from the central region.

In Fig. A.1, we present the evolution of eccentricity e of the SMBH pairs for five realisations of the numerical model with 540 k particles: for S1+S2 and S+N. As we would expect due to the ZLK mechanism, we can see the significant oscillation of the S1+S2 binary SMBH eccentricity after the formation of bound S and N systems. The exact range of eccentricity oscillations depends on the particle randomisation, but in all cases, it is very significant and reaches nearly ≈1 for each case.

The evolution of the pericenter, rp (Fig. A.2), also shows the periodicity for inner S1+S2 binary. The inner binary pericenter oscillations are mainly defined by the eccentricity oscillations. The inverse semimajor axis of the S1+S2 binary during the ZLK phase of simulation is extending almost linearly (i.e. the semimajor axis of the inner binary is constantly shrinking); however, this happens on a much longer timescale compared to the rp oscillation period, which clearly reflects the eccentricity oscillations of the inner binary (see Fig. A.1).

In Fig. A.3, we present the mutual inclination angle between the angular momentum of the inner S1+S2 and outer S+N binaries (Naoz 2016, see Fig. 1). As we can see, the inclination angle evolution is in very good agreement with the expected theoretical limits of the Kozai angles: cos i min = ± 3 / 5 ≈ 39.2 ° $ \cos i_{\mathrm{min}} = {\pm}\sqrt{3/5} \approx 39.2^\circ $.

Here, we can note that our triple SMBH behaviour is more complex compared to the pure three-body ZLK effect. Together with the three SMBHs, our dynamical system contains a few hundred thousand field particles (see the video on YouTube6). Of course, the presented above spatial configuration is only one of the possible representations of the observed data set. In the present work, for simplicity and a quick set of the modelling, we limited ourselves to this single interpretation. In some sense, the chosen physical model works in favour of the ZLK oscillations in the triple SMBH system and it is enough to demonstrate this possible effect.

As an additional run for 540 k in the first randomisation, we generated the system with a smaller initial orbital inclination angle (≈75° instead of 90°) and with a slightly larger initial orbital eccentricity ≈0.6 instead of 0.5 of the S+N outer binary. The results of this extra run is presented in the top-right panel (with orange colour) in Figs. A.3–A.5.

For a detailed analysis of the ZLK oscillations in our triple systems, we ran separate extremely high time-resolution simulations with the data output time interval of ≈1.6 yr. These simulations were run for the 540 k for two intervals of the time around 18 Myr (in a moment of ZLK formation) and around 22 Myr (see middle and right panels in Fig. A.4). These time intervals are also marked as vertical black dotted lines in the left panels in Fig. A.4. Here, we see the stable ZLK effect with the almost constant periodicity and inclination angle amplitude.

As an extra test of the possible S1+S2 binary post-Newtonian merger attributed to the high orbital eccentricity of the ZLK effect, we started a set of restricted three-body SMBH simulations, including the combination of different post-Newtonian terms (PN0 – Newtonian, PN1 ∼ 1/c2, PN2 ∼ 1/c4, and PN2.5 ∼ 1/c5). The resulting plots of the eccentricity evolution of the different models are presented in Fig. B.1. The initial positions and velocities were chosen from the full scale N-body model, with a N = 540 k first randomisation at three moments of time in the first ZLK oscillation period (around ∼18.5 Myr), namely: when orbital eccentricity for the inner binary S1+S2 was 0.5, 0.7, and 0.9.

We individually turned on PN0 (“Newtonian”), PN0+PN1, PN0+PN1+PN2, only PN0+PN2.5, and PN-full = PN0+PN1+PN2+PN2.5 (i.e. with “full”) terms for three-body post-Newtonian code (Sobolenko et al. 2017, 2022). As we can see, the post-Newtonian terms (mainly the PN1) significantly change the whole ZLK effect evolution. With the PN-Full regime, we see the possible gravitational wave merger only for a very high initial eccentricity (≳0.9). The starting points with the lower initial eccentricities (0.5 and 0.7) simply show the non-growing, fixed eccentricity, evolution (i.e. basically no ZLK effect), with no subsequent gravitational wave merging for the S1+S2 inner binary in NGC 6240. In this sense, our set of restricted three-body post-Newtonian simulations support the ideas already presented in the works of Mannerkoski et al. (2021), Bonetti et al. (2016), and Blaes et al. (2002) – ultimately demonstrating that the PN1 term may suppress the ZLK effects.

4. Conclusions

We carried out a dynamical modelling of the evolution for the triple SMBH system in the dense stellar environment in the core of the NGC 6240 galaxy. Our physical model was constrained based on the latest available observational data (Kollatschny et al. 2020) and information from the literature (Komossa et al. 2003; Fabbiano et al. 2020) for this object. We run a set of numerical simulations using the well-tested dynamical N-body φ-GPU code. Based on the general physical model, we constructed three numerical models with different numbers of particles (135 k, 270 k, and 540 k), where each model had five independent sets of randomisations.

According to our simulation results, for all models and randomisations, the bound system from S1+S2 components formed at ≈3.6 Myr. The formation of a HTS triple bound system S+N occurs at ≈18 Myr. During the subsequent pure-Newtonian simulation of the evolution of triple SMBHs and the surrounding environment in NGC 6240, a ZLK oscillation was detected for the binary components S (the centre of mass of the S1+S2) and N systems. For these runs, the inclination angle between the orbital angular momentum of binary components aptly coincides with the theoretical calculations of the ZLK mechanism.

In contrast, our extra few-body post-Newtonian runs (for one randomisation case) show it is only for the large inner binary initial eccentricity (in our case ≳0.9) that we have a possibility of the inner binary merging due to the post-Newtonian energy radiation effects. For the lower eccentricity cases, the test runs show no sign of any possible merging and no hint of any ZLK oscillations in the system.

2

Triple SMBH evolution with the ZLK oscillation: https://www.youtube.com/watch?v=Jy919R6f61g

5

Triple SMBH evolution with the ZLK oscillation: https://www.youtube.com/watch?v=Jy919R6f61g

6

Triple SMBH evolution with the ZLK oscillation: https://www.youtube.com/watch?v=Jy919R6f61g

Acknowledgments

The authors thank the anonymous referee for a very constructive report and suggestions that helped significantly improve the quality of the manuscript. P.B. and M.S. thanks 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. P.B. and M.I. acknowledge the support of the Science Committee of the Ministry of Education and Science of the Republic of Kazakhstan (Grant No. AP14870501). M.S. acknowledges the support under the Fellowship of the President of Ukraine for young scientists 2022–2024. P.B. and M.S. gratefully acknowledge the Gauss Centre for Supercomputing e.V. (https://www.gauss-centre.eu) for funding this project by providing computing time through the John von Neumann Institute for Computing (NIC) on the GCS Supercomputer JUWELS at Julich Supercomputing Centre (JSC).

References

  1. Amaro-Seoane, P., Sesana, A., Hoffman, L., et al. 2010, MNRAS, 402, 2308 [NASA ADS] [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. Berczik, P., Merritt, D., & Spurzem, R. 2005, ApJ, 633, 680 [Google Scholar]
  5. Berczik, P., Merritt, D., Spurzem, R., & Bischof, H.-P. 2006, ApJ, 642, L21 [Google Scholar]
  6. Berczik, P., Nitadori, K., Zhong, S., et al. 2011, International Conference on High Performance Computing, 8 [Google Scholar]
  7. Berczik, P., Spurzem, R., Wang, L., Zhong, S., & Huang, S. 2013, Third International Conference “High Performance Computing”, HPC-UA 2013, 52 [Google Scholar]
  8. Berczik, P., Arca Sedda, M., Sobolenko, M., et al. 2022, A&A, 665, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Beswick, R. J., Pedlar, A., Mundell, C. G., & Gallimore, J. F. 2001, MNRAS, 325, 151 [NASA ADS] [CrossRef] [Google Scholar]
  10. Blaes, O., Lee, M. H., & Socrates, A. 2002, ApJ, 578, 775 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bonetti, M., Haardt, F., Sesana, A., & Barausse, E. 2016, MNRAS, 461, 4419 [NASA ADS] [CrossRef] [Google Scholar]
  12. Bonetti, M., Sesana, A., Haardt, F., Barausse, E., & Colpi, M. 2019, MNRAS, 486, 4044 [NASA ADS] [CrossRef] [Google Scholar]
  13. Dehnen, W. 1993, MNRAS, 265, 250 [NASA ADS] [CrossRef] [Google Scholar]
  14. Engel, H., Davies, R. I., Genzel, R., et al. 2010, A&A, 524, A56 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Fabbiano, G., Paggi, A., Karovska, M., et al. 2020, ApJ, 902, 49 [NASA ADS] [CrossRef] [Google Scholar]
  16. Gallimore, J. F., & Beswick, R. 2004, AJ, 127, 239 [NASA ADS] [CrossRef] [Google Scholar]
  17. Gerssen, J., van der Marel, R. P., Axon, D., et al. 2004, AJ, 127, 75 [NASA ADS] [CrossRef] [Google Scholar]
  18. Gualandris, A., Khan, F. M., Bortolas, E., et al. 2022, MNRAS, 511, 4753 [NASA ADS] [CrossRef] [Google Scholar]
  19. Haehnelt, M. G., & Kauffmann, G. 2000, MNRAS, 318, L35 [NASA ADS] [CrossRef] [Google Scholar]
  20. Hagiwara, Y., Diamond, P. J., & Miyoshi, M. 2002, A&A, 383, 65 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Hagiwara, Y., Diamond, P. J., & Miyoshi, M. 2003, A&A, 400, 457 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Hagiwara, Y., Baan, W. A., & Klöckner, H.-R. 2011, AJ, 142, 17 [NASA ADS] [CrossRef] [Google Scholar]
  23. Hernquist, L. 1990, ApJ, 356, 359 [Google Scholar]
  24. Hoffman, L., & Loeb, A. 2007, MNRAS, 377, 957 [NASA ADS] [CrossRef] [Google Scholar]
  25. Iono, D., Wilson, C. D., Takakuwa, S., et al. 2007, ApJ, 659, 283 [NASA ADS] [CrossRef] [Google Scholar]
  26. Ito, T., & Ohtsuka, K. 2019, Monogr. Environ. Earth Planets, 7, 1 [Google Scholar]
  27. Jaffe, W. 1983, MNRAS, 202, 995 [NASA ADS] [Google Scholar]
  28. Kauffmann, G., & Haehnelt, M. 2000, MNRAS, 311, 576 [Google Scholar]
  29. Khan, F. M., Holley-Bockelmann, K., Berczik, P., & Just, A. 2013, ApJ, 773, 100 [Google Scholar]
  30. Khan, F. M., Capelo, P. R., Mayer, L., & Berczik, P. 2018, ApJ, 868, 97 [Google Scholar]
  31. Koehn, H., Just, A., Berczik, P., & Tremmel, M. 2023, A&A, 678, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Kollatschny, W., Weilbacher, P. M., Ochmann, M. W., et al. 2020, A&A, 633, A79 [EDP Sciences] [Google Scholar]
  33. Komossa, S., Burwitz, V., Hasinger, G., et al. 2003, ApJ, 582, L15 [NASA ADS] [CrossRef] [Google Scholar]
  34. Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511 [Google Scholar]
  35. Kozai, Y. 1962, AJ, 67, 591 [Google Scholar]
  36. Lidov, M. L. 1961, Artif. Earth Satellite, 8, 5 [Google Scholar]
  37. Lidov, M. L. 1962, Planet. Space Sci., 9, 719 [Google Scholar]
  38. Makino, J., & Aarseth, S. J. 1992, PASJ, 44, 141 [NASA ADS] [Google Scholar]
  39. Makino, J., & Ebisuzaki, T. 1994, ApJ, 436, 607 [NASA ADS] [CrossRef] [Google Scholar]
  40. Makino, J., & Ebisuzaki, T. 1996, ApJ, 465, 527 [CrossRef] [Google Scholar]
  41. Mannerkoski, M., Johansson, P. H., Rantala, A., Naab, T., & Liao, S. 2021, ApJ, 912, L20 [NASA ADS] [CrossRef] [Google Scholar]
  42. Medling, A. M., Ammons, S. M., Max, C. E., et al. 2011, ApJ, 743, 32 [NASA ADS] [CrossRef] [Google Scholar]
  43. Milosavljević, M., & Merritt, D. 2003a, ApJ, 596, 860 [Google Scholar]
  44. Milosavljević, M., & Merritt, D. 2003b, in The Astrophysics of Gravitational Wave Sources, ed. J. M. Centrella, AIP Conf. Ser., 686, 201 [CrossRef] [Google Scholar]
  45. Mirza, M. A., Tahir, A., Khan, F. M., et al. 2017, MNRAS, 470, 940 [Google Scholar]
  46. Naoz, S. 2016, ARA&A, 54, 441 [Google Scholar]
  47. Nitadori, K., & Makino, J. 2008, New A, 13, 498 [NASA ADS] [CrossRef] [Google Scholar]
  48. Nixon, C. J., King, A. R., & Pringle, J. E. 2011, MNRAS, 417, L66 [NASA ADS] [CrossRef] [Google Scholar]
  49. Pasquali, A., de Grijs, R., & Gallagher, J. S. 2003, MNRAS, 345, 161 [NASA ADS] [CrossRef] [Google Scholar]
  50. Plummer, H. C. 1911, MNRAS, 71, 460 [Google Scholar]
  51. Preto, M., Berentzen, I., Berczik, P., & Spurzem, R. 2011, ApJ, 732, L26 [Google Scholar]
  52. Puccetti, S., Comastri, A., Bauer, F. E., et al. 2016, A&A, 585, A157 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Ryu, T., Perna, R., Haiman, Z., Ostriker, J. P., & Stone, N. C. 2018, MNRAS, 473, 3410 [Google Scholar]
  54. Sanders, D. B., Mazzarella, J. M., Kim, D. C., Surace, J. A., & Soifer, B. T. 2003, AJ, 126, 1607 [Google Scholar]
  55. Sobolenko, M., Berczik, P., Spurzem, R., & Kupi, G. 2017, Kinemat. Phys. Celest. Bodies, 33, 21 [NASA ADS] [Google Scholar]
  56. Sobolenko, M., Berczik, P., & Spurzem, R. 2021, A&A, 652, A134 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Sobolenko, M., Kompaniiets, O., Berczik, P., et al. 2022, MNRAS, 517, 1791 [NASA ADS] [CrossRef] [Google Scholar]
  58. Tacconi, L. J., Genzel, R., Tecza, M., et al. 1999, ApJ, 524, 732 [CrossRef] [Google Scholar]
  59. Tecza, M., Genzel, R., Tacconi, L. J., et al. 2000, ApJ, 537, 178 [CrossRef] [Google Scholar]
  60. Treister, E., Messias, H., Privon, G. C., et al. 2020, ApJ, 890, 149 [NASA ADS] [CrossRef] [Google Scholar]
  61. Tremaine, S., Gebhardt, K., Bender, R., et al. 2002, ApJ, 574, 740 [NASA ADS] [CrossRef] [Google Scholar]
  62. Valtonen, M. J., Mikkola, S., Heinamaki, P., & Valtonen, H. 1994, ApJS, 95, 69 [CrossRef] [Google Scholar]
  63. Zhu, L.-G., & Chen, X. 2024, ApJ, 960, 43 [NASA ADS] [CrossRef] [Google Scholar]
  64. von Zeipel, H. 1910, Astron. Nachr., 183, 345 [Google Scholar]

Appendix A: Evolution of the orbital parameters for the runs

thumbnail Fig. A.1.

Evolution of eccentricity e for S1+S2 and S+N systems for numerical model with particle number 540k and for the same five different randomisation seeds.

thumbnail Fig. A.2.

Evolution of pericenter, rpm for S1+S2 and S+N systems for a numerical model with the particle number 540k and for the five different randomisation seeds. The dotted grey line is the level of the ten Schwarzschild radii (as in Fig. A.5)

thumbnail Fig. A.3.

Evolution (ZLK oscillation) of the angle, i, between the S1+S2 and S+N SMBH binary systems for a numerical model with particle number 540k and for the five different randomisation seeds. The dotted grey lines represent the limiting theoretical Kozai angles.

thumbnail Fig. A.4.

Overall detail dynamical evolution of the selected 540k model. Left column: Long evolution till ≈30 Myr, where the vertical dotted black lines show the two intervals of the time for the further frequent output runs. Middle and right columns: Higher time resolution views of the ZLK effect for the simulations. Blue and light blue lines represents the parameters evolution for S1+S2 and S+N SMBHs from top to bottom: separation, ΔR, eccentricity, e, pericenter, rp, angle, i, between the S1+S2 and S+N systems. On the separation and pericenter plots, the dashed grey lines represent ten Schwarzschild radii for SMBHs with summarised mass. On the angle plot, the dashed grey line represents the theoretical limiting Kozai angles, imin.

thumbnail Fig. A.5.

Evolution of the separations, ΔR, between SMBHs in pairs S1+S2 and S+N for the numerical model with particle number 540k. Five plots show five different initial randomisation seeds for the particles distributions. The dotted grey line is the level of the ten Schwarzschild radii for SMBHs with summarised mass: M = MS1 + MS2 + MN.

thumbnail Fig. A.6.

Evolution of separations ΔR between SMBHs in pairs S1+S2 and S+N for numerical model with particle number 135k. Five plots show five different initial randomisation seeds for the particles distributions. The dotted grey line is the level of the ten Schwarzschild radii for SMBHs with summarised masses M = MS1 + MS2 + MN.

thumbnail Fig. A.7.

Evolution of eccentricity, e, for S1+S2 and S+N systems for numerical model with particle number 135k and for the same five different randomisation seeds.

thumbnail Fig. A.8.

Evolution of pericenter rp for S1+S2 and S+N systems for numerical model with the particle number 135k and for the five different randomisation seeds. The dotted grey line is the level of the ten Schwarzschild radii (as in Fig. A.5)

thumbnail Fig. A.9.

Evolution (ZLK oscillation) of the angle i between S1+S2 and S+N SMBH binary system for numerical model with particle number 135k and for the five different randomisation seeds. The dotted grey lines represent the limiting theoretical Kozai angles.

thumbnail Fig. A.10.

Same as Fig. A.6 but for N = 270k.

thumbnail Fig. A.11.

Same as Fig. A.7 but for N = 270k.

thumbnail Fig. A.12.

Same as Fig. A.8 but for N = 270k.

thumbnail Fig. A.13.

Same as Fig. A.9 but for N = 270k.

Appendix B: Evolution of the orbital parameters for the runs with PN terms

thumbnail Fig. B.1.

Evolution of the eccentricity for the 3-body PN integration model. Initial positions and velocities were chosen from the full scale N-body model with N = 540k first randomisation at three moments of time in the first ZLK oscillation period from top to bottom: when orbital eccentricity for the inner binary S1+S2 was 0.5, 0.7 and 0.9. We individually turned “on” the different combination of PN terms in our 3-body PN code: PN0 (“Newtonian”), PN0+PN1, PN0+PN1+PN2, only PN0+PN2.5 and PN-Full = PN0+PN1+PN2+PN2.5 (i.e. with “Full”).

All Tables

Table 1.

NGC 6240 parameters of the physical model.

In the text
Table 2.

Simulation results of the formation or non-formation of hierarchical triple SMBH systems and follow-up ZLK oscillations.

In the text

All Figures

thumbnail Fig. 1.

NGC 6240 initial configuration for our triple SMBH physical model. The colour spheres represent the stellar bulges S1, S2, and N with central SMBHs (circles with dash edges) and stellar bulge B without a central black hole. SMBH initial velocities are noted VS1, S2, N. Initial separation for pair S1+S2 is noted as ΔR12 and for pair S+N is noted as ΔR23. The initial velocities are calculated in the assumption of initial orbital eccentricity e0 = 0.5 for both binary pairs. The radii of the spheres are scaled according to their bulge sizes (cutting the radius of the Plummer spheres).
In the text
thumbnail Fig. A.1.

Evolution of eccentricity e for S1+S2 and S+N systems for numerical model with particle number 540k and for the same five different randomisation seeds.
In the text
thumbnail Fig. A.2.

Evolution of pericenter, rpm for S1+S2 and S+N systems for a numerical model with the particle number 540k and for the five different randomisation seeds. The dotted grey line is the level of the ten Schwarzschild radii (as in Fig. A.5)
In the text
thumbnail Fig. A.3.

Evolution (ZLK oscillation) of the angle, i, between the S1+S2 and S+N SMBH binary systems for a numerical model with particle number 540k and for the five different randomisation seeds. The dotted grey lines represent the limiting theoretical Kozai angles.
In the text
thumbnail Fig. A.4.

Overall detail dynamical evolution of the selected 540k model. Left column: Long evolution till ≈30 Myr, where the vertical dotted black lines show the two intervals of the time for the further frequent output runs. Middle and right columns: Higher time resolution views of the ZLK effect for the simulations. Blue and light blue lines represents the parameters evolution for S1+S2 and S+N SMBHs from top to bottom: separation, ΔR, eccentricity, e, pericenter, rp, angle, i, between the S1+S2 and S+N systems. On the separation and pericenter plots, the dashed grey lines represent ten Schwarzschild radii for SMBHs with summarised mass. On the angle plot, the dashed grey line represents the theoretical limiting Kozai angles, imin.
In the text
thumbnail Fig. A.5.

Evolution of the separations, ΔR, between SMBHs in pairs S1+S2 and S+N for the numerical model with particle number 540k. Five plots show five different initial randomisation seeds for the particles distributions. The dotted grey line is the level of the ten Schwarzschild radii for SMBHs with summarised mass: M = MS1 + MS2 + MN.
In the text
thumbnail Fig. A.6.

Evolution of separations ΔR between SMBHs in pairs S1+S2 and S+N for numerical model with particle number 135k. Five plots show five different initial randomisation seeds for the particles distributions. The dotted grey line is the level of the ten Schwarzschild radii for SMBHs with summarised masses M = MS1 + MS2 + MN.
In the text
thumbnail Fig. A.7.

Evolution of eccentricity, e, for S1+S2 and S+N systems for numerical model with particle number 135k and for the same five different randomisation seeds.
In the text
thumbnail Fig. A.8.

Evolution of pericenter rp for S1+S2 and S+N systems for numerical model with the particle number 135k and for the five different randomisation seeds. The dotted grey line is the level of the ten Schwarzschild radii (as in Fig. A.5)
In the text
thumbnail Fig. A.9.

Evolution (ZLK oscillation) of the angle i between S1+S2 and S+N SMBH binary system for numerical model with particle number 135k and for the five different randomisation seeds. The dotted grey lines represent the limiting theoretical Kozai angles.
In the text
thumbnail Fig. A.10.

Same as Fig. A.6 but for N = 270k.
In the text
thumbnail Fig. A.11.

Same as Fig. A.7 but for N = 270k.
In the text
thumbnail Fig. A.12.

Same as Fig. A.8 but for N = 270k.
In the text
thumbnail Fig. A.13.

Same as Fig. A.9 but for N = 270k.
In the text
thumbnail Fig. B.1.

Evolution of the eccentricity for the 3-body PN integration model. Initial positions and velocities were chosen from the full scale N-body model with N = 540k first randomisation at three moments of time in the first ZLK oscillation period from top to bottom: when orbital eccentricity for the inner binary S1+S2 was 0.5, 0.7 and 0.9. We individually turned “on” the different combination of PN terms in our 3-body PN code: PN0 (“Newtonian”), PN0+PN1, PN0+PN1+PN2, only PN0+PN2.5 and PN-Full = PN0+PN1+PN2+PN2.5 (i.e. with “Full”).
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.