The three-dimensional general relativistic Poynting-Robertson effect I: radial radiation field

TL;DR

This paper extends the understanding of the 3D motion of particles under the Poynting-Robertson effect in general relativity, specifically in Kerr spacetime with a radial radiation field, revealing new dynamical behaviors and equilibrium conditions.

Contribution

It develops the 3D equations of motion for particles in Kerr spacetime under radial radiation, extending previous 2D models, and analyzes the resulting particle dynamics and equilibrium surfaces.

Findings

Particles reach a critical hypersurface where radiation balances gravity.

Particles with angular momentum drift latitudinally towards the equator.

Particles without angular momentum stay at fixed latitude on the hypersurface.

Abstract

In this paper we investigate the three-dimensional (3D) motion of a test particle in a stationary, axially symmetric spacetime around a central compact object, under the influence of a radiation field. To this aim we extend the two-dimensional (2D) version of the Poynting-Robertson effect in General Relativity (GR) that was developed in previous studies. The radiation flux is modeled by photons which travel along null geodesics in the 3D space of a Kerr background and are purely radial with respect to the zero angular momentum observer (ZAMO) frames. The 3D general relativistic equations of motion that we derive are consistent with the classical (i.e. non-GR) description of the Poynting-Robertson effect in 3D. The resulting dynamical system admits a critical hypersurface, on which radiation force balances gravity. Selected test particle orbits are calculated and displayed, and their properties described. It is found that test particles approaching the critical hypersurface at a finite latitude and with non-zero angular moment are subject to a latitudinal drift and asymptotically reach a circular orbit on the equator of the critical hypersurface, where they remain at rest with respect to the ZAMO. On the contrary, test particles that have lost all their angular momentum by the time they reach the critical hypersurface do not experience this latitudinal drift and stay at rest with respects to the ZAMO at fixed non-zero latitude.

Tables

Table 1. Table 1: Explicit expressions of metric and ZAMO kinematical quantities for the Kerr metric.

Metric quantity	Explicit expression
$N={\left(-g^{tt}\right)}^{-1/2}$	${\left[\Delta /\rho \right]}^{1/2}$
$N^{\phi }=g_{t\phi }/g_{\phi \phi }$	$-2Mar/\left[\Sigma \rho \right]$
ZAMO quantity	Explicit expression
Radial components
$a{\left(n\right)}^{\hat{r}}$	$M/\left[\rho \sqrt{{\Sigma }^5\Delta }\right]\left\{{\Sigma }^2\left(r^2-a^2\right)+a^2{\sin }^2\theta \left[r^2\left(3r^2-4Mr+a^2\right)+a^2{\cos }^2\theta \left(r^2-a^2\right)\right]\right\}$
$\theta {\left(n\right)}^{\hat{r}}{}_{\hat{\phi }}$	$aM\sin \theta \left[\left(r^2+a^2\right)\left(\Sigma -2r^2\right)-2r^2\Sigma \right]/\left[\rho \sqrt{{\Sigma }^5}\right]$
$k_{\left(\mathrm{Lie}\right)}{\left(n\right)}^{\hat{r}}$	$-\sqrt{\Delta /{\Sigma }^5}\left[r{\Sigma }^2+a^{2}M{\sin }^2\theta \left(\Sigma -2r^2\right)\right]/\rho$
Polar components
$a{\left(n\right)}^{\hat{\theta }}$	$-a^{2}rM\sin \left(2\theta \right)\left[r^2+a^2\right]/\left[\rho \sqrt{{\Sigma }^5}\right]$
$\theta {\left(n\right)}^{\hat{\theta }}{}_{\hat{\phi }}$	$a^{2}rM\sin \left(2\theta \right)\sin \theta \sqrt{\Delta }/\left[\rho \sqrt{{\Sigma }^5}\right]$
$k_{\left(\mathrm{Lie}\right)}{\left(n\right)}^{\hat{\theta }}$	$-\sin \left(2\theta \right)\left[\left(r^2+a^2\right)\left(2a^{2}rM{\sin }^2\theta +{\Sigma }^2\right)+2a^{2}rM\Sigma {\sin }^2\theta \right]/\left[2\rho \sqrt{{\Sigma }^5}{\sin }^2\theta \right]$

Equations

$$\mathrm{d}{s}^2$$

$$+\frac{\Sigma }{\Delta }\mathrm{d}{r}^2+\Sigma \mathrm{d}{\theta }^2+\rho {\sin }^2\theta \mathrm{d}{\phi }^2,$$

$$\mathbf{n}=\frac{1}{N}({\partial }_{\mathbf{t}}-N^{\phi }{\partial }_{\mathbf{\phi }}),$$

$${\mathbf{e}}_{\hat{\mathbf{t}}}=\mathbf{n},{\mathbf{e}}_{\hat{\mathbf{r}}}=\frac{1}{\sqrt{g_{rr}}}{\partial }_{\mathbf{r}},$$

$${\mathbf{e}}_{\hat{\mathbf{\theta }}}=\frac{1}{\sqrt{g_{\theta \theta }}}{\partial }_{\mathbf{\theta }},{\mathbf{e}}_{\hat{\mathbf{\phi }}}=\frac{1}{\sqrt{g_{\phi \phi }}}{\partial }_{\mathbf{\phi }}.$$

$${\mathbf{\omega }}^{\hat{\mathbf{t}}}=Ndt,{\mathbf{\omega }}^{\hat{\mathbf{r}}}=\sqrt{g_{rr}}dr,$$

$${\mathbf{\omega }}^{\hat{\mathbf{\theta }}}=\sqrt{g_{\theta \theta }}d\theta ,{\mathbf{\omega }}^{\hat{\mathbf{\phi }}}=\sqrt{g_{\phi \phi }}(d\phi +N^{\phi }dt).$$

$$\mathbf{a}(n)$$

$$={\partial }_{\hat{r}}(\ln N){\mathbf{e}}_{\hat{\mathbf{r}}}+{\partial }_{\hat{\theta }}(\ln N){\mathbf{e}}_{\hat{\mathbf{\theta }}},$$

$${\mathbf{\theta }}_{\hat{\mathbf{\phi }}}(n)$$

$$=-\frac{\sqrt{g_{\phi \phi }}}{2N}({\partial }_{\hat{r}}{N}^{\phi }{\mathbf{e}}_{\hat{\mathbf{r}}}+{\partial }_{\hat{\theta }}{N}^{\phi }{\mathbf{e}}_{\hat{\mathbf{\theta }}}),$$

