hermitian matrix $\bm{M}_{X}{\hskip 0.85358pt}\bm{M}_{X}^{\dagger}$ by means
of the unitary transformation

$\displaystyle\bm{U}_{X}^{\dagger}{\hskip 0.85358pt}\bm{M}_{X}{\hskip
0.85358pt}\bm{M}_{X}^{\dagger}{\hskip
0.85358pt}\bm{U}_{X}=\text{diag}(M_{X}^{2},0,0,0)\,,$ (B.9)

where

$\displaystyle\bm{U}_{X}=\left(\begin{array}[]{ccccccc}c_{X}&&-s_{X}&&0&&0\\\\[2.84526pt]
s_{X}{\hskip 0.85358pt}\hat{\bm{Y}}_{X}&&c_{X}{\hskip
0.85358pt}\hat{\bm{Y}}_{X}&&\bm{n}_{X}^{1}&&\bm{n}_{X}^{2}\end{array}\right).$
(B.12)

We have defined the quantities

$\displaystyle c_{X}=\frac{m_{X}}{M_{X}}\,,\qquad
s_{X}=\frac{u}{\sqrt{2}{\hskip
0.85358pt}M_{X}}\big{(}\bm{Y}_{X}^{\dagger}\bm{Y}_{X}\big{)}^{1/2}\,,\qquad\hat{\bm{Y}}_{X}=\frac{\bm{Y}_{X}}{\big{(}\bm{Y}_{X}^{\dagger}\bm{Y}_{X}\big{)}^{1/2}}\,,$
(B.13)

which satisfy $c_{X}^{2}+s_{X}^{2}=1$ and $\hat{\bm{Y}}_{X}^{\dagger}{\hskip
0.85358pt}\hat{\bm{Y}}_{X}=1$. The complex unit vectors $\bm{n}_{X}^{i}$ with
$i=1,2$ are defined such that $\hat{\bm{Y}}_{X}^{\dagger}{\hskip
0.85358pt}\bm{n}_{X}^{i}=0$ and $\bm{n}_{X}^{i{\hskip
0.85358pt}\dagger}{\hskip 0.85358pt}\bm{n}_{X}^{j}=\delta^{ij}$. The three
unit vectors $\\{\hat{\bm{Y}}_{X},\bm{n}_{X}^{1},\bm{n}_{X}^{2}\\}$ form an
orthonormal basis in generation space. The mass eigenstates of the left-handed
$SU(2)_{L}$ doublets and right-handed $SU(2)_{L}$ singlets are related to the
interaction states by the matrices $\bm{U}_{Q}^{\dagger}$ and
$\bm{U}_{U,D}^{\dagger}$, respectively. Written out in components, this gives

$\displaystyle\left(\mathbb{Q}_{L}\right)_{\rm
mass}=c_{Q}\,\mathbb{Q}_{L}+s_{Q}{\hskip
0.85358pt}\hat{\bm{Y}}_{Q}^{\dagger}\,Q_{L}\,,\qquad\left(Q_{L}\right)_{\rm
mass}=\left(\begin{array}[]{c}c_{Q}{\hskip
0.85358pt}\hat{\bm{Y}}_{Q}^{\dagger}\,Q_{L}-s_{Q}\,\mathbb{Q}_{L}\\\
\bm{n}_{Q}^{1\,\dagger}\,Q_{L}\\\
\bm{n}_{Q}^{2\,\dagger}\,Q_{L}\end{array}\right),$ (B.17)

and similarly for the other cases. The second relation, in particular, implies
$\big{(}\hat{\bm{Y}}_{Q}^{\dagger}\,Q_{L}\big{)}_{\rm mass}=c_{Q}{\hskip
0.85358pt}\hat{\bm{Y}}_{Q}^{\dagger}\,Q_{L}-s_{Q}\,\mathbb{Q}_{L}$. It then
follows that

$\displaystyle\mathbb{Q}_{L}=c_{Q}\left(\mathbb{Q}_{L}\right)_{\rm
mass}-s_{Q}{\hskip
0.85358pt}\big{(}\hat{\bm{Y}}_{Q}^{\dagger}\,Q_{L}\big{)}_{\rm mass}\,.$
(B.18)

Analogous relations hold for the right-handed components $\mathbb{U}_{R}$ and
$\mathbb{D}_{R}$ of the $SU(2)_{L}$-singlet VLQs.

## References

  * [1] C. W. Bauer, S. Fleming, D. Pirjol and I. W. Stewart, Phys. Rev. D 63, 114020 (2001) [arXiv:hep-ph/0011336 [hep-ph]].
  * [2] C. W. Bauer, D. Pirjol and I. W. Stewart, Phys. Rev. D 65, 054022 (2002) [arXiv:hep-ph/0109045 [hep-ph]].
  * [3] C. W. Bauer, S. Fleming, D. Pirjol, I. Z. Rothstein and I. W. Stewart, Phys. Rev. D 66, 014017 (2002) [arXiv:hep-ph/0202088 [hep-ph]].
  * [4] M. Beneke, A. P. Chapovsky, M. Diehl and T. Feldmann, Nucl. Phys. B 643, 431 (2002) [hep-ph/0206152].
  * [5] S. Alte, M. König and M. Neubert, JHEP 08, 095 (2018) [arXiv:1806.01278 [hep-ph]].
  * [6] S. Alte, M. König and M. Neubert, Eur. Phys. J. C 79, no.4, 352 (2019) [arXiv:1902.04593 [hep-ph]].
  * [7] E. Eichten and B. R. Hill, Phys. Lett. B 234, 511-516 (1990).
  * [8] H. Georgi, Phys. Lett. B 240, 447-450 (1990).
  * [9] For a review, see: M. Neubert, Phys. Rept. 245, 259-396 (1994) [arXiv:hep-ph/9306320 [hep-ph]].
  * [10] P. Langacker, Rev. Mod. Phys. 81, 1199-1228 (2009) [arXiv:0801.1345 [hep-ph]].
  * [11] A. J. Buras and J. Girrbach, JHEP 12, 009 (2013) [arXiv:1309.2466 [hep-ph]].
  * [12] W. Altmannshofer, S. Gori, M. Pospelov and I. Yavin, Phys. Rev. D 89, 095033 (2014) [arXiv:1403.1269 [hep-ph]].
  * [13] L. Calibbi, A. Crivellin, F. Kirk, C. A. Manzari and L. Vernazza, Phys. Rev. D 101, no.9, 095003 (2020) [arXiv:1910.00014 [hep-ph]].
  * [14] L. Di Luzio, A. Greljo and M. Nardecchia, Phys. Rev. D 96, no.11, 115011 (2017) [arXiv:1708.08450 [hep-ph]].
  * [15] M. Bordone, C. Cornella, J. Fuentes-Martin and G. Isidori, Phys. Lett. B 779, 317-323 (2018) [arXiv:1712.01368 [hep-ph]].
  * [16] M. J. Baker, J. Fuentes-Martín, G. Isidori and M. König, Eur. Phys. J. C 79, no.4, 334 (2019) [arXiv:1901.10480 [hep-ph]].
  * [17] R. J. Hill and M. Neubert, Nucl. Phys. B 657, 229-256 (2003) [arXiv:hep-ph/0211018 [hep-ph]].
  * [18] R. J. Hill, T. Becher, S. J. Lee and M. Neubert, JHEP 07, 081 (2004) [arXiv:hep-ph/0404217 [hep-ph]].
  * [19] L. D. Landau, Dokl. Akad. Nauk SSSR 60, no.2, 207-209 (1948).
  * [20] C. N. Yang, Phys. Rev. 77, 242-245 (1950).
  * [21] C. Marcantonini and I. W. Stewart, Phys. Rev. D 79, 065028 (2009) [arXiv:0809.1093 [hep-ph]].
  * [22] W. Y. Keung, I. Low and J. Shu, Phys. Rev. Lett. 101, 091802 (2008) [arXiv:0806.2864 [hep-ph]].
  * [23] T. Becher and M. Neubert, Phys. Rev. Lett. 102, 162001 (2009) [Erratum: Phys. Rev. Lett. 111, no.19, 199905 (2013)] [arXiv:0901.0722 [hep-ph]].
  * [24] T. Becher and M. Neubert, JHEP 06, 081 (2009) [Erratum: JHEP 11, 024 (2013)] [arXiv:0903.1126 [hep-ph]].
  * [25] G. P. Korchemsky and A. V. Radyushkin, Nucl. Phys. B 283, 342-364 (1987).
  * [26] I. A. Korchemskaya and G. P. Korchemsky, Phys. Lett. B 287, 169-175 (1992).
  * [27] B. Jantzen, J. H. Kuhn, A. A. Penin and V. A. Smirnov, Nucl. Phys. B 731, 188-212 (2005) [Erratum: Nucl. Phys. B 752, 327-328 (2006)] [arXiv:hep-ph/0509157 [hep-ph]].
  * [28] L. N. Mihaila, J. Salomon and M. Steinhauser, Phys. Rev. D 86, 096008 (2012) [arXiv:1208.3357 [hep-ph]].
  * [29] B. Grzadkowski and M. Lindner, Phys. Lett. B 193, 71 (1987).
  * [30] T. Becher and M. Neubert, Phys. Rev. D 79, 125004 (2009) [Erratum: Phys. Rev. D 80, 109901 (2009)] [arXiv:0904.1021 [hep-ph]].
  * [31] T. Becher and M. Neubert, JHEP 01, 025 (2020) [arXiv:1908.11379 [hep-ph]].
  * [32] A. Ferroglia, M. Neubert, B. D. Pecjak and L. L. Yang, Phys. Rev. Lett. 103, 201601 (2009) [arXiv:0907.4791 [hep-ph]].
  * [33] A. Ferroglia, M. Neubert, B. D. Pecjak and L. L. Yang, JHEP 11, 062 (2009) [arXiv:0908.3676 [hep-ph]].
  * [34] J. Martin Camalich, M. Pospelov, P. N. H. Vuong, R. Ziegler and J. Zupan, Phys. Rev. D 102, no.1, 015023 (2020) [arXiv:2002.04623 [hep-ph]].


