Brownian motion of supermassive black holes in galaxy cores

TL;DR

This paper models the Brownian motion of supermassive black holes in galaxy cores using a semi-analytic Langevin equation approach, validated against simulations, revealing typical displacements of a few parsecs.

Contribution

It introduces a semi-analytic Langevin-based model for SMBH dynamics in galactic cores, incorporating dynamical friction and stochastic effects, validated with N-body simulations.

Findings

SMBHs exhibit Brownian-like motion with displacements of a few parsecs.

Model aligns well with N-body simulations for intermediate-mass black holes.

Results are consistent with parameters of galaxy M87.

Abstract

We investigate the dynamics of supermassive black holes (SMBHs) in galactic cores by means of a semi-analytic model based on the Langevin equation, including dynamical friction and stochastic noise accounting for the gravitational interactions with stars. The model is validated against direct $N$-body simulations of intermediate-mass black holes in stellar clusters where a realistic number of particles is accessible. For the galactic case, we find that the SMBH experiences a Brownian-like motion with a typical displacement from the geometric center of the Galaxy of a few parsecs, for system parameters compatible with M87. \keywords{stellar dynamics, black hole physics, methods: n-body simulations, methods: statistical.

Equations

$$\ddot{\mathbf{r}}=-\nabla {\Phi }_{\mathrm{tot}}(\mathbf{r})-\eta (\mathbf{r},\mathbf{v})\mathbf{v}+\mathbf{F}(\mathbf{r}),$$

$$\rho (r)=\frac{\gamma -3}{4\pi }\frac{M{r}_c}{r^{\gamma }(r+r_c{)}^{4-\gamma }},$$

$$\rho (r)=\frac{3}{4\pi }\frac{M{r}_c^2}{(r_c^2+r^2{)}^{5/2}}.$$

$$\eta =4\pi G^2{n}_{\ast }{m}_{\ast }(M_{\mathrm{BH}}+m_{\ast })\ln \Lambda \frac{\Psi (v)}{v^3},$$

$$\Psi (v)=4\pi {\int }_0^{v}f(v_{\ast })v_{\ast }^2\mathrm{d}{v}_{\ast },$$

$$\Psi (v)=\mathrm{Erf}\left(\frac{v}{\sqrt{2}\sigma }\right)-\frac{2v\exp (-v^2/2{\sigma }^2)}{\sqrt{2\pi }\sigma },$$

$$H(F)=\frac{2}{\pi F}{\int }_0^{\infty }\exp \left[-\alpha (\xi /F{)}^{3/2}\right]\xi \sin (\xi )\mathrm{d}\xi ;\alpha =(4/15)(2\pi G{m}_{\ast }{)}^{3/2}{n}_{\ast },$$

$$\tilde{H}(F)\sim 2\pi n_{\ast }(G{m}_{\ast }{)}^{3/2}{F}^{-5/2}.$$

$$v(t+\Delta t)=c_2\left[c_{1}v(t)+\Delta t\nabla \Phi (x^{\prime })+d_1\tilde{F}(x^{\prime })\right]x(t+\Delta t)=x^{\prime }+\frac{\Delta t}{2}v(t+\Delta t),$$

Full text

††volume: 351††journal: Star Clusters: From the Milky Way to the Early Universe††editors: A. Bragaglia, M.B. Davies, A. Sills & E. Vesperini, eds.

1 Introduction

Here we report the preliminary results of our exploration of the dynamics of supermassive black holes (SMBHs) under the effect of gravitational scattering with the stars in galaxy cores.

Due to the high computational cost of direct $N-$body simulations, and the difficulty to attain realistic black hole-to-star mass ratios $M_{\rm BH}/m_{*}$ in the latter, we adopt a stochastic differential equation approach to simulate the effect of multiple gravitational “collisions” with field stars. At variance with previous studies involving the solution of the Fokker-Planck equation for the probability distribution function of position $\mathbf{r}$ and velocity $\mathbf{v}=\dot{\mathbf{r}}$ of the black hole (BH) [e.g. see [Chatterjee et al. (2002), Merritt (2015)] and references therein], here we use instead the Langevin equation

[TABLE]

where $\Phi_{\rm tot}$ is the smooth potential generated by the stellar and dark matter distributions, $\eta(\mathbf{r},\mathbf{v})$ is the dynamical friction [[Chandrasekhar (1943), Chandrasekhar (1949)]] coefficient, and $\mathbf{F}(\mathbf{r})$ a fluctuating force component (per unit mass).

2 Methods

We model the density profiles $\rho_{*}$ and $\rho_{DM}$ of stars and dark matter, respectively, using the family of spherical $\gamma-$models

[TABLE]

where $M$ is the total mass, $r_{c}$ the core radius, and $\gamma$ the central logarithmic density slope. In addition, for the comparisons with direct $N-$body simulations, we use the Plummer model

[TABLE]

The dynamical friction coefficient evaluated for a single component model with number density $n_{*}(\mathbf{x})$, particle mass $m_{*}$, mass density $\rho=n_{*}m_{*}$ and velocity distribution $f(v_{*})$ is

[TABLE]

where $G$ is the gravitational constant, $\ln\Lambda$ is the Coulomb logarithm, $v=||\mathbf{v}||$, and

[TABLE]

is the fractional velocity volume function. We assume that $\eta$ depends on $\mathbf{x}$ through $n_{*}$ [see e.g [Alessandrini et al. (2014), Arca-Sedda & Capuzzo-Dolcetta (2014a), Arca-Sedda & Capuzzo-Dolcetta (2014b)]] and that $f$ is approximated by a Maxwellian with velocity dispersion $\sigma(r)$ obtained solving the Jeans equations for the given density profile (2) or (3). In this case

[TABLE]

so $\Psi(\infty)=1$, [Binney & Tremaine (2008)].

The norm $F$ of the stochastic acceleration term in Equation (1) is sampled from the [Holtsmark (1919)] distribution

[TABLE]

introduced originally in the context of plasma physics, and used in stellar dynamics by [Chandrasekhar & von Neumann (1942), Chandrasekhar & von Neumann (1943)] to study the fluctuations of the gravitational field acting on a test star.

As Equation (7) cannot be written explicitly in terms of simple functions, one either solves numerically the integral, or expands in series the integrand up to the desired order integrating separately the terms of the sum [see e.g. [Hummer (1986)]]. In the limit of large $F$, Equation (7) is well approximated [see [Kandrup (1980), Petrovskaya (1986), Gabrielli et al. (1999), Bottaccio et al. (2002)]] by

[TABLE]

We integrate Eq. (1) with the quasi-symplectic scheme of [Mannella (2004)], which for the one dimensional case reads

[TABLE]

where $x^{\prime}=x(t)+{\Delta t}/{2}v(t)$. In the equations above $\Delta t$ is the constant time-step, $\tilde{F}$ the adimensional stochastic force, $c_{1}=1-{\eta\Delta t}/{2}$, $c_{2}=(1+\eta\Delta t/2)^{-1}$, $d_{1}=\sqrt{2\tau\eta\Delta t}$ and $\tau$ is fixed by the standard deviation of the distribution of $F$, $\langle F(x,t)F(x,t_{*})\rangle=2\eta\tau\delta(t-t_{*})$.

Since the Holtsmark distribution has by construction singular standard deviation, we are forced to pose a cut-off large $F$ in order to make it re-normalizable and usable in the Mannella scheme. Note that, for vanishing $\eta$ and $\tau$, Equations (9) yield back the standard symplectic leapfrog method.

3 Results and discussion

Equation (1) can be solved with any desired mass ratio $M_{\rm BH}/m_{*}$ and density profile $\rho$. In order to test the validity of our model we have compared the orbit of a BH in a system of $3\times 10^{4}$ stars distributed with a Plummer profile with mass ratio $M_{\rm BH}/m_{*}=100$ obtained with a standard direct $N-$body solver and with our stochastic method. In the case of the stochastic method we explored two different forms of the noise term, where $F$ is sampled from a Gaussian and a truncated Holtsmark distribution.

In Figure 1 (top panels) we show from left to right the orbit projections in the $(x,y)$ plane for the BH propagated in a direct $N-$body simulation and with the Langevin solver with Gaussian and Holtsmark noise terms. It appears that the case using the Holtsmak distribution better approaches the results of the $N-$body simulation (at least in term of radial displacement). In Figure 1 (bottom panels), we show the distribution $f(x)$ of displacement along $x$ and the associated velocity component distribution $f(v_{x})$, revealing that the $N-$body system has a somewhat intermediate behaviour between the stochastic models with Holtsmark and the Gaussian noise.

Once established the convergence between our model and the direct collisional $N-$body dynamics we have investigated the possibility that certain SMBHs that appear to be off-centered in the parent galaxy

[see e.g. the case of M87; [Batcheldor et al. (2010)]; but see also [Gebhardt et al. (2011)]] are effectively diffusing due to several dynamical encounters with stars.

In Figure 2, we show the evolution of the radial coordinate $r$ measured from the centre of the parent model (i.e. the centre of mass of the galaxy) for a $M_{\rm BH}=6\times 10^{9}M_{\odot}$ SMBH starting at rest in a galaxy with stellar and dark matter density distributions (2.1), total mass $M_{\rm gal}=3\times 10^{12}M_{\odot}$, $\gamma=1.2$ (for both stellar and dark components), $r_{c}=3$ kpc, and a dark to visible matter ratio of $\approx 6$, [parameters roughly corresponding to the case of M87, e.g. [Wu & Tremaine (2006)], [Event Horizon telescope collaboration (2019)]]. We observe that over a time of 10 Gyrs the SMBH reaches radii of the order of $\approx 6$ pc, that is compatible with the off-centre displacement claimed for the SMBH of M87, only due to multiple dynamical collisions with stars.

We suggest that the advantages (for instance in terms of computational time) given by the stochastic models can be crucial in the study of the processes involving SMBHs and stellar encounters.