$${\mathbf{k}}_{(\mathrm{Lie})}(n)$$

$$=-[{\partial }_{\hat{r}}(\ln \sqrt{g_{\phi \phi }}){\mathbf{e}}_{\hat{\mathbf{r}}}+{\partial }_{\hat{\theta }}(\ln \sqrt{g_{\phi \phi }}){\mathbf{e}}_{\hat{\mathbf{\theta }}}].$$

$$T^{\mu \nu }={\Phi }^2{k}^{\mu }{k}^{\nu },k^{\mu }{k}_{\mu }=0,k^{\mu }{\nabla }_{\mu }{k}^{\nu }=0,$$

$$\mathbf{k}=E(n)[\mathbf{n}+\hat{\mathbf{\nu }}(k,n)],$$

$$\hat{\mathbf{\nu }}(k,n)=\sin \zeta \sin \beta {\mathbf{e}}_{\hat{\mathbf{r}}}+\cos \zeta {\mathbf{e}}_{\hat{\mathbf{\theta }}}+\sin \zeta \cos \beta {\mathbf{e}}_{\hat{\mathbf{\phi }}},$$

$$E(n)$$

$$=\frac{E+L_z{N}^{\phi }}{N}=\frac{E}{N}(1+b{N}^{\phi }),$$

$$E(n)\cos \beta \sin \zeta$$

$$=\mathbf{k}\cdot \frac{{\partial }_{\mathbf{\phi }}}{\sqrt{g_{\phi \phi }}}=\frac{L_z}{\sqrt{g_{\phi \phi }}}$$

$$\cos \beta$$

$$=\frac{bN}{\sin \zeta \sqrt{g_{\phi \phi }}(1+b{N}^{\phi })}.$$

$$k^t$$

$$k^r$$

$$k^{\theta }$$

$$k^{\phi }$$

$$R_{b,q}\left(r\right)$$

$${\Theta }_{b,q}\left(\theta \right)$$

$$q\equiv {\left(\frac{k_{\theta }}{k_t}\right)}^2+{\left[b\tan \left(\frac{\pi }{2}-\theta \right)\right]}^2-a^2{\cos }^2\theta .$$

$$\mathbf{k}={\partial }_{\mathbf{t}}+{\partial }_{\mathbf{r}}.$$

$${\Theta }_{b=0,q}(\theta )=0,$$

$$q=-a^2{\cos }^2\theta .$$

$$k^{\hat{\mu }}={\omega }_{\alpha }^{\hat{\mu }}{k}^{\alpha },$$

$$\cos \zeta =-\frac{k^{\hat{\theta }}}{k^{\hat{t}}}.$$

$$k^{\hat{\theta }}={\omega }_{\theta }^{\hat{\theta }}{k}^{\theta }=0.$$

$$0$$

$$={\partial }_r(\sqrt{-g}{\Phi }^2{k}^r).$$

$$\sqrt{-g}{\Phi }^2{k}^r$$

$$=const=E{\Phi }_0^2,$$

Full text

††preprint: APS/123-QED

1 Introduction

Matter motion in a gravitational field may be affected, among others, by radiation forces. The case in which both the gravity and radiation fields originate in the same body has been discussed extensively in astrophysical context for decades. The Eddington argument, for instance, describes the outward radial force exerted by momentum transfer by radiation from a star and determines the conditions under which such force balances the inward gravitational force. The corresponding critical luminosity of a star separates the regime of radial infall from that of radial escape. It has long been known also that radiation can remove angular momentum from matter in motion around a star. Radiation absorbed by a particle is in general re-emitted isotropically in the reference frame of the particle. However re-emission in the star reference frame is (slightly) beamed along the direction of motion owing to relativistic aberration, such that the particle recoils opposite to its velocity and a fraction of its angular momentum is transferred to the radiation field. This effect was first studied by Poynting in 1903 Poynting (1903), based on a classical treatment, and extended to special relativity by Robertson in 1937 Robertson (1937); it has since then been termed Poynting-Robertson (PR) effect (or radiation drag). Both authors had in mind applications to the motion of comets, dusts, and cm-size bodies in the solar system, for which Newtonian gravity suffices 111Note however that Robertson mentions general relativistic corrections to inertial terms in describing the perihelion shift in the quasi-Newtonian orbits Robertson (1937).. The Sun is treated as a point-like source of gravity and radiation in a flat Minkowski spacetime, with straight light rays propagating outwards; the test particle moves in a planar orbit and reradiates energy at the same rate at which it receives it from the radiation field. Special relativistic equations are written in the test particle reference frame and then transformed to the reference system of the Sun Robertson (1937). More advanced descriptions of test particle motion under the influence of PR effects were discussed in a some classical works in the 50 – 60’s (Wyatt and Whipple, 1950; Guess, 1962). A detailed description of the PR effect in relation to other radiation forces acting on small particles in the solar system was given in Burns et al. (1979).

The range of applications of PR drag to astrophysical problems has grown steadily since 80’s coming to encompass also compact objects, especially neutron stars (NSs) and black holes (BHs) that accrete matter down to the very strong gravitational fields in their vicinity. For instance Walker et al. Walker and Meszaros (1989); Walker (1992) studied the increase in mass accretion rate that is caused by PR drag when a bright thermonuclear flash occurs on the surface of a NS. This has motivated theory developments involving the PR effect in which GR is taken into account. The radial motion of test particles under the influence of a central isotropically emitting star was first investigated in the full GR in Abramowicz et al. (1990), and detailed calculations in the Schwarzschild metric presented. General relativistic equations of motion for a fluid in an arbitrary radiation field were formulated in Carroll (1990). Miller and collaborators Miller and Lamb (1996); Miller et al. (1998) investigated the velocity field of accreting matter as affected by the PR effect in the vicinity of a rotating NS, by carrying out approximate calculations in the Kerr metric.