# Post-Newtonian properties of EMRI with Power Law Potential

Chinmay N. Gandevikar<EMAIL_ADDRESS>BITS Pilani K.K. Birla Goa
Campus, Sancoale, Goa 403726, India  Divyesh N. Solanki
<EMAIL_ADDRESS>Sardar Vallabhbhai National Institute of
Technology, Surat GUJ 395007, India  Dipanjan Dey<EMAIL_ADDRESS>International Center for Cosmology, Charusat University, Gujarat 388421, India

###### Abstract

There are many astrophysical scenarios where extreme mass ratio inspiral
(EMRI) binaries can be surrounded by matter distribution. The distribution of
mass can affect the dynamical properties (e.g. orbital frequency, average
energy radiation rate, etc.) of the EMRI. In this matter distribution, instead
of Kepler-Newton potential, one may consider a more general form of potential
i.e. power law potential. Moreover, due to the power law potential, the
velocity profile of test particles does not fall as much as that predicted by
Kepler-Newton potential and this feature of the velocity profile may be
observationally important. In this study, we have obtained the first post-
Newtonian (1PN) expressions for dynamical quantities and the average energy
radiation rate from the circular orbit EMRI which is surrounded by a matter
distribution. We show that the energy radiation rate and orbital frequency of
EMRI can be significantly different in the presence of power law potential as
compared to that in the Kepler-Newton potential, signatures of which may be
observed in gravitational waves from EMRI.

Keywords : Power law potential, Gravitational energy radiation, Post-Newtonian
approximation, EMRI.

##  I Introduction

In a galaxy, the study of dynamics of the bodies around the central object has
been of great interest for the last four decades. According to Kepler-Newton
potential, the velocities of the stars about a central supermassive object in
a galaxy should drop drastically with increasing radial distance. However,
observations suggest otherwise. The Kepler-Newton expression (also denoted as
KN potential) of gravitational potential ($\frac{Gm}{r}$) known to us can be
obtained from the effective potential of the Schwarzschild spacetime at an
infinite distance from the central body. However, the Schwarzschild metric is
a vacuum solution of the Einstein’s field equation. But at the galactic scale,
due to dark matter and baryonic matter, the system is non-vacuum. Due to this
limitation, the Schwarzschild metric or the KN potential cannot theoretically
explain galactic dynamics.

In the late 1970s and early 1980s, it was confirmed by the extensive study
reported in Sofue and Rubin in 2001 [1] that most of the stars in a galaxy
rotate at nearly the same velocities. This means that the velocities of these
stars are independent of their distance from the central mass. Investigation
of this novel observation has been a hot topic for several decades. There have
been multiple theories like modified Newtonian dynamic (MOND) [2], the
existence of dark matter, modified gravity, etc. [3, 4, 5]. In 1933, Zwicky
[6] suggested “missing mass” to account for the orbital velocities of galaxies
in clusters. Later, V Rubin and their group studied the possibility of dark
matter around galaxies [7, 8]. In [9, 10], authors study the mass density
distribution in galaxies and compare them with the observations. The density
profile of dark matter famously known as the Navvaro-Frenk-White (NFW) profile
[11] in the galaxy justifies the nature of the curve. Dey et al. have worked
on developing a general relativistic approach to understand galactic dynamics.
Bertrand Spacetimes is one of the non-vacuum spacetime which can be thought to
be seeded by dark matter [14, 15, 12, 13]. In the Newtonian limit, this
spacetime can effectively give the NFW dark matter density profile and hence
can explain the observed behaviour of rotational velocity curves. In 1983,
Milgrom suggested modifying Newtonian dynamics (MOND) which is an alternate
way to justify the curves. According to the MOND, Newton’s laws of gravity get
modified at a large distance scale. The theoretical prediction of MOND
satisfactorily fit the observed velocity profiles of stars for a large number
of galaxies [16]. Later, relativistic study of MOND was done by Bekenstein
[17], Moffat [19, 18] and Brownstein [20].

Logarithmic potential and appropriate power law potential can lead to the
rotational velocities (matching the observations) higher than that in the case
of KN potential. The velocity curve plotted using the logarithmic potential is
completely flat. Recently, Munera and Delgado-Correal [21] have shown
derivation of logarithmic potential as a solution of non-homogeneous Laplace
equation. Power law forces are mentioned in classical books like Goldstein
[22] and Danby [23]. The power law potentials (also denoted as PL potential)
have been extensively studied in the non-relativistic case. In [24, 26, 25],
the nature of orbits in the presence of power law potential is studied
extensively.

The KN potential can be written as an approximate potential derived from the
Schwarzschild solution of Einstein field equations which is only valid in
vacuum region. Hence, it is not entirely applicable to study binary systems in
galaxies or systems with dark matter distribution around them. Power law
potential can incorporate the effect of this mass distribution and hence also
can justify the behaviour of the velocity profile. In this paper, we consider
an extreme mass ratio inspiral (EMRI) binary of a stellar object which is
orbiting the central supermassive object in power law potential. The orbital
motion makes the EMRI radiate its energy which propagates as the gravitational
waves. We consider that this stellar object is in the faraway region from the
central supermassive object. This region allows the usability of Post-
Newtonian (PN) theory to obtain more precise expressions for the quantities
that can be described for the system, e.g., orbital velocity of the object,
the orbital radius, the acceleration of the object etc. Subsequently, we
predict the average energy radiation rate from the system. The upcoming
advanced gravitational-wave detectors like LISA may be capable of detecting
the gravitational radiations from EMRI binaries [27, 28, 29, 30].

We organize the paper in the following way. In Sec. (II), we derive the
acceleration of a body up to 1PN approximation using power law potential
instead of KN potential. The binary system and the PN corrected expressions
for the dynamical quantities are discussed in Sec. (III). In Sec. (IV), we
discuss the PN corrected average energy radiation rate. In Sec. (V), the
significant dynamical variables like the velocity and the orbital frequency,
and energy radiation in the case of power law potential are compared with that
in the case of KN potential system. In Sec. (VI), we discuss the important
results and possible future works.

##  II Acceleration of the reference body in the presence of power law
potential

The first order approximate expression for acceleration vector ($a_{P}^{j}$)
of a body named as body-P can be written as [31],

$\displaystyle
a^{j}_{P}=\partial^{j}U_{P}+\frac{1}{c^{2}}[(v_{P}^{2}-4U_{P})\partial^{j}U_{P}-4v_{P}^{j}v_{P}^{k}\partial_{k}U_{P}$
$\displaystyle-3v_{P}^{j}\dot{U_{P}}+\partial^{j}\psi_{P}+\frac{1}{2}\partial^{j}\ddot{X_{P}}+4\dot{U_{P}}^{j}$
$\displaystyle-4(\partial^{j}U_{P}^{k}-\partial^{k}U_{P}^{j})v_{P,k}]$
$\displaystyle+O(c^{-4}),$ (1)

where $U_{P}$, $U^{j}_{P}$, $\psi^{j}_{P}$ and $X_{P}$ are the external scalar
gravitational potential, vector gravitational potential, post-Newtonian
correction to scalar gravitational potential and the super-potential
respectively.

The components of the EMRI under consideration are a supermassive compact
object-$Q$ of mass $M_{Q}$ and a stellar mass object-$P$ of mass $M_{P}$
(which is much less than that of body-$Q$). We use the following notations
throughout the paper: $\mathbf{y}_{Q}\text{ and }\mathbf{y}_{P}$ are the
position vectors of the bodies Q and P in the center of mass frame, $r_{PQ}$
is the distance between the two bodies, also given as
$|\mathbf{y}_{P}-\mathbf{y}_{Q}|$. $\mathbf{v}_{Q}$ and $\mathbf{v}_{P}$ are
the absolute velocities of the bodies, $\mathbf{a}_{Q}$ and $\mathbf{a}_{P}$
are the absolute accelerations of the bodies, $\mathbf{x}_{S}$ is the position
vector of any field point ‘$S$’ where potential is being evaluated. The unit
vector in $\mathbf{x}_{SQ}$ direction and the unit vector in the direction of
$\mathbf{x}_{PQ}$ can be written as
$\mathbf{\hat{n}}_{SQ}=\frac{\mathbf{x_{SQ}}}{|\mathbf{x_{SQ}}|}=\frac{\mathbf{x}_{S}-\mathbf{y}_{Q}}{|\mathbf{x}_{S}-\mathbf{y}_{Q}|}$
and
$\mathbf{\hat{n}}_{PQ}=\frac{\mathbf{x_{PQ}}}{|\mathbf{x_{PQ}}|}=\frac{\mathbf{y}_{P}-\mathbf{y}_{Q}}{|\mathbf{y}_{P}-\mathbf{y}_{Q}|}$
respectively.