A fully general relativistic treatment of the PR effect in the context of stationary and axially symmetric spacetimes was developed by Bini and collaborators Bini et al. (2009, 2011a). Similar to the classic model of Robertson these authors consider a compact object radiating as a point-like source, with photons traveling along null geodesics of the background spacetime (Schwarzschild or Kerr), and test particles moving in the equatorial plane around the compact object. Equations of motion are written in the ZAMO frame and then transformed to the rest reference frame of the compact object. The relativity of observer splitting formalism is adopted, which permits to clearly distinguish between gravitational and inertial contributions (see Jantzen et al. (1992); Bini et al. (1997a, b, 1999a, 1999b); De Falco et al. (2018), for further details). The equations are solved numerically and test particle trajectories and motion analysed for both a purely radial photon field (null impact parameter) Bini et al. (2009) and a photon field endowed with angular momentum (non-null impact parameter) Bini et al. (2011a).

Recent theoretical works on the extension of PR drag in GR have included studies of: test particle motion in the Vaidya spacetime Bini et al. (2011b) and around a slowly rotating relativistic star emitting isotropic radiation Oh et al. (2010); the general relativistic PR effect on a spinning test particle Bini and Geralico (2010); finite size effects Sok Oh et al. (2011), and the Lagrangian formulation of the general relativistic PR effect De Falco et al. (2018). More astrophysically-oriented studies of PR effect in strong gravitational fields have concentrated on: the development of the Eddington capture sphere concept around luminous stars, the surface where gravity, radiation, and PR forces balance Wielgus et al. (2012); Stahl et al. (2012, 2013); Wielgus (2016); Wielgus et al. (2016); the cosmic battery model in astrophysical accretion discs Koutsantoniou and Contopoulos (2014); Contopoulos et al. (2015); the dynamical evolution of accretion discs suddenly invested by a constant radiation filed (Bakala et al., 2018, A&A submitted, Lančová et al. (2017)). Research in this area has acquired further momentum from the growing body of observational evidence for PR effect in matter motion around compact objects, especially accreting NSs undergoing thermonuclear flashes Ballantyne and Strohmayer (2004); Ballantyne and Everett (2005); Worpel et al. (2013); Ji et al. (2014); Keek et al. (2014); Worpel et al. (2015); Keek et al. (2018).

Virtually all previous works on the general relativistic properties of the PR effect have been based upon a 2D model of the effect, i.e. planar (and arbitrarily oriented) orbits in spherically symmetric spacetimes (e.g. Schwarzschild’s) and equatorial orbits in the (axially symmetric) Kerr metric. A necessary improvement consists in developing the 3D theory of the PR effect in GR. That would allow to investigate the motion of test particles immersed in non-spherically symmetric radiation fields (e.g. latitude-dependent fields) and/or orbiting away from the equatorial plane of the Kerr metric. That is the aim of the present study, which builds on the formalism developed in Refs. Bini et al. (2009, 2011a). Our paper is structured as follows: in Sec. 2 we generalise to the 3D case the previous 2D equations for the PR effect in a stationary and axially symmetric general relativistic spacetimes. We adopt a simple prescription for the radiation, namely a field with zero angular momentum. In Sec. 3 we define the critical hypersurface on which radiation force balances gravity and discuss its salient features. In Sec. 4 we present calculations of selected orbits in the Schwarzschild and Kerr spacetimes; our concluding remarks are in Sec. 5.

2 Scenario and spacetime geometry

Our scenario for the description of the interaction between the radiation field and the motion of a test particle in the extreme gravitational field of a BH, or a NS, is constituted as follows: we consider the radiation field coming from an emitting region, located outside of the event horizon. The test particle motion is determined by its position in spherical coordinates and its velocity field in the ZAMO frame. The photon four-momentum is described by a pair of polar coordinates (see Fig. 1). In order to derive such set of equations we compute first the quantities in the ZAMO frame and then we transform them in the static observer frame. To deal with the relative motion of two non-inertial observers in GR we use the relativity of observer splitting formalism.

We consider a central compact object (BH or NS), whose outside spacetime is described by the Kerr metric with signature $(-,+,+,+)$ Kerr (1963). In geometrical units ($c=G=1$), the line element of the Kerr spacetime, $ds^{2}=g_{\alpha\beta}dx^{\alpha}dx^{\beta}$, in Boyer-Lindquist coordinates, parameterized by mass $M$ and spin $a$, reads as Boyer and Lindquist (1967)

[TABLE]

where $\Sigma\equiv r^{2}+a^{2}\cos^{2}\theta$, $\Delta\equiv r^{2}-2Mr+a^{2}$, and $\rho\equiv r^{2}+a^{2}+2Ma^{2}r\sin^{2}\theta/\Sigma$. The determinant of the Kerr metric is $g=\sqrt{(\Delta/\rho)g_{rr}g_{\varphi\varphi}g_{\theta\theta}}\equiv-\Sigma^{2}\sin^{2}\theta$.

2.1 ZAMO frame

ZAMOs are dragged by the rotation of the spacetime with angular velocity $\Omega_{\mathrm{ZAMO}}=-g_{\phi t}/g_{\phi\phi}$, while their radial coordinate remains constant. The four-velocity $\bm{n}$ of ZAMOs is the future-pointing unit normal to the spatial hypersurfaces, i.e. Bini et al. (1997a, b, 2009, 2011a),

[TABLE]

where $N=(-g^{tt})^{-1/2}$ is the time lapse function, $g^{tt}=g_{\varphi\varphi}/(g_{tt}g_{\varphi\varphi}-g_{t\varphi}^{2})$, and $N^{\varphi}=g_{t\varphi}/g_{\varphi\varphi}$ the spatial shift vector field. We focus our attention on the region outside the event horizon, where the time coordinate hypersurfaces are spacelike, i.e., $g^{tt}<0$. An orthonormal frame adapted to the ZAMOs is given by Bardeen et al. (1972)

[TABLE]

The relative dual tetrad of one-forms is given by

[TABLE]

All the indices associated to the ZAMO frame will be labeled by a hat, instead all the quantities measured in the ZAMO frame will be followed by $(n)$.

2.2 ZAMO kinematical quantities

Since the accelerated ZAMOs are locally nonrotating, their vorticity vector $\bm{\omega}(n)$ vanishes, but they have a nonzero expansion tensor $\bm{\theta}(n)$. For this reason it is more convenient to use the Lie transport (see Bini et al. (1997a); De Falco et al. (2018), for further details). The nonzero ZAMO kinematical quantities (i.e., acceleration $\bm{a}(n)=\nabla_{\bm{n}}\bm{n}$, expansion tensor along the $\hat{\varphi}$-direction $\bm{\theta_{\hat{\varphi}}}(n)$, also termed shear vector, and the relative Lie curvature vector $\bm{k}_{(\rm Lie)}(n)$) have only nonzero components in the $\hat{r}-\hat{\theta}$ plane of the tangent space Bini et al. (2009, 2011a); De Falco et al. (2018):

[TABLE]

In Table 1 we summarize the expressions of such quantities for the Kerr spacetime.

2.3 Radiation field

The stress-energy tensor, describing the radiation field, is modeled as a coherent flux of photons traveling along null geodesics in the Kerr geometry and acting on the test particle in the following manner Robertson (1937); Bini et al. (2009, 2011a)

[TABLE]

where parameter $\Phi$ is related to the intensity of the radiation field and $\bm{k}$ is the four-momentum field describing the null geodesics. The photon four-momentum, $\bm{k}$, and the photon spatial unit relative velocity with respect to the ZAMOs, $\bm{\hat{\nu}}(k,n)$, are respectively given by

[TABLE]

where $\beta$ and $\zeta$ are the two angles in the azimuthal and polar direction, respectively (see Fig. 1). The case of $\sin\beta>0$ corresponds to an outgoing photon beam (increasing $r$) while the case of $\sin\beta<0$ corresponds to an incoming photon beam (decreasing $r$, see Fig. 1). The photon four-momentum in the background Kerr geometry is identified by two impact parameters $(b,q)$, which are associated with two emission angles $(\beta,\zeta)$, respectively.

Using Eq. (7), the photon energy with respect to the ZAMO, $E(n)$, is expressed in the frame of a distant static observer by

[TABLE]

where $E=-k_{t}>0$ is the conserved photon energy, $L_{z}=k_{\varphi}$ is the conserved angular momentum along the polar $z$ axis orthogonal to the equatorial plane, and $b\equiv-k_{\phi}k_{t}=L_{z}/E$ is the first (azimuthal) photon impact parameter (constant of motion) Carter (1968); note that all these quantities are measured by a distant static observer Bini et al. (2011a).

This impact parameter is associated with the relative azimuthal angle $\beta$, measured in the ZAMO frame Bini et al. (2011a) (see Fig. 1). The angular momentum along the polar $\hat{\theta}$-axis in the ZAMO frame, $L_{z}(n)$, is expressed in the distant static observer frame by

[TABLE]

From such equation, we obtain

[TABLE]

An equation for $\zeta$ is needed to completely determine $\beta$.

The photon specific four-momentum components in the Kerr geometry are given by Chandrasekhar (1992)

[TABLE]

where $P\equiv r^{2}+a^{2}-b\,a$, and the pair of signs $s_{r}$, $s_{\theta}$ describes the orientation of the radial and latitudinal evolution, respectively Carter (1968). The radial and latitudinal effective potentials are respectively Chandrasekhar (1992):

[TABLE]

Here, $q$ is the second (latitudinal) photon impact parameter (constant of motion) related to the covariant components of the photon four-momentum through the relation Chandrasekhar (1992)

[TABLE]

2.3.1 Impact parameters

We consider here a radiation field which consists of photons moving in a purely radial direction at the infinity (this physically admissible because of the asymptotic flatness of the Kerr spacetime). In this case we have

[TABLE]

From Eq. (15) we note that the azimuthal component of the four-momentum, $k^{\varphi}$, vanishes at infinity. Moreover, in order to simplify the calculations we assume that the azimuthal impact parameter of the radiation field, $b$, takes null value, i.e., $b=0$. The latitudinal impact parameter, $q$, can be calculated from the condition

[TABLE]

which results from the absence of latitudinal photon motion ($k^{\theta}=0$). From Eqs. (13) and (16) we can express $q$ as a function of the polar angle $\theta$:

[TABLE]

This is possible because the latitudinal potential, Eq. (13), is independent of the radial coordinate and therefore the polar angle $\theta$ along a given photon trajectory is conserved. Photons with a given value of $q$ move only in the radial and azimuthal directions on the surface of the cone with the vertex located in the coordinates origin and with the vertex angle $\theta$ given by Eq. (17). Note that the azimuthal motion on finite values of the radial coordinate is caused only by frame dragging.

The above-defined radiation field significantly simplifies the integration of test particle trajectories in that only a single photon beam, described by the constants of motion $b=0,\ q=-a^{2}\cos^{2}\theta$, must be considered at the test particle position. In such case, the radial potential, Eq. (12), is always positive above the event horizon: this proves that the radiation field reaches every positions (for all $r$ and $\theta$) above the event horizon. The second constant of motion $q$ ranges in the interval $[-a^{2},0]$. The value $q=0$ corresponds to the motion of photons in the equatorial plane, while the value of $q=-a^{2}$ corresponds to the motion of photons along the polar axis on the south or north directions. We note also that since $q\leq 0$ radiation field photons can never cross the equatorial plane.

The local components of the photon four-momentum in the ZAMO frame are obtained through the following transformation:

[TABLE]

where $\omega^{\hat{\mu}}_{\,\,\,\alpha}$ represents the transformation matrix from the holonomic basis $\bm{\partial_{\alpha}}$ to the anholomic (tetaraed) basis $\bm{e_{\hat{\alpha}}}$, see Eq. (3) for determining its components. The local polar direction $\zeta$ of the photon four-momentum is given by (see Fig. 1)

[TABLE]

For the considered radiation field ($q=-a^{2}\cos^{2}\theta\,,k^{\theta}=0$) we simply obtain

[TABLE]

Consequently from Eqs. (19) and (20), the local polar direction of the radiation field photons in the ZAMO frame is always $\zeta=\pi/2$. From Eq. (9) and $b=0$, the local azimuthal direction $\beta$ of the photon four-momentum is $\cos\beta=0$ (see Fig. 1). Therefore, the local azimuthal angle of the test field photons in the ZAMO frame always take the value of $\beta=\pi/2$. We can conclude that in all ZAMO frames radiation field photons move in a purely radial direction. The source of the radiation field can thus be considered as centered in the coordinate origin, (differentially) rotating with a latitude-dependent angular velocity $\Omega_{\mathrm{ZAMO}}$ and emitting photons only along the radial direction in the appropriate, locally comoving ZAMO frame.

2.3.2 Intensity parameter