With increasing radial distance from the central supermassive object,
according to Kepler’s laws, the rotational velocities are expected to fall
down to nearly zero. However, in order to have higher values of rotational
velocities (which is generally observed in galactic scales [1, 8, 7]), one can
use specific forms of power law potential [33, 32]. A general form of the
expression of the Newtonian order (i.e. 0 PN order) gravitational acceleration
$\mathbf{a}_{P,N}(r)$ for the power law potential can be written as,

$\mathbf{a}_{P,N}(r_{PQ})=\partial^{j}U_{P}(r_{PQ})=-\frac{GM^{*}_{Q}r_{PQ}}{(\epsilon^{2}+r_{PQ}^{2})^{\delta+1}}\mathbf{\hat{n}}_{PQ}.$
(2)

Here, $M^{*}_{Q}$ is the constant scale mass of the potential due to body-Q,
$\epsilon$ is the scale radius of the source and $\delta$ can have any real
value [34]. Note that the subscript in the acceleration is “$P,N$” denoting
acceleration of body-P at Newtonian order, which is different from “$PN$”
which denotes the post-Newtonian corrections. Two interesting special cases
are

  1. 1.

The circular velocity is constant ($v_{circ}=\sqrt{|\mathbf{a}|r}=$ const)
i.e. $\delta=0$ and

  2. 2.

The circular velocity is inversely proportional to the square root of the
radial distance ($v_{circ}=\sqrt{|\mathbf{a}|r}\propto r^{-1/2}$) i.e.
$\delta=\frac{1}{2}$.

Above two cases correspond to flat and Keplerian rotation curves,
respectively. The exactly flat rotation curves can be obtained at $\delta=0$.
Using Eq. (2), it can be seen that the corresponding potential is the
logarithmic potential. In this study, the general power law potential is taken
under consideration except for the logarithmic potential case. The potential
($U_{P}$) corresponding to the above mentioned acceleration (Eq. (2)) at point
P due to super-massive body-Q [33] is

$U_{P}=\frac{GM^{*}_{Q}}{2\delta(\epsilon^{2}+r_{PQ}^{2})^{\delta}}.$ (3)

Corresponding vector potential ($U_{P}^{j}$) [31] is accordingly obtained

$U^{j}_{P}=\frac{GM^{*}_{Q}}{2\delta(\epsilon^{2}+r_{PQ}^{2})^{\delta}}{v}^{j}_{Q}.$
(4)

###  II.1 External potentials and their derivatives

In this sub-section, we derive the post-Newtonian potentials and their
derivatives as introduced in Eq. (II). It is considered that the distance
between the two bodies is much larger than the sizes of the two bodies. From
here on, external potentials are referred to as just “potentials”.

Before introducing other external potentials and their derivatives, we
introduce the following tools which are used to obtain the same

$\displaystyle\partial^{j}r_{SQ}$ $\displaystyle=$ $\displaystyle n_{SQ}^{j},$
(5) $\displaystyle\partial_{j}\mathbf{x}^{k}_{SQ}$ $\displaystyle=$
$\displaystyle\delta_{j}^{k},$ (6) $\displaystyle\dot{\mathbf{x}}_{SQ}$
$\displaystyle=$ $\displaystyle-\mathbf{v}_{Q},$ (7)
$\displaystyle\dot{{r}}_{SQ}$ $\displaystyle=$
$\displaystyle-\mathbf{\hat{n}}_{SQ}.\mathbf{v}_{Q},$ (8)
$\displaystyle\partial^{j}n_{SQ}^{k}$ $\displaystyle=$
$\displaystyle\frac{\delta^{jk}}{r_{SQ}}-\frac{n_{SQ}^{j}n_{SQ}^{k}}{r_{SQ}},$
(9) $\displaystyle\partial^{j}v_{Q}^{k}$ $\displaystyle=$ $\displaystyle 0.$
(10)

The over-dot represents the time derivative of the quantity. Following are the
expressions for the derivatives of the two potentials mentioned in Eq. (3) and
Eq. (4) [31].

$\displaystyle\dot{U}_{P}$ $\displaystyle=$
$\displaystyle\frac{GM^{*}_{Q}r_{PQ}}{(r_{PQ}^{2}+\epsilon^{2})^{\delta+1}}(\mathbf{\hat{n}}_{PQ}.\mathbf{v}_{Q}),$
(11) $\displaystyle\dot{U}^{j}_{P}$ $\displaystyle=$
$\displaystyle\frac{GM^{*}_{Q}r_{PQ}}{2\delta(r_{PQ}^{2}+\epsilon^{2})^{\delta+1}}\bigg{[}2\delta(\mathbf{\hat{n}}_{PQ}.\mathbf{v}_{Q})v_{Q}^{j}$
(12) $\displaystyle\text{
}+\frac{GM_{P}^{*}}{(r_{PQ}^{2}+\epsilon^{2})^{\delta}}{n}^{j}_{PQ}\bigg{]},$
$\displaystyle\partial^{j}U_{P}$ $\displaystyle=$
$\displaystyle-\frac{GM_{Q}^{*}r_{PQ}}{(r_{PQ}^{2}+\epsilon^{2})^{\delta+1}}n^{j}_{PQ},$
(13) $\displaystyle\partial^{j}U^{k}_{P}$ $\displaystyle=$
$\displaystyle-\frac{GM^{*}_{Q}r_{PQ}}{(r_{PQ}^{2}+\epsilon^{2})^{\delta+1}}{n}^{j}_{PQ}{v}^{k}_{Q}.$
(14)

Post-Newtonian correction ($\psi_{P}$) to the Newtonian Potential.

$\psi_{P}=U_{P}\left(\frac{3}{2}v_{Q}^{2}-U_{Q}\right)\,\,,$ (15)

where $U_{Q}$ is the potential at body-Q. Spatial derivative of $\psi_{P}$ is

$\displaystyle\partial^{j}\psi_{P}$ $\displaystyle=$
$\displaystyle\partial^{j}\left[U_{P}\left(\frac{3}{2}v_{Q}^{2}-U_{Q}\right)\right]$
(16) $\displaystyle=$
$\displaystyle-\frac{3}{2}\frac{GM^{*}_{Q}r_{PQ}}{(r_{PQ}^{2}+\epsilon^{2})^{\delta+1}}v_{Q}^{2}n_{PQ}^{j}.$

Mathematically, potential is the quantity obtained from Poisson’s equation
using a source term. We need mass density to get the gravitational potential
($U$). Using this Newtonian gravitational potential as source, a higher order
of potential called super-potential ($X$) can be obtained. In other words
super-potential is the potential generated by external potential. Here we will
briefly explain the derivation:
The Poisson’s equation for super-potential ($X$) is given as

$\nabla^{2}X=2U.$ (17)

The source term contains the potential $U$. The potential due to a mass-$Q$ at
a field point ($\mathbf{x^{\prime}}$) is a function of the distance between
the field point and the source mass is

$U(\mathbf{x^{\prime}})=U(\mathbf{x^{\prime}}-\mathbf{y}_{Q}).$ (18)

The corresponding super-potential ($X(\mathbf{x}_{S})$) at another field point
($\mathbf{x}_{S}$) due to the potential $U(\mathbf{x^{\prime}})$ is the
solution of the Poisson’s equation,
$\nabla^{2}X(\mathbf{x}_{S})=2U(\mathbf{x}_{S}-\mathbf{x^{\prime}})$.
$X(\mathbf{x}_{S})$ is spherically symmetric hence it’s dependence on vector
($\mathbf{x}_{S}-\mathbf{x^{\prime}}$) can be simplified into the dependence
over the magnitude of this vector, which we denote as ‘$r$’. The derivative of
this super-potential as required in Eq. (II) is formalized as

$\displaystyle\partial^{j}\ddot{X}(\mathbf{x_{S}})\ $ $\displaystyle=$
$\displaystyle\partial^{j}(\partial_{tt}X(r))$ $\displaystyle=$
$\displaystyle\partial^{j}\left(X^{\prime\prime}(r)\dot{r}^{2}+X^{\prime}(r)\ddot{r}\right)$
$\displaystyle=$ $\displaystyle
X^{\prime\prime\prime}(r)\dot{r}^{2}\partial^{j}r+2X^{\prime\prime}_{P}(r)\dot{r}\partial^{j}\dot{r}$
(19)
$\displaystyle+X^{\prime\prime}(r)\ddot{r}\partial^{j}r+X^{\prime}(r)\partial^{j}\ddot{r}.$

Here, the primes denote radial derivative and the over-dots are for time
derivatives. The general solution to the Poisson’s equation is given as

$X(\mathbf{x}_{S})=-\frac{1}{2\pi}\int\frac{U(\mathbf{x}^{\prime})}{|\mathbf{x}_{S}-\mathbf{x^{\prime}}|}d^{3}x^{\prime}+X_{0}(\mathbf{x}_{S}).$
(20)