Since the photon four-momentum $\bm{k}$ is completely determined by $(b\,,q)$, the coordinate dependence of $\Phi$ then follows from the conservation equations $\nabla_{\beta}T^{\alpha\beta}=0$. Exploiting the absence of photon latitudinal motion ($k^{\theta}=0$) and symmetries of the Kerr spacetime, these can be written as

[TABLE]

Therefore, we have

[TABLE]

where $\Phi_{0}$ is a new constant related to the intensity of the radiation field at the emitting surface. This equation, however, does not fix the intensity parameter unambiguously. In fact, the conservation equations will be fulfilled, even if we multiply the constant expression $E\Phi_{0}^{2}$ by an arbitrary function of the $\theta$ coordinate. Thus this condition determines the class of radiating fields that differ from one another by the latitudinal dependence of the intensity. A radiation field which is independent of latitude (and whole intensity parameter is thus independent of $\theta$) is a natural choice, especially in the Schwarzschild limit because of its spherical symmetry. This can easily achieved by multiplying (22) by a factor $\sin\theta$, such that the intensity parameter becomes

[TABLE]

where we have used Eqs. (9) and (10), together with the fact that $N|b\tan\beta|=\sin\zeta\sqrt{g_{\varphi\varphi}}$ for $b=0$. In a Schwarzschild spacetime limit Eq. (23) thus reads

[TABLE]

which matches the spacetime spherical symmetry.

2.4 Test particle motion

We consider a test particle moving in the 3D space, with four-velocity $\mathbb{U}$ and spatial three-velocity with respect to the ZAMOs, $\bm{\nu}(U,n)$:

[TABLE]

where $\gamma(U,n)=1/\sqrt{1-||\bm{\nu}(U,n)||^{2}}$ is the Lorentz factor (see Fig. 1). We use the following abbreviated notations $\nu^{\hat{\alpha}}=\nu(U,n)^{\hat{\alpha}}$, $\nu=||\bm{\nu}(U,n)||\geq 0$, $\gamma(U,n)=\gamma$ throughout this paper. We have that $\nu$ represents the magnitude of the test particle spatial velocity $\bm{\nu}(U,n)$, $\alpha$ is the azimuthal angle of the vector $\bm{\nu}(U,n)$ measured clockwise from the positive $\hat{\varphi}$ direction in the $\hat{r}-\hat{\varphi}$ tangent plane in the ZAMO frame, and $\psi$ is the polar angle of the vector $\bm{\nu}(U,n)$ measured from the axis orthogonal to the $\hat{r}-\hat{\varphi}$ tangent plane of the ZAMO frame (see Fig. 1). The explicit expression for the test particle velocity components with respect to the ZAMOs are Bini et al. (2009, 2011a):

[TABLE]

where $\tau$ is the proper time parameter along $\mathbb{U}$.

2.4.1 Relativity of observer splitting formalism

The acceleration of the test particle relative to the ZAMO congruence, $\bm{a}(U)=\nabla_{\mathbb{U}}\mathbb{U}$, is given by the formula (see Eq. (29) in De Falco et al. (2018) and references therein) 222A complementary approach to the relativity of observer splitting formalism is the general relativistic Lagrangian formulation of the PR effect De Falco et al. (2018).:

[TABLE]

where $\alpha,\beta,\gamma=\hat{r},\hat{\theta},\hat{\varphi}$ run on the spatial indices of the metric coordinates 333 Terms $C_{\rm(Lie)}(n)^{\alpha}_{\beta\gamma},\ C_{\rm(Lie)}(n)^{\alpha}_{\beta}$, representing respectively the temporal and spatial constant structures, are missing in Eq. (28), because they vanish in a stationary and axially-symmetric spacetime (see Bini et al. (1997a, b); De Falco et al. (2018), for details).. Calculating the Christoffel symbols $\Gamma(n)^{\alpha}{}_{\beta\gamma}$, we have Bini et al. (1997a, b); De Falco et al. (2018)

[TABLE]

Therefore, Eqs. (28) in explicit form become

[TABLE]

From the orthogonality between $(\mathbb{a}(U),\mathbb{U})$ we have:

[TABLE]

2.5 Radiation test particle interaction

We assume that the interaction between the test particle and the radiation files takes place through Thomson scattering, characterized by a constant $\sigma$, independent of direction and frequency of the radiation field. The radiation force is Abramowicz et al. (1990); Bini et al. (2009, 2011a)

[TABLE]

where $P(U)^{\alpha}{}_{\beta}=\delta^{\alpha}_{\beta}+U^{\alpha}U_{\beta}$ projects a vector orthogonally to $\mathbb{U}$, namely on the spatial hypersurfaces or local rest spaces. The test particle equations of motion then become $m\mathbb{a}(U)=\bm{{\mathcal{F}}_{\rm(rad)}}(U)$, where $m$ is the test particle mass. By definition the radiation force lies in the local rest space of the test particle; to calculated it we decompose the photon four-momentum $\mathbb{k}$ first with respect to the four-velocity of the test particle, $\mathbb{U}$, and then to the previous ZAMO decomposition, $\mathbb{n}$, i.e. Bini et al. (2009, 2011a),

[TABLE]

By projecting $\mathbb{k}$ with respect to the test particle four-velocity, $\mathbb{U}$ we get

[TABLE]

Using Eq. (36) in Eq. (34) we obtain

[TABLE]

In this way the test particle acceleration is aligned with the photon relative velocity in the test particle local rest space, i.e.,

[TABLE]

where $\tilde{\sigma}=\sigma/m$. Hereafter we use the simplified notation $\bm{\hat{\mathcal{V}}}(k,U)=\bm{\hat{\mathcal{V}}}$. Multiplying scalarly Eq. (35) by $\mathbb{U}$ and using Eqs. (7) (i.e. the decomposition of $\mathbb{k}$ in the ZAMO frame), and (25) – (2.4) (i.e. the decomposition of $\mathbb{U}$ in the ZAMO frame), we find

[TABLE]

where we have used Eqs. (9) and (10) and the value of the assumed local angles. Such procedure is very useful for determining the spatial velocity $\bm{\hat{\mathcal{V}}}$:

[TABLE]

The frame components of $\bm{\hat{\mathcal{V}}}=\hat{\mathcal{V}}^{t}\bm{n}+\hat{\mathcal{V}}^{r}\bm{e_{\hat{r}}}+\hat{\mathcal{V}}^{\theta}\bm{e_{\hat{\theta}}}+\hat{\mathcal{V}}^{\varphi}\bm{e_{\hat{\varphi}}}$ are therefore