Here, $X_{0}$ is the constant of integral, which is the solution to the
Laplace equation $\nabla^{2}X_{0}=0$. As prescribed [31] (in the box 7.3) the
domain of integration is truncated to spatial near zone. This changes the form
of Poisson’s equation to $\nabla^{2}X=2U\Theta(R-r)$, where $\Theta(R-r)$ is
the Heaviside step function and $R$ is the radius of the boundary of the near
zone (value of which is immaterial and this fact will be used later).

According to our assumption, the body-Q is heavy enough compared to the body-P
such that the effect due to the body-P is negligible, and we consider only the
effect due to body-Q as

$U(\mathbf{x^{\prime}})\approx
U(\mathbf{x^{\prime}}-\mathbf{y}_{Q})=\frac{GM^{*}_{Q}}{2\delta(r_{1Q}^{2}+\epsilon^{2})^{\delta}},$
(21)

here, $r_{1Q}=|\mathbf{x^{\prime}}-\mathbf{y}_{Q}|$.
Substituting Eq. (21) into Eq. (20), we obtain

$X(\mathbf{x}_{S})=-\dfrac{GM^{*}_{Q}}{4\pi\delta}\int\frac{\Theta(R-r^{\prime})}{|\mathbf{x}_{S}-\mathbf{x^{\prime}}|(r^{2}_{1Q}+\epsilon^{2})^{\delta}}d^{3}x^{\prime}+X_{0}(\mathbf{x}_{S}).$
(22)

From here on, $r^{\prime}=|\mathbf{x^{\prime}}|$. If the distance between
points, $r_{1Q}$, is very large compared to the scale length $\epsilon$, i.e.
$\epsilon<<r_{1Q}$, we can approximate the potential as

$U(\mathbf{x^{\prime}-y}_{Q})\approx\frac{GM^{*}_{Q}}{2\delta(r_{1Q})^{2\delta}}.$
(23)

Therefore,

$X(\mathbf{x}_{S})=-\dfrac{GM^{*}_{Q}}{4\pi\delta}\int\frac{\Theta(R-r^{\prime})}{|\mathbf{x}_{S}-\mathbf{x^{\prime}}|(r_{1Q})^{2\delta}}d^{3}x^{\prime}+X_{0}(\mathbf{x}_{S}).$
(24)

For the mathematical simplicity, we take $\mathbf{y}_{Q}=0$.
At the end of this section, we generalize our result by introducing
$\mathbf{y}_{Q}$. We set the Cartesian coordinates such that the Z-axis aligns
with the position vector $\mathbf{x_{S}}$.
Therefore, $\mathbf{x_{S}}=r\mathbf{\hat{z}}$ ,
$\mathbf{x_{S}}.\mathbf{x^{\prime}}=rr^{\prime}\cos\theta^{\prime}$ , and
$|\mathbf{x}_{S}-\mathbf{x^{\prime}}|=\sqrt{r^{2}+r^{\prime
2}-2rr^{\prime}\cos\theta^{\prime}}$. Hence

$X(\mathbf{x}_{S})=-\dfrac{GM^{*}_{Q}}{4\pi\delta}\int\frac{\Theta(R-r^{\prime})r^{\prime
2}\sin\theta^{\prime}dr^{\prime}d\theta^{\prime}d\phi^{\prime}}{r^{\prime
2\delta}\sqrt{r^{2}+r^{\prime 2}-2rr^{\prime}\cos\theta^{\prime}}}+X_{0}.$
(25)

The integrand is independent of $\phi^{\prime}$, thus the $\phi$ integral will
give a factor $2\pi$. Now, we integrate $\theta^{\prime}$ part first. Let us
define it as $I(r,r^{\prime})$.

$\displaystyle I(r,r^{\prime})$ $\displaystyle=$
$\displaystyle\int_{0}^{\pi}(r^{2}+r^{\prime
2}-2rr^{\prime}\cos\theta^{\prime})^{-1/2}\sin\theta^{\prime}d\theta^{\prime}$
(26) $\displaystyle=$
$\displaystyle\dfrac{1}{rr^{\prime}}\\{(r+r^{\prime})-|r-r^{\prime}|\\}$ (27)

Now we integrate over $r^{\prime}$.

$\displaystyle X(\mathbf{x}_{S})$ $\displaystyle=$
$\displaystyle-\dfrac{GM^{*}_{Q}}{2\delta}\int_{0}^{R}r^{\prime
2(1-\delta)}I(r,r^{\prime})dr^{\prime}+X_{0}$ (28) $\displaystyle=$
$\displaystyle-\dfrac{GM^{*}_{Q}}{2\delta}\bigg{(}\int_{0}^{r}r^{\prime
2(1-\delta)}\dfrac{2}{r}dr^{\prime}+\int_{r}^{R}r^{\prime
2(1-\delta)}\dfrac{2}{r^{\prime}}dr^{\prime}\bigg{)}+X_{0}$ $\displaystyle=$
$\displaystyle\dfrac{GM^{*}_{Q}}{2\delta(1-\delta)}\bigg{(}\dfrac{r^{2(1-\delta)}}{3-2\delta}-R^{2(1-\delta)}\bigg{)}+X_{0}.$

Now, one should consider
$X_{0}=\dfrac{GM^{*}_{Q}}{2\delta(1-\delta)}R^{2(1-\delta)}$ to get rid of the
terms dependent on the boundary of near zone. One can get the following
expression of the super-potential

$X(\mathbf{x}_{S})=\dfrac{GM^{*}_{Q}}{2\delta(1-\delta)(3-2\delta)}(r_{SQ})^{2(1-\delta)}\,\,,$
(29)

where we replace $\mathbf{x}$ by $\mathbf{x-y}_{Q}$ to generalize the result
(Eq. (28)) and we drop the assumption $\mathbf{y}_{Q}=0$. Placing the body-P
at point S, it can be shown that

$\partial^{j}\ddot{X}_{P}=\dfrac{2GM_{Q}^{*}}{(3-2\delta)(r_{PQ})^{2\delta+1}}\Bigg{[}2(\delta+1)(\mathbf{v}_{Q}.\mathbf{\hat{n}}_{PQ})^{2}n^{j}_{PQ}-v_{Q}^{2}n^{j}_{PQ}-2(\mathbf{v}_{Q}.\mathbf{\hat{n}}_{PQ})v_{Q}^{j}+\dfrac{2GM^{*}_{P}}{(r_{PQ})^{2\delta}}n^{j}_{PQ}\Bigg{]}.$
(30)

Acceleration of body-P can be obtained by substituting the expressions for the
derivatives of the potentials as derived in Eq. (11), (12), (13), (14), (16),
and (30) in Eq. (II). As we previously mentioned, we assume that the body-Q is
much heavier than the body-P. Therefore, for mathematical simplicity, we can
consider the body-Q is always stationary which implies $\mathbf{v}_{Q}\approx
0$ and $\mathbf{a}_{Q}\approx 0$. For simplicity in notations, we use
$r_{PQ}=r$, $\mathbf{\hat{n}}_{PQ}=\mathbf{\hat{n}}$,
$\mathbf{v}_{P}=\mathbf{v}$, and $\dfrac{GM_{Q}^{*}}{2}=K_{Q}$. Therefore, the
potential and it’s differentials reduce to the expressions listed as

$\displaystyle U_{P}$ $\displaystyle=$ $\displaystyle\dfrac{K_{Q}}{\delta
r^{2\delta}},$ (31) $\displaystyle\partial^{j}U_{P}$ $\displaystyle=$
$\displaystyle\dfrac{-2K_{Q}}{r^{2\delta+1}}n^{j},$ (32) $\displaystyle
U_{P}^{j}$ $\displaystyle=$ $\displaystyle\dfrac{K_{Q}v_{Q}^{j}}{\delta
r^{2\delta}}=0,$ (33) $\displaystyle\dot{U}_{P}$ $\displaystyle=$
$\displaystyle\dfrac{2K_{Q}(\mathbf{v}_{Q}.\mathbf{\hat{n}})}{r^{2\delta+1}}=0,$
(34) $\displaystyle\partial^{j}\psi_{P}$ $\displaystyle=$
$\displaystyle\dfrac{-3K_{Q}v_{Q}^{2}}{r^{2\delta+1}}n^{j}=0,$ (35)
$\displaystyle\dot{U}_{P}^{j}$ $\displaystyle=$
$\displaystyle\dfrac{2K_{P}K_{Q}}{\delta
r^{4\delta+1}}n^{j}+\dfrac{2K_{Q}(\mathbf{v}_{Q}.\mathbf{\hat{n}})}{r^{2\delta+1}}v_{Q}^{j}$
(36) $\displaystyle=$ $\displaystyle\dfrac{2K_{P}K_{Q}}{\delta
r^{4\delta+1}}n^{j},$ $\displaystyle\partial^{j}U_{P}^{i}$ $\displaystyle=$
$\displaystyle\dfrac{-2K_{Q}}{r^{2\delta+1}}n^{j}v_{Q}^{i}=0.$ (37)

One can substitute the above expressions (Eqs. (31) to (37)) along with the
expression for the derivative of the super-potential as given in Eq. (30) into
Eq. (II) to obtain the acceleration ($a_{P}^{j}$) as follows.

$a_{P}^{j}=-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{r^{2\delta+1}}n^{j}-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{c^{2}r^{2\delta+1}}\bigg{[}-4(\mathbf{v}_{P}.\mathbf{\hat{n}})v_{P}^{j}+\bigg{(}v_{P}^{2}-\dfrac{2GM_{Q}\epsilon^{2\delta-1}}{\delta
r^{2\delta}}\bigg{)}n^{j}\bigg{]},$ (38)