[TABLE]

where the second equality of Eq. (44) is due to the orthogonality of the $(\bm{\hat{\mathcal{V}}},\bm{U})$ pair and we have simplified the components of $\bm{\hat{\mathcal{V}}}$ of the radiation field.

2.5.1 General relativistic equations of motion

In order to make the equations of motion for the test particle moving in a 3D space explicit, Eqs. (38), we consider the ZAMO frame components of the test particle acceleration $\mathbb{a}(U)$, Eqs. (2.4.1) – (33), and the ZAMO frame components of the radiation force field $\bm{{\mathcal{F}}}_{\rm(rad)}(U)$, Eqs. (39), (41) – (44). The motion of the test particle is completely defined by the following six parameters $(r,\theta,\varphi,\nu,\psi,\alpha)$, the first three describing the position and the last three the velocity field. The displacement field is simply described by $(U^{r},U^{\theta},U^{\varphi})\equiv(dr/d\tau,d\theta/d\tau,d\varphi/d\tau)$. Instead the velocity field is connected to Eqs. (38) for determining $(d\nu/d\tau,d\psi/d\tau,d\alpha/d\tau)$. We note that using Eq. (33), it is possible to isolate $d\nu/d\tau$, indeed $a(U)^{\hat{t}}$ is the energy balance equation (see discussions in De Falco et al. (2018)). Then by using the expression of $d\nu/d\tau$ in $a(U)^{\hat{\theta}}$, Eq. (2.4.1), it is possible to determine $d\psi/d\tau$. Finally using the expressions of $d\nu/d\tau$ and $d\psi/d\tau$ in $a(U)^{\hat{r}}$, Eq. (2.4.1) yields $d\alpha/d\tau$.

Therefore, the general relativistic equations in the Kerr metric for the 3D motion of a test particle immersed in the radiation field defined in Sec. 2.3 and 2.5 are the following six coupled ordinary differential equations of the first order

[TABLE]

where $\tilde{\sigma}=\sigma/m$ and the two angles $\beta$ and $\zeta$ are calculated in terms of the two impact parameters $b$ and $q$. For $\psi=\zeta=\pi/2$ the equations of motion reduce to the 2D case Bini et al. (2009). Such set of equations reduce also to the classical 3D case in the weak field limit (see Appendix A).

Following Abramowicz et al. (1990); Bini et al. (2009, 2011a) we define the relative luminosity of the radiation field as

[TABLE]

Eq. (51) can be recast in the terms of the relative luminosity $A=L_{\infty}/L_{\rm EDD}$, taking thus the values in $[0,1]$, where $L_{\infty}$ is the luminosity of the central source as seen by an observer at infinity and $L_{\rm EDD}=4\pi Mm/\sigma$ is the Eddington luminosity at infinity. Then for the investigated radiation field with zero angular momenta ($b=0$, $\beta=\pi/2$) and without latitudinal photon motion ($q=-a^{2}\cos^{2}\theta$, $\zeta=\pi/2$), the term $\tilde{\sigma}[\Phi E(U)]^{2}$ becomes

[TABLE]

3 Critical hypersurface

The system of six differential equations (45) – (50) admits a critical solution of radial equilibrium, which corresponds to the axially-symmetric hypersurface where radiation pressure balances the attraction of the gravitational field. Let us consider a test particle moving purely radially with respect to the ZAMO frame ($\alpha=\psi=\pm\pi/2$). Then, at the critical radius $r_{\rm(crit)}$, where the test particle is in rest with respect to the ZAMO frame ($\nu=0,\gamma=1$), the first equation of motion, Eq. (45), takes the form

[TABLE]

In the case of pure radial motion ($\cos\alpha=0,\frac{d\alpha}{d\tau}=0$), the third equation of motion, Eq. (47), is automatically fulfilled. If we multiply the second equation of motion, Eq. (46), by the term $\nu^{2}$, thus removing its divergence, one can easily see that it is fulfilled in the radial equilibrium case ($\cos\psi=0,\frac{d\psi}{d\tau}=0,\nu=0$). For $\theta=\pi/2$ (i.e. in the equatorial plane) relation (53) corresponds to the equilibrium condition Eq. (2.33) derived in Bini et al. (2009), which gives the values $r_{\rm(crit)}$ of the radial coordinate where the test particle comoves with the ZAMOs in the equatorial circular orbit. However, relation (53) generalizes this condition also for the case of test particles with arbitrary polar angle $\theta$ and therefore describes a critical hypersurface which envelops the central compact object and where the test particles comoves with the local ZAMOs in a bound quasi-circular orbits. 444A different mechanism that leads to the formation of similar off-equatorial circular orbits is the interaction of charged test particles with the magnetic field of a neutron star Kovář et al. (2008).

In the case of a non-zero spin, the critical radius given by Eq. (53) is function of the polar angle of $r_{\rm(crit)}=r_{\rm(crit)}(A,\theta)$ (in addition to the relative luminosity $A$). The radial equilibrium therefore occurs at the axially symmetric hypersurface, whose shorter axis lies in the equatorial plane and longer axis in the polar direction. This is due to the properties of frame-dragging, as photons (and test particles) are dragged maximally in the azimuthal direction in the equatorial plane $\theta=\pi/2$. Therefore the radial component of the photon four-velocity reaches a maximum (and thus the radial momentum transfer is largest) along the polar axis and decreases for increasing polar angles; that is the reason why the critical hypersurface is elongated along the polar axis. In the case of zero spin (Schwarzschild spacetime), the critical hypersurface turns into a sphere with radius corresponding to the value given by Eq. (2.33) in Bini et al. (2009). The left panel of Fig. 2 compares the shape of the critical hypersurfaces for a high-spin Kerr spacetime $a=0.9995$ and for a Schwarzschild spacetime with $a=0$, where the relative luminosity of the radiating field is in both cases set to the value of $A=0.8$. In the high-spin case the critical radius is $r^{\rm eq}_{\rm(crit)}\sim 5.52M$ in the equatorial plane and $r^{\rm pole}_{\rm(crit)}\sim 6.56M$ at the poles. In the case of a Schwarzschild spacetime, the radius of the critical sphere is $r_{\rm(crit)}\sim 5.56M$. The right panel of Fig. 2 illustrates the shape of the critical hypersurfaces for the values of the relative luminosity in the interval $0.5-0.9$ and for a constant value of the spin $a=0.9995$.