where the mass of the body-Q ($M_{Q}$) is related to the scale mass
$M_{Q}^{*}$ by $M^{*}_{Q}=M_{Q}\epsilon^{2\delta-1}$. In the next section, we
derive the dynamical variables corresponding to the binary.

##  III THE BINARY SYSTEM AND THE DYNAMICAL QUANTITIES

The EMRI under consideration comprises of a supermassive compact object
(body-Q) and a stellar mass object (body-P) moving in a circular orbit around
body-Q. The position $\mathbf{x}$ of the body-P is written as

$\mathbf{x}=r(\cos{\phi}\mathbf{\hat{i}}+\sin{\phi}\mathbf{\hat{j}}),$ (39)

where $\phi$ is the azimuthal angle, $r$ is the radius. Since we don’t need to
work with the dynamics of any body other than body-P, we drop the
subscript-“P”, i.e. the quantities $r_{P}$, $\mathbf{x}_{P}$, $\mathbf{v}_{P}$
and $\mathbf{a}_{P}$ will be denoted as $r$, $\mathbf{x}$, $\mathbf{v}$ and
$\mathbf{a}$ and the subscript “PQ” in the quantities that are measured with
reference to body-Q are also dropped. i.e. $r_{PQ}$, $\mathbf{x}_{PQ}$ and
$\mathbf{\hat{n}}_{PQ}$ are henceforth denoted as $r$, $\mathbf{x}$,
$\mathbf{\hat{n}}$. In this study, we have used the following sign convention:
$U(r)=+\dfrac{GM_{Q}\epsilon^{2\delta-1}}{\delta r^{2\delta}}$, and
$a(r)=+\dfrac{dU}{dr}$. The differential orbit equation (for the body-P) can
be written as [22, 32, 33, 35]

$\dfrac{d^{2}u}{d\phi^{2}}+u=\dfrac{M_{P}}{L^{2}}\dfrac{dU(1/u)}{du},$ (40)

where $u=1/r$, $L$ is conserved angular momentum of the body-P. Therefore

$\dfrac{d^{2}u}{d\phi^{2}}+u=\dfrac{GM_{Q}M_{P}}{L^{2}}\epsilon^{2\delta-1}u^{2\delta-1}.$
(41)

For the circular orbits, the solution of the above differential orbit equation
is $u=1/p$. Therefore, we obtain

$p^{2(\delta-1)}=\dfrac{GM_{Q}M_{P}}{L^{2}}\epsilon^{2\delta-1}.$ (42)

###  III.1 Dynamical quantities at Newtonian order

Henceforth, the dynamical quantities, say F, will appear in the form of
$F=F_{N}+\frac{1}{c^{2}}F_{PN}$, where $F_{N}$ denotes the Newtonian order or
the 0-PN term which has been derived in this sub-section (III.1) and $F_{PN}$
denotes the 1PN correction to the quantity F which will be derived in the next
sub-section (III.2). For example, the Newtonian order and post-Newtonian order
expression for velocity of body-P are denoted as $\mathbf{v}_{N}$ and
$\mathbf{v}_{PN}$ respectively while the complete PN-corrected term is denoted
as $\mathbf{v}$.

We derive the Kepler’s third law for the power law potential by balancing the
centrifugal acceleration with the radial attractive acceleration due to
gravity,

$\dfrac{v_{N}^{2}}{r}=\dfrac{GM_{Q}}{r^{2\delta+1}}\epsilon^{2\delta-1}.$ (43)

Here, $v_{N}=\dfrac{2\pi r}{T}$ is the 0PN or the Newtonian order expression
of the velocity of the body-P and $T$ is the time period of its orbit.
Therefore,

$\dfrac{4\pi^{2}}{T^{2}}=\dfrac{GM_{Q}}{r^{2(\delta+1)}}\epsilon^{2\delta-1}.$
(44)

The Kepler’s third law in the case of power law potential becomes
$T^{2}\propto r^{2(\delta+1)}$. We have restricted our attention to the
circular orbits, therefore $r=p$ (constant). Thus the angular velocity
$\dot{\phi}_{N}$ becomes

$\dot{\phi_{N}}=\dfrac{2\pi}{T}=\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2(\delta+1)}}}.$
(45)

The 0PN velocity $\mathbf{v}_{N}$ of the body-P using Eq. (43) becomes

$\mathbf{v}_{N}=\dfrac{d\mathbf{x}}{dt}=\dfrac{d\mathbf{x}}{d\phi}\dot{\phi}=p\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2(\delta+1)}}}\bm{\hat{\nu}}$
(46)

where
$\bm{\hat{\nu}}=\dfrac{d\mathbf{\hat{n}}}{d\phi}=-\sin{\phi}\mathbf{\hat{i}}+\cos{\phi}\mathbf{\hat{j}}$.
The Newtonian order term of the acceleration is

$\mathbf{a}_{N}=\dfrac{\partial U}{\partial
r}\bigg{|}_{r=p}=-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta+1}}\mathbf{\hat{n}}\,\,.$
(47)

###  III.2 Dynamical quantities with Post-Newtonian correction

We use the 1PN corrected acceleration (38) to determine the 1PN correction in
other dynamical variables, i.e. velocity, radial coordinate and angular
frequency. In the case of circular orbits,
$\mathbf{v}_{P}.\mathbf{\hat{n}}=0$, and the acceleration (38) reduces to

$\mathbf{a}=-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta+1}}\mathbf{\hat{n}}+\dfrac{1}{c^{2}}\dfrac{G^{2}M_{Q}^{2}\epsilon^{2(2\delta-1)}}{\delta
p^{4\delta+1}}(2-\delta)\mathbf{\hat{n}},$ (48)

which can be written as

$\mathbf{a}=\mathbf{a}_{N}+\dfrac{1}{c^{2}}\mathbf{a}_{PN},$ (49)

where the 1PN correction in the acceleration is

$\mathbf{a}_{PN}=\dfrac{G^{2}M_{Q}^{2}\epsilon^{2(2\delta-1)}}{\delta
p^{4\delta+1}}(2-\delta)\mathbf{\hat{n}}.$ (50)

The expression for velocity at the Newtonian order (0PN) can be given as
$v_{N}^{2}=\dfrac{GM_{Q}}{p^{2\delta}}\epsilon^{2\delta-1}$ as derived in Eq.
(46). After substituting this expression, the post-Newtonian parameter
($v_{N}^{2}/c^{2}$) becomes evidently apparent. Hence we can simplify Eq.
(48),

$\mathbf{a}=-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta+1}}\bigg{[}1-\dfrac{v_{N}^{2}}{c^{2}}\dfrac{(2-\delta)}{\delta}\bigg{]}\mathbf{\hat{n}}.$
(51)

In Kepler-Newton (KN) case, $\delta=1/2$, and the acceleration
($\mathbf{a}_{KN}$) reduces to ([36], [37])

$\mathbf{a}_{KN}=-\dfrac{GM_{Q}}{p^{2}}\bigg{[}1-3\dfrac{v_{KN}^{2}}{c^{2}}\bigg{]}\mathbf{\hat{n}}\,\,,$
(52)

where $v_{KN}^{2}=\dfrac{GM_{Q}}{p}$.

Post-Newtonian corrections are only valid as long as $(\frac{v_{N}}{c})^{2}\ll
1$ and the PN correction is much smaller than the Newtonian order part. This
can be simplified to give the following condition

$1\gg\big{(}\frac{v_{N}}{c}\big{)}^{2}\big{(}\frac{2-\delta}{\delta}\big{)}=\frac{GM_{Q}\epsilon^{(2\delta-1)}}{p^{2\delta}c^{2}}\frac{2-\delta}{\delta}\implies\delta>0.$
(53)

We shall consider the cases where $\delta>0$.

Now, we proceed to determine the PN corrections in the other dynamical
variables. Let the PN corrected angular momentum per unit mass
($r^{2}\dot{\phi}$) and PN corrected velocity ($\mathbf{v}_{P}$) be written as

$\displaystyle r^{2}\dot{\phi}$ $\displaystyle=$
$\displaystyle|\mathbf{r}\times\mathbf{v}|=p\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}(1+\delta
h)\,,$ (54) $\displaystyle\mathbf{v}$ $\displaystyle=$
$\displaystyle\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\bm{\hat{\nu}}+\delta\mathbf{v}\,,$
(55)

where $p\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\delta h$ and
$\delta\mathbf{v}$ are the PN correction terms. As in this paper, we consider
circular orbits, for a particular orbit, the angular momentum of the body is
constant. Therefore, Eq. (54) becomes [37],

$r^{2}\dot{\phi}=p\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}},$
(56)

here $r$ and $\dot{\phi}$ are corrected up to 1-PN.

Since circular orbits are considered, the radial component of the velocity is
zero. Hence the velocity can be written as $v=r\dot{\phi}$. Also the relation
between acceleration $\mathbf{a}$ and the velocity $\mathbf{v}$ can be given
as:

$\mathbf{a}=\dfrac{d\mathbf{v}}{dt}=\dfrac{d\mathbf{v}}{d\phi}\dfrac{d\phi}{dt}=\dfrac{d\mathbf{v}}{d\phi}\dot{\phi}.$
(57)