4 Test particle orbits

We have developed the 3D PRtrajectories code to integrate the test particles trajectories described by equations (45) – (50). The integration of the equations of motion in three spatial dimensions turns out to be substantially more sensitive to integration errors than the the 2D case. Therefore we adapted the highly-accurate core for the integration of photon trajectories used in LSDCode+ Bakala et al. (2015) to the case of massive particles. The code implements the Runge-Kutta method of the eighth order (the Dorman – Prince method) Press et al. (2002) with an adaptive step. Successful integration of the 3D trajectory of test particles influenced by the radiation field (especially in the latitudinal direction) requires advanced monitoring of integration errors. In the 3D PRtrajectories code, the PI stepsize control algorithm (see Press et al. (2002) for details) is implemented, which easily attains an average relative accuracy of $\sim 10^{-14}$. Such a value allow precise and consistent integration of 3D trajectories even in the most sensitive parts, the vicinity of turning points.

We integrated equations (45) – (50) for a set of different boundary conditions and model parameters. Our results show that the main qualitative features of the 2D case examined in Bini et al. (2009, 2011a) remain the same for the trajectories in three spatial dimensions. Similarly to the 2D case, we can divide the orbits into two distinct classes depending on the initial radial position $r(0)=r_{0}$: inside and outside the critical hypersurface. Also in the 3D case, a test particle trajectory can have only two possible ends, either $(i)$ it goes to infinity or $(ii)$ it reaches the critical hypersurface. Moreover in the presence of (an outgoing) radiation field, the test particle cannot cross the event horizon. We compared representative trajectories of test particles with a polarly and azimuthally-oriented initial velocity for the case of a Schwarzschild spacetime, for the case of the Kerr metrics with a small spin (a=0.05) that approximates the spacetime in the vicinity of NSs and, finally, for the case of the Kerr metrics with very high spin ($a=0.9995$) that corresponds to the spacetime in the vicinity of almost extreme BHs.

In the case of the Schwarzschild metric our results fully agree with those from earlier analyses of the 2D case in which motions are confined to the equatorial plane Bini et al. (2009, 2011a) (note however, that this can be chosen arbitrarily for spherically symmetric metrics and radiation fields). The left panel of Fig. 3 shows examples of 3D trajectories which reflect the spherically-symmetric limit of equations (45 - 50) in the case of zero spin. The trajectories of test particles starting from the same location but with initial velocity oriented polarly and azimuthally are identical except for the different orientation of the plane on which they lie. In such a case, the 3D trajectories are easily transformed to the corresponding 2D trajectories through coordinate rotation. The left panel of Fig.3 then shows that in a Schwarzschild spacetime, once the trajectory of a test particle reaches the spherical critical hypersurface it stops precisely there.

Note that the presence of even a very small spin ($a=0.05$) breaks the spherical symmetry of the spacetime geometry and radiation field and introduces qualitatively-new features in test particle trajectories, owing to frame dragging effects. In the Kerr case, the test particle, once captured on the critical hypersurface, gets dragged azimuthally at the angular velocity $\Omega_{\mathrm{ZAMO}}$ and undergoes a latitudinal drift towards the equatorial plane (see Sec. 4.1, for a detailed explanation). Hence the test particle spirals on the critical hypersurface as it shifts to lower and lower latitudes. The results of our numerical integrations show that, besides the angular velocity $\Omega_{\mathrm{ZAMO}}$, the velocity of the latitudinal drift increases for increasing spins. In fact in for small spin values, test particles caught on the critical hypersurface encircle multiple spirals before attaining the final purely circular equatorial trajectory. Such behavior is exhibited by test particles whose motion starts from the inside (see right panel of Fig. 3 and left panel of Fig. 4), as well as the outside of the critical hypersurface (see Fig. 5 and right panel of Fig. 4). The behaviour of a test particles beginning its motion near the polar axis in a outgoing, purely radial direction is illustrated in Fig. 6: the test particle initially travels outward, reaches the turning point and then falls back along a nearly identical trajectory; after being captured near the pole of the critical hypersurface, it drifts toward the equator in tight spiralling trajectory that spans over most of northern hemisphere of the critical hypersurface.

For a value of the spin ($a=0.9995$) close to that of an extreme Kerr BH, frame dragging is faster and leads to a faster latitudinal drift, besides a higher $\Omega_{\mathrm{ZAMO}}$. Therefore test particles captured on the critical hypersurface at any value of the polar coordinate $\theta$ are dragged quickly to the equatorial plane where they attain a purely circular trajectory (see Figs 7 and 8).

4.1 Orbits bound to the critical hypersurface

In this section, we investigate in greater detail test particle trajectories bound to the critical hypersurface and their latitudinal drift towards the equatorial plane. We first emphasise that the condition for the radial balance of (outward) radiation force and gravitational attraction given by the equation (53) is satisfied also when the test particle reaches the critical hypersurface with non-zero angular momentum (with its space velocity vector thus forming an arbitrary angle $\alpha$ in the azimuthal direction; $\alpha\neq\pm\pi/2,\ \psi=\pm\pi/2,\ \nu=0,\ \gamma=1$). The orbits of the test particles, reaching the critical hypersurface, can be divided into two classes with qualitatively different behavior.

• (I)

Test particles with zero angular momentum achieve a complete balance of all forces acting at the critical hypersurface. Such case corresponds to test particle trajectories which satisfy the condition $\nu=0$ at any $r(0)=r_{0}$ and $\theta(0)=\theta_{0}$ and are thus carried around by frame dragging in the azimuthal direction, along with photons of the radiation field. At the critical hypersurface, such test particles then move along the purely off-equatorial circular orbits at constant latitude with angular velocity $\Omega_{\mathrm{ZAMO}}$ (see Fig. 9 and Fig. 10) remaining at rest relative to the appropriate ZAMO frame.

• (II)