Using Eq. (55) and Eq. (51), the above equation can be re-written as,

$\displaystyle-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta+1}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\frac{(2-\delta)}{\delta}\bigg{]}\mathbf{\hat{n}}$
$\displaystyle=$
$\displaystyle\bigg{(}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\dfrac{d\bm{\hat{\nu}}}{d\phi}+\dfrac{d\delta\mathbf{v}}{d\phi}\bigg{)}\dot{\phi}$
(58) $\displaystyle=$
$\displaystyle\bigg{(}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\dfrac{d\bm{\hat{\nu}}}{d\phi}+\dfrac{d\delta\mathbf{v}}{d\phi}\bigg{)}\dfrac{1}{p}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}$
$\displaystyle=$
$\displaystyle-\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta+1}}\mathbf{\hat{n}}+\dfrac{1}{p}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\dfrac{d\delta\mathbf{v}}{d\phi}.$

Equating the 1-PN terms, we get

$\displaystyle\dfrac{2-\delta}{c^{2}}\dfrac{G^{2}M_{Q}^{2}\epsilon^{2(2\delta-1)}}{\delta
p^{4\delta+1}}\mathbf{\hat{n}}$ $\displaystyle=$
$\displaystyle\dfrac{1}{p}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\dfrac{d\delta\mathbf{v}}{d\phi}.$
(59)

Therefore,

$\displaystyle\dfrac{d\delta\mathbf{v}}{d\phi}$ $\displaystyle=$
$\displaystyle\dfrac{1}{c^{2}}\dfrac{(2-\delta)}{\delta}\bigg{[}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\bigg{]}^{3/2}\mathbf{\hat{n}},$
(60) $\displaystyle\therefore\delta\mathbf{v}$ $\displaystyle=$
$\displaystyle-\dfrac{1}{c^{2}}\dfrac{(2-\delta)}{\delta}\bigg{[}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\bigg{]}^{3/2}\bm{\hat{\nu}}.$
(61)

Hence, the velocity $\mathbf{v}_{P}$ of the body-P corrected up to the first
order Post-Newtonian term is

$\therefore\mathbf{v}_{P}=\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\frac{(2-\delta)}{\delta}\bigg{]}\bm{\hat{\nu}}.$
(62)

Substituting Eq. (62) in the following relation

$\mathbf{v}_{P}=\dfrac{d\mathbf{x}}{dt}=\dfrac{d\mathbf{x}}{d\phi}\dfrac{d\phi}{dt}=r\dot{\phi}\bm{\hat{\nu}}=\dfrac{r^{2}\dot{\phi}}{r}\bm{\hat{\nu}},$
(63)

and using Eq. (56), we obtain

$\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\frac{(2-\delta)}{\delta}\bigg{]}\bm{\hat{\nu}}\\\
=\dfrac{r^{2}\dot{\phi}}{r}\bm{\hat{\nu}}=\dfrac{p}{r}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\bm{\hat{\nu}},$
(64)

which implies

$\dfrac{p}{r}=\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\frac{(2-\delta)}{\delta}\bigg{]}.$
(65)

Therefore, the radial coordinate ($r$) corrected up to the first order Post-
Newtonian term is

$r=p\bigg{[}1+\dfrac{1}{c^{2}}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\frac{(2-\delta)}{\delta}\bigg{]}+O(c^{-4}).$
(66)

Next, we derive the expression for the angular velocity $\dot{\phi}$.
Substituting Eq. (66) into Eq. (56) and simplifying for $\dot{\phi}$, we
obtain

$\dot{\phi}=\dfrac{1}{p}\sqrt{\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{2GM_{Q}\epsilon^{2\delta-1}}{\delta
p^{2\delta}}(2-\delta)\bigg{]}.$ (67)

In the next section, we discuss the average energy radiation rate of EMRI
using the expressions of PN corrected quantities derived here.

##  IV Average Energy radiation rate from an EMRI

The average energy radiation rate $\langle\frac{dE}{dt}\rangle$ from a binary
is given by the following expression [37, 36, 31].

$\bigg{\langle}\frac{dE}{dt}\bigg{\rangle}=-\frac{G}{45c^{5}}\frac{1}{T}\int^{T}_{0}(\dddot{D}^{\alpha\beta})^{2}dt,$
(68)

where the mass quadrupole tensor $D^{\alpha\beta}$ is defined as [36]

$D^{\alpha\beta}=M_{P}(3x^{\alpha}x^{\beta}-|\mathbf{x}|^{2}\delta^{\alpha\beta}).$
(69)

Therefore,

$\displaystyle\dddot{D}^{\alpha\beta}$ $\displaystyle=$
$\displaystyle\dfrac{d^{3}}{dt^{3}}\big{[}M_{P}(3x^{\alpha}x^{\beta}-|\mathbf{x}|^{2}\delta^{\alpha\beta})\big{]}$
(70) $\displaystyle=$ $\displaystyle
M_{P}\dfrac{d^{2}}{dt^{2}}\big{[}3\dot{x}^{\alpha}x^{\beta}+3x^{\alpha}\dot{x}^{\beta}-2\mathbf{x}.\dot{\mathbf{x}}\delta^{\alpha\beta}\big{]}$
$\displaystyle=$ $\displaystyle
M_{P}\dfrac{d}{dt}\big{[}3\ddot{x}^{\alpha}x^{\beta}+6\dot{x}^{\alpha}\dot{x}^{\beta}+3x^{\alpha}\ddot{x}^{\beta}-2\dot{x}^{2}\delta^{\alpha\beta}-2\mathbf{x}.\ddot{\mathbf{x}}\delta^{\alpha\beta}\big{]}$
$\displaystyle=$ $\displaystyle
M_{P}\big{[}3\dddot{x}^{\alpha}x^{\beta}+9\ddot{x}^{\alpha}\dot{x}^{\beta}+9\dot{x}^{\alpha}\ddot{x}^{\beta}+3x^{\alpha}\dddot{x}^{\beta}-6\dot{\mathbf{x}}.\ddot{\mathbf{x}}\delta^{\alpha\beta}-2\mathbf{x}.\dddot{\mathbf{x}}\delta^{\alpha\beta}\big{]}$

Hence, $(\dddot{D}^{\alpha\beta})^{2}$ can be written as

$\displaystyle(\dddot{D}^{\alpha\beta})^{2}$ $\displaystyle=$
$\displaystyle(\dddot{D}^{\alpha\beta})(\dddot{D}_{\alpha\beta})$ (71)
$\displaystyle=$ $\displaystyle
6M_{P}^{2}\big{[}3\dddot{x}^{2}x^{2}+18(\dddot{\mathbf{x}}.\ddot{\mathbf{x}})(\dot{\mathbf{x}}.\mathbf{x})+18(\dddot{\mathbf{x}}.\dot{\mathbf{x}})(\ddot{\mathbf{x}}.\mathbf{x})+(\mathbf{x}.\dddot{\mathbf{x}})^{2}+27\dot{\mathbf{x}}^{2}\ddot{\mathbf{x}}^{2}+9(\dot{\mathbf{x}}.\ddot{\mathbf{x}})^{2}-12(\dddot{\mathbf{x}}.\mathbf{x})(\ddot{\mathbf{x}}.\dot{\mathbf{x}})\big{]},$

where the acceleration terms ($\ddot{\mathbf{x}}$) can be written as the
gradient of the potential or the magnitude of the force corresponding to this
acceleration is given as $\mathbf{F_{a}}=M_{P}|\mathbf{a}|=\dfrac{\partial
V}{\partial r}$, where $V$ is the potential energy. Therefore, the
acceleration of the body-P becomes

$\mathbf{a}=\frac{1}{M_{P}}\dfrac{\partial V}{\partial
r}\mathbf{\hat{n}}+\dfrac{1}{c^{2}}\mathbf{a}_{PN}.$ (72)

Substituting the above expression of $\mathbf{a}$ into Eq. (71) and after some
simplification, we obtain the following expression for
$(\dddot{D}^{\alpha\beta})^{2}$ ([38, 36]):

$\displaystyle(\dddot{D}^{\alpha\beta})^{2}=24\Bigg{[}\bigg{(}r\dot{r}\dfrac{\partial^{2}V}{\partial
r^{2}}+3\dot{r}\dfrac{\partial V}{\partial r}\bigg{)}^{2}$ $\displaystyle+$
$\displaystyle 12r^{2}\dot{\phi}^{2}\bigg{(}\dfrac{\partial V}{\partial
r}\bigg{)}^{2}\Bigg{]}+\dfrac{24M_{P}}{c^{2}}\Bigg{[}2r\bigg{(}\dot{r}\dfrac{\partial^{2}V}{\partial
r^{2}}\mathbf{x}+3\dfrac{\partial V}{\partial
r}\dot{\mathbf{x}}\bigg{)}.\dot{\mathbf{a}}_{PN}$ (73) $\displaystyle+$
$\displaystyle\bigg{\\{}9\dfrac{\partial V}{\partial
r}\big{(}\dot{r}^{2}+2r^{2}\dot{\phi}^{2}\big{)}\mathbf{\hat{n}}+3\dot{r}\dfrac{\partial^{2}V}{\partial
r^{2}}\big{(}3\dot{r}\mathbf{x}-r\dot{\mathbf{x}}\big{)}+9\dot{r}\dfrac{\partial
V}{\partial r}\dot{\mathbf{x}}\bigg{\\}}.\mathbf{a}_{PN}\Bigg{]},$

For circular orbits $\dot{r}=0$, therefore Eq. (73) reduces to the following
expression.

$\displaystyle(\dddot{D}^{\alpha\beta})^{2}$ $\displaystyle=$ $\displaystyle
24\Bigg{[}12\bigg{(}\dfrac{\partial V}{\partial
r}\bigg{)}^{2}r^{2}\dot{\phi}^{2}$ $\displaystyle+$
$\displaystyle\dfrac{M_{P}}{c^{2}}\bigg{\\{}18^{2}\dot{\phi}^{2}\dfrac{\partial
V}{\partial r}\mathbf{\hat{n}}.\mathbf{a}_{PN}+6r\dfrac{\partial V}{\partial
r}\mathbf{v}_{P}.\dot{\mathbf{a}}_{PN}\bigg{\\}}\Bigg{]}.$

Here, the time derivative of the PN correction term $\dot{\mathbf{a}}_{PN}$ of
the acceleration is

$\dot{\mathbf{a}}_{PN}=\dfrac{d\mathbf{a}_{PN}}{dt}=\dfrac{d\mathbf{a}_{PN}}{d\phi}\dfrac{d\phi}{dt}=\dfrac{d\mathbf{a}_{PN}}{d\phi}\dot{\phi}.$
(75)

Therefore,

$\dot{\mathbf{a}}_{PN}=\bigg{[}\dfrac{GM_{Q}\epsilon^{2\delta-1}}{p^{2\delta}}\bigg{]}^{5/2}\dfrac{(2-\delta)}{\delta
p^{2}}\bm{\hat{\nu}}+O(c^{-2}),$ (76)

where the higher order terms are ignored as they will contribute to the 2PN or
higher orders, which is out of the scope of this current study. We can also
write,

$\displaystyle V(r)$ $\displaystyle=$
$\displaystyle\dfrac{GM_{Q}M_{P}}{2\delta r^{2\delta}}\epsilon^{2\delta-1},$
(77) $\displaystyle\mathbf{v}_{P}.\dot{\mathbf{a}}_{PN}$ $\displaystyle=$
$\displaystyle\dfrac{G^{3}M_{Q}^{3}\epsilon^{3(2\delta-1)}}{\delta
p^{2(3\delta+1)}}(2-\delta)+O(c^{-2}),$ (78)
$\displaystyle\mathbf{\hat{n}}.\mathbf{a}_{PN}$ $\displaystyle=$
$\displaystyle\dfrac{G^{2}M_{Q}^{2}\epsilon^{2(2\delta-1)}}{\delta
p^{4\delta+1}}(2-\delta)+O(c^{-2}).$ (79)

Now substituting Eq. (66), (67), and the above expressions into Eq.
(LABEL:D^2_for_circular_orbits), we obtain

$(\dddot{D}^{\alpha\beta})^{2}=\dfrac{288G^{3}M_{Q}^{3}M_{P}^{2}\epsilon^{3(2\delta-1)}}{p^{2(3\delta+1)}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{2GM_{Q}\epsilon^{2\delta-1}}{\delta
p^{2\delta}}(2-\delta)(3+2\delta)\bigg{]}.$ (80)

The average rate of the energy radiation mentioned in Eq. (68) can be written
as,

$\displaystyle\bigg{\langle}\frac{dE}{dt}\bigg{\rangle}=-\frac{G}{45c^{5}}\frac{1}{2\pi}\int^{2\pi}_{0}(\dddot{D}^{\alpha\beta})^{2}d\phi\,\,,$
(81)

where $\dot{\phi}=2\pi/T$.

(a) For $\delta=0.27\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\
\delta<\delta_{0})$

(b) For $\delta=0.01\text{ i.e. }\approx 0$

(c) For $\delta=0.535\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\
\delta>0.5)$

Figure 1: Figure shows velocity of the smaller object in the extreme mass
ratio binary in power law (PL) potential (Orange curve) and in Kepler-Newton
(KN) potential (Blue curve) with respect to radial distance between the two
bodies.

In the case of circular orbits, the integrand is independent of the variable
$\phi$. Therefore, the above expression of the average energy radiation rate
reduces to

$\bigg{\langle}\frac{dE}{dt}\bigg{\rangle}=-\dfrac{G}{45c^{5}}(\dddot{D}^{\alpha\beta})^{2}.$
(82)

Thus

$\bigg{\langle}\frac{dE}{dt}\bigg{\rangle}=-\dfrac{32G^{4}M_{Q}^{3}M_{P}^{2}\epsilon^{3(2\delta-1)}}{5c^{5}p^{2(3\delta+1)}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{2GM_{Q}\epsilon^{2\delta-1}}{\delta
p^{2\delta}}(2-\delta)(3+2\delta)\bigg{]}.$ (83)

The condition for using PN approximations to evaluate average energy radiation
rate for power law potential is

$1\gg\frac{GM_{Q}}{c^{2}p^{2\delta}}\frac{(2-\delta)(3+2\delta)}{\delta}.$
(84)

For a given EMRI and in a range of radial distance, this condition can be
narrowed down to a cut-off $\delta_{0}$ below which PN approximations cannot
be used to work out the average energy radiation rate. This cut-off depends on
the properties of the binary under consideration.

(a) $\delta=0.27\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\
\delta<\delta_{0})$

(b) $\delta=0.395\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\
\delta_{0}<\delta<0.5)$

(c) $\delta=0.535\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\ \delta>0.5)$

Figure 2: Figure shows variation of 1-PN corrected rate of energy radiation
from an extreme mass ratio binary in power law (PL) potential (Orange curve)
and in Kepler-Newton (KN) potential (Blue curve) with respect to radial
distance ($p$) between the two bodies.

(a) $\delta=0.395\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\
\delta_{0}<\delta<0.5)$, $\Delta=0.105$

(b) $\delta=0.535\leavevmode\nobreak\ (i.e.\leavevmode\nobreak\
\delta>0.5),\leavevmode\nobreak\ \Delta=-0.035$

Figure 3: Figure shows the comparison of the orbital frequency of the extreme
mass ratio binary in power law (PL) potential (Orange curve) and in Kepler-
Newton (KN) potential (Blue curve) with respect to radial distance ($p$)
between the two bodies.

Eq. (83) shows the energy radiation from the binaries following power law
potential with the assumptions that (i) the scale length ($\epsilon$) is very
small compared to the distance between two bodies i.e. $r_{PQ}$, (ii) the
body-Q is much heavier than the body-P, therefore the motion of the body-Q is
negligible, and it can be treated as a stationary body and (iii) the body-P is
moving in the circular orbit around Q.

##  V  Comparison of power Law potential with the Kepler-Newton Potential

Consider an example of an EMRI system of a supermassive black hole, and a star
as massive as the sun.
Let the properties of the supermassive compact object be: mass, $M_{Q}$ =
$4.1\times 10^{6}M_{\odot}=8.2\times 10^{36}kg$, the scale radius, $\epsilon$=
$5.8\times 10^{-4}pc$ and the properties of the star or the stellar mass
object be: mass, $M_{P}$ = $1M_{\odot}$ = $2\times 10^{30}kg$.

In this section, we compare the velocity curves, the average energy radiation
rate and the orbital frequency of the binary in the PL potential and in the KN
Potential. The PN corrections are much smaller than the Newtonian values of
the evaluated quantities. Hence, in the magnitude scales considered below, the
graphs of the Newtonian and the PN corrected quantities will apparently
overlap.

###  V.1 Velocity profile

The comparison of the velocity curves for power law potential (Eq. (62)) and
Kepler-Newton potential is shown in this subsection. For $0<\delta<0.5$, the
velocity falls more slowly in the case of power law (PL) potential as compared
to that in KN potential case which is usually observed in astrophysical
scenarios (Fig. 1(a)). In fact, for $\delta=0.01\text{ i.e. }\approx 0$, the
velocity curve is nearly flat (Fig. 1(b)). For the case of $0.5<\delta$, the
velocity for PL potential drops faster with increasing radial distance (Fig.
1(c)). Therefore, $0.5<\delta$ is astrophysically improbable.

###  V.2 Average energy radiation rate

Next, we compare the average energy radiation rate from an EMRI in PL
potential and that in KN Potential.

Average energy radiation rate from a system in KN potential can be determined
by substituting $\delta=1/2$ in our final result (Eq. (83))

$\bigg{\langle}\frac{dE}{dt}\bigg{\rangle}_{KNP}=-\frac{32G^{4}M_{Q}^{3}M_{P}^{2}}{5c^{5}p^{5}}\bigg{(}1-\dfrac{1}{c^{2}}\dfrac{24GM_{Q}}{p}\bigg{)}.$
(85)

The comparisons of the two energy radiation rates, Eq. (85), and Eq. (83) with
$\delta\neq 0.5$, are shown in the Fig. (2) for different $\delta$-s.