Test particles that reach the critical hypersurface while still endowed with residual (non-zero) angular momentum (not coaligned with the spin axis, $\alpha\neq\pm\pi/2,\ \psi=\pm\pi/2,\theta\neq\pi/2$) attain radial balance, but the PR effect still operates on them because the radiation field is not yet directed in the radial direction in the test particle frame Bini et al. (2011a). Such particles exhibit a latitudinal drift on the critical hypersurface under the influence of the polar components of acceleration and consequently experience a polarly-oriented dissipative force originating from the interaction with the radiation field (see Eq. 46). In the latitudinal drift the residual angular momentum of the test particle is progressively removed. Then in accordance with the reflection symmetry of the Kerr spacetime, full equilibrium ($\alpha=\psi=\pi/2,\nu=0,\gamma=1$) is attained in the equatorial plane, where latitudinal drift stops, the motion stabilises in a circular orbit and the angular momentum of the test particle is completely removed (see right panel of Fig. 3, Figs 4 - 8 and Fig. 10).

We note that in the Schwarzschild case, the spin and polar acceleration are absent, and thus latitudinal drift does not occur (see left panel of Fig. 3).

5 Conclusions

We developed a fully general relativistic treatment of the 3D PR effect in the Kerr geometry, therefore extending previous works describing 2D PR motion in the equatorial plane of relativistic compact objects. The outgoing radiation field we adopted assumes that photons propagate radially with respect to the ZAMO frames. Such a boundary condition implies a purely radial propagation of the photons in any local ZAMO frame and may be considered as a simple approximation of the radiation field from a static emitting source very close to the horizon of a Kerr BH.

The resulting equations of motion for a test particle moving in the 3D space consist of a system of six coupled ordinary, highly nonlinear differential equations of first order. The non-linearity arises because of the general relativistic environment, further complicated by the PR effect which is a dissipative process and thus always entails nonlinearity. This set of equations is consistent with the previous 2D case for both test particles and photons moving in the equatorial plane Bini et al. (2009).

Our analytical and numerical calculations in both the Schwarzschild and Kerr metric, show that 3D PR orbits are strongly affected by general relativistic effects, including frame-dragging. We have demonstrate the existence of a critical hypersurface, where the attraction of gravity is balanced by the outgoing radiation forces. In the case of the Schwarzschild geometry, the critical hypersurface is a sphere, on which the test particles are captured and remain at rest. In the case of the Kerr spacetime (with non-zero spin), the critical hypersurface is elongated in the polar direction. Test particles that are captured by it are dragged at an azimuthal angular velocity $\Omega_{\mathrm{ZAMO}}$, and, if still endowed with a residual and offset angular momentum, they exhibit a latitudinal drift that lead to spiraling towards the equatorial plane. Analysis of the $\nu$ profile shows that the test test particle spatial velocity attain that of the local ZAMO in an infinite time (see Figs. 5 – 6 – 7). The test particle approaches the equatorial plane and radius of the hypersurface asymptotically. In future works we plan to relax some of the simple assumptions of the present study (e.g. by adopting more realistic radiation fields) and to investigate some possible astrophysical applications.

Acknowledgements

V.D.F. thanks the International Space Science Institute in Bern for support and the Max Planck Institute für Radioastronomie in Bonn for hospitality, since part of this work has been carried out there. V.D.F. thanks the Silesian University in Opava for partially funding this work. V.D.F. acknowledges useful discussions with Dr. Andrea Geralico on the relativity of observer splitting formalism. V.D.F. is grateful to Prof. Antonio Romano for the valuable discussions. P.B. and D.L. acknowledge the Czech Science Foundation (GAČR) grant GAČR 17-16287S. P.B. and E.B. thank the International Space Science Institute in Bern for the hospitality to carry out part of this work. D.L. acknowledges SU SGS/15/2016 and MSK 03788/2017/RRC grants. L.S. acknowledges financial contributions from ASI-INAF agreements 2017-14-H.O and I/037/12/0 and from iPeska research grant (P.I. Andrea Possenti) funded under the INAF call PRIN-SKA/CTA (resolution 70/2016).

Appendix A Classical 3D Poynting-Robertson effect

The classical radiation drag force was described and introduced by Poynting (1903) Poynting (1903) and Robertson (1937) Robertson (1937) in the 2D case. We extend the planar motion to the 3D case, written in spherical coordinates, $(r,\theta,\varphi)$. Noting that the classical drag force can be seen as a viscous effect depending linearly on the test particle velocity Poynting (1903); Robertson (1937); De Falco et al. (2018), and assuming that the radiation propagates radially in the whole 3D space, the test particle equations of motion read

[TABLE]

where the dot means the derivative with respect to the time, $G$ is the gravitational constant, $M$ the mass of the central object, $c$ the speed of the light, and $A=Sd^{2}/(6c^{2}\rho a)$ is the luminosity parameter with $S$ being the surface luminosity density of the compact object and $d$ the distance Earth – compact object. The term $-Ac/r^{2}$ represents the radiation pressure and $-2A\dot{r}/r^{2}$ is the specific angular momentum removed from the test particle due to the PR drag force.

1.1 Weak field approximation of the general relativistic equations

We show here the way in which the 3D general relativistic equations of motion $ma(U)^{\alpha}=F_{\rm(rad)}(U)^{\alpha}$, Eqs. (45) – (50), reduce to the classical 3D case, Eqs. (54) – (56), in the weak field limit ($a\to 0,\ r\to\infty,\ \nu/c\to 0$). Eqs. (48) – (50) are by definition

[TABLE]

The radial components of the ZAMO kinematical quantities reduce to

[TABLE]

expressed in geometrical units $G=c=1$, and where the relative Lie radial curvature reduces to the curvature of the osculating sphere (see Bini et al. (2009, 2011a); De Falco et al. (2018), for comparisons). Instead for the polar components of the ZAMO kinematical quantities we have

[TABLE]

where the relative Lie polar curvature describes the longitudinal Euler acceleration Bini et al. (1997a, b); De Falco et al. (2018). Now it is easy to see how the test particle acceleration, $a(U)^{\alpha}$, reduces to the left members of Eqs. (54) – (56). Approximating the radiation force, $F_{\rm(rad)}(U)^{\alpha}$, through linear terms in the velocity field, we have (see De Falco et al. (2018), for comparisons)

[TABLE]

which reduce to the right members of Eqs. (54) – (56). We note that the time component of the equations of motion, $ma(U)^{t}=F_{\rm(rad)}(U)^{t}$, reduces to De Falco et al. (2018)

[TABLE]

which represents the energy conservation equation. Indeed, the left term represents the total mechanical energy, while the right term corresponds to the dissipated energy.