In the EMRI considered above, we focus on the region where the radial distance
$p\approx 10^{16}$ meter and the corresponding $\delta_{0}\approx 0.38$. For
$0<\delta<0.38$, the 1-PN term becomes comparable and much higher than the
0-PN (or the Newtonian) part. Therefore, the 1-PN part blows up and it looks
as if the EMRI is absorbing energy instead of radiating (Fig. 2(a)). Hence, as
far as the energy radiation rate is concerned, we have to limit our study to
$\delta>\delta_{0}:=0.38$. For $0.38<\delta<0.5$, the energy radiation rate is
greater than that from KN potential (Fig. 2(b)). For $\delta>0.5$, the energy
radiation rate is smaller than that from KN potential (Fig. 2(c)). However, as
we mentioned before, $\delta>0.5$ is astrophysically improbable.

###  V.3 Orbital frequency

Presence of matter distribution around a binary affects the dynamics of the
binary. The effect of dynamical friction as studied in [39, 40, 41] is ignored
here. Let’s consider a power law (where the power is $\delta$) such that it
has a small deviation ($\Delta$) in the power from the Kepler-Newton potential
($\delta_{KN}=1/2$) which could be because of the presence of matter
distribution around the EMRI. Accordingly, from the Eq. (67), the orbital
frequency $\dot{\phi}$ can be written as

$\dot{\phi}=\dfrac{1}{p}\sqrt{\dfrac{GM_{Q}\epsilon^{1-2\Delta}}{p^{(1-2\Delta)}}}\bigg{[}1-\dfrac{1}{c^{2}}\dfrac{2GM_{Q}\epsilon^{-2\Delta}}{(1/2-\Delta)p^{(1-2\Delta)}}(3/2+\Delta)\bigg{]}.$
(86)

When $\Delta>0$, which is also the case where we expect non-vacuum
environment, the orbital frequency is increased (Fig. 3(a)). While in the case
where $\Delta<0$ we get orbital frequency lower than that of an EMRI in KN
potential but this case astrophysically improbable.

From the observation data of the orbital frequency of EMRI, the value of
$\Delta$ can be obtained which gives the precise power law potential which can
be used to get the matter density profile around the central supermassive
compact object.

##  VI Conclusion

In this paper, we consider an extreme mass ratio inspiral (EMRI) in a system
where the supermassive object lies in a non vacuum region such that we can use
power law potential. The main motivation to consider power law potential is to
get nearly flat profile of orbital velocity which is observed in astrophysical
scenarios (e.g. galactic rotation curve). We derive the 1-PN corrected
dynamical variables of the orbiting stellar mass object. We consider circular
orbits for the ease of calculations. Using the Newtonian gravitational
potential ($U$) and the PN potentials ($U^{j},$ $\Psi$ and $X$), we derive the
acceleration of a body in the PL potential. Using this, we derive other PN
corrected dynamical variables. Subsequently, we obtain the mass quadrupole
tensor which is ultimately used to obtain the average energy radiation rate.

Next, we discuss an example of an EMRI. We calculate the average energy
radiation rate from this EMRI in power law potential at different $\delta$-s
and compared them with that from the same EMRI at $\delta=0.5$ i.e. the KN
potential case. When $\delta$ is less than a particular cut-off $\delta_{0}$,
we cannot use PN approximations because it violates the conditions that the
value of 1-PN part should be much less than 0-PN part. The important results
we get here are,

  * •

The average energy radiation rate that we have derived is applicable for every
EMRI which is surrounded by a matter distribution that can give power law
potential in the region around the supermassive compact object. We show that
the average energy radiation rate is higher in general PL potential case as
compared to that in KN potential.

  * •

We use the comparison of orbital frequency of an EMRI in PL and KN potentials
to show how matter distribution (e.g. dark matter distribution) causes change
in the signals. This can be used to further investigate the quantities like
frequencies of the gravitational waves radiated from such an EMRI. The
deviation of the observed orbital frequencies from the frequency of a similar
binary in KN potential can be used to evaluate the value of $\delta$ and that
can be used to find out the effective potential and subsequently, can
calculate the mass density profile around the binary.

##  VII Acknowledgement

C.N.G. would like to acknowledge the support of the International Center for
Cosmology, CHARUSAT, India for funding the work.

## References

  * [1] Sofue. Y and Rubin. V, Annual Review of Astronomy and Astrophysics 39, no. 1, 137 (2001).
  * [2] Milgrom, M., Astrophysical Journal, Vol. 270, p. 365-370 (1983).
  * [3] Evstigneeva, E. A., Yu.A. Nikolaev publisher,p. 180-196, (2004).
  * [4] Rodrigues, D. C., Letelier, P. S., & Shapiro, I. L., Journal of Cosmology and Astroparticle Physics 2010, no. 04, 020 (2010).
  * [5] Boris E. Meierovich, Phys. Rev. D 87, 103510 (2013).
  * [6] Zwicky, F., General Relativity and Gravitation, 41, no. 1, pp.207-224 (2009).
  * [7] Rubin, V. C., Ford, W. K., & Thonnard, N., Astrophysical Journal, 238, p. 471-487 (1980).
  * [8] Rubin, V. C. & Ford, W. K., Astrophysical Journal, 159, p.379 (1970).
  * [9] Díaz E., Zandivarez A., Merchán M. E., et al. Astrophys. J. , 629, 158, (2005).
  * [10] Donato, F., Gentile, G., Salucci, P. Mon. Not. Roy. Astron. Soc., 353, L17, (2004).
  * [11] Navarro, J. F., Frenk, C. S., & White, S. D. M., Astrophysical Journal 462, p. 563 (1996).
  * [12] Banik, U., Dey, D., Bhattacharya, K., et al. , General Relativity and Gravitation, 49, 116, (2017).
  * [13] Dey, D., Bhattacharya, K., & Sarkar, T., General Relativity and Gravitation, 47, 103 (2015).
  * [14] Dey, D., Bhattacharya, K., & Sarkar, T., Physical Review D, 87, no. 10, 103505 (2013).
  * [15] Dey, D., Bhattacharya, K., & Sarkar, T., Physical Review D, 88, no. 8, 083532 (2013).
  * [16] The, L. S. & White, S. D. M., Astronomical Journal 95, p. 1642 (1988).
  * [17] Bekenstein, J. D., Physical Review D, 70, no. 8, 083509 (2004).
  * [18] J W Moffat, Journal of Cosmology and Astroparticle Physics, 2005, no. 5, (2005).
  * [19] Moffat, J. W., Journal of Cosmology and Astroparticle Physics, 2006, no. 03, 004 (2006).
  * [20] Brownstein, J. R. and Moffat, J. W., The Astrophysical Journal, 636, no. 2, (2006)
  * [21] Munera, H. A. & Delgado-Correal, C., Research Notes of the AAS, 4, no. 9, 151 (2020).
  * [22] Goldstein, H, Addison-Wesley, (1950).
  * [23] Danby, J. M. A., Willmann-Bell, 2nd ed., rev. & enl (1988). .
  * [24] Valluri, S. R. and Wiegert, P. A. and Drozd, J. and Da Silva, M.,  Mon. Not. Roy. Astron. Soc. 427, no. 3, 2392 (2012).
  * [25] Lynden-Bell, D. and Jin, S.,  Mon. Not. Roy. Astron. Soc. 386, no. 1, 245 (2008).
  * [26] Valluri, S. R. and Yu, P. and Smith, G. E. and Wiegert, P. A.,  Mon. Not. Roy. Astron. Soc. 358, no. 4, 1273 (2005).
  * [27] Barack Leor, Cutler Curt, Physical Review D, 69, no. 8, 1550-2368 (2004).
  * [28] Amaro-Seoane, P., Aoudia, S., Babak, S., et al. , Classical and Quantum Gravity, 29, 124016, (2012).
  * [29] Glampedakis Kostas, Classical and Quantum Gravity, 22, no. 15, S605–S659 (2005).
  * [30] Gair, J. R., Vallisneri, M., Larson, S. L. et al., Living Reviews in Relativity, 16, no. 1, 1433-8351 (2013).
  * [31] Poisson. E. and Will, C. M., Cambridge press (2014).
  * [32] Struck, Curtis,  Mon. Not. Roy. Astron. Soc. 446, no. 3, 3139 (2015).
  * [33] Struck, C, The Astron. Journ. 131, no. 3, 1347 (2006).
  * [34] J. Binney and S. Tremaine, Princeton, N.J. : Princeton University Press, c (1987).
  * [35] P. Bambhaniya, A. B. Joshi, D. Dey and P. S. Joshi, Phys. Rev. D 100, no.12, 124020 (2019).
  * [36] Dionysiou D. D., Astrophysics and Space Science 125, no. 1, pp.115-125 (1986).
  * [37] Wagoner, R. V. and Will, C. M., Astrophysical Journal, vol. 210, Dec. 15, pt. 1, p. 764-775, (1976).
  * [38] Dehnen, H. and Ghaboussi, F., Nuovo Cimento B, Serie 11 (ISSN 0369-3554), vol. 89 B, Oct. 11, p. 131-138, (1985).
  * [39] Bekenstein, J. D. & Zamir, R. Astrophys. J. , 359, 427, (1990).
  * [40] Gómez, L. G. & Rueda, J. A. Phys. Rev. D, 96, 063001, (2017).
  * [41] A. Caputo, J. Zavala and D. Blas, Phys. Dark Univ. 19, 1-11 (2018) [arXiv:1709.03991 [astro-ph.HE]].