Contact
CoCalc Logo Icon
StoreFeaturesDocsShareSupport News AboutSign UpSign In
| Download

5-dimensional Kerr-AdS spacetime with a Nambu-Goto string

Project: KerrAdS
Views: 147
License: GPL3
Image: ubuntu2004
Kernel: SageMath 9.3

5D Kerr-AdS spacetime with a Nambu-Goto string

Case b = n a

This SageMath notebook is relative to the article Heavy quarks in rotating plasma via holography by Anastasia A. Golubtsova, Eric Gourgoulhon and Marina K. Usova, arXiv:2107.11672.

The involved differential geometry computations are based on tools developed through the SageManifolds project.

NB: a version of SageMath at least equal to 9.1 is required to run this notebook:

version()
'SageMath version 9.3, Release Date: 2021-05-09'

First we set up the notebook to display mathematical objects using LaTeX rendering:

%display latex

Since some computations are quite long, we ask for running them in parallel on 8 cores:

Parallelism().set(nproc=8)

Spacetime manifold

We declare the Kerr-AdS spacetime as a 5-dimensional Lorentzian manifold:

M = Manifold(5, 'M', r'\mathcal{M}', structure='Lorentzian', metric_name='G') print(M)
5-dimensional Lorentzian manifold M

Let us define Boyer-Lindquist-type coordinates (rational polynomial version) on M\mathcal{M}, via the method chart(), the argument of which is a string expressing the coordinates names, their ranges (the default is (,+)(-\infty,+\infty)) and their LaTeX symbols:

BL.<t,r,mu,ph,ps> = M.chart(r't r:(0,+oo) mu:(0,1):\mu ph:(0,2*pi):\phi ps:(0,2*pi):\psi') BL
(M,(t,r,μ,ϕ,ψ))\renewcommand{\Bold}[1]{\mathbf{#1}}\left(\mathcal{M},(t, r, {\mu}, {\phi}, {\psi})\right)

The coordinate μ\mu is related to the standard Boyer-Lindquist coordinate θ\theta by μ=cosθ \mu = \cos\theta

The coordinate ranges are

BL.coord_range()
t: (,+);r: (0,+);μ: (0,1);ϕ: (0,2π);ψ: (0,2π)\renewcommand{\Bold}[1]{\mathbf{#1}}t :\ \left( -\infty, +\infty \right) ;\quad r :\ \left( 0 , +\infty \right) ;\quad {\mu} :\ \left( 0 , 1 \right) ;\quad {\phi} :\ \left( 0 , 2 \, \pi \right) ;\quad {\psi} :\ \left( 0 , 2 \, \pi \right)

Note that contrary to the 4-dimensional case, the range of μ\mu is (0,1)(0,1) only (cf. Sec. 1.2 of R.C. Myers, arXiv:1111.1903 or Sec. 2 of G.W. Gibbons, H. Lüb, Don N. Page, C.N. Pope, J. Geom. Phys. 53, 49 (2005)). In other words, the range of θ\theta is (0,π2)\left(0, \frac{\pi}{2}\right) only.

Metric tensor

The 4 parameters mm, aa, bb and \ell of the Kerr-AdS spacetime are declared as symbolic variables, aa and bb being the two angular momentum parameters and \ell being related to the cosmological constant by Λ=62\Lambda = - 6 \ell^2:

var('m a b', domain='real')
(m,a,b)\renewcommand{\Bold}[1]{\mathbf{#1}}\left(m, a, b\right)
var('l', domain='real', latex_name=r'\ell')
\renewcommand{\Bold}[1]{\mathbf{#1}}{\ell}

We assume that b=nab = n a:

n = var('n', domain='real') b = n*a

Some auxiliary functions:

keep_Delta = True # change to False to provide explicit expression for Delta_r, Xi_a, etc...
sig = (1 + r^2*l^2)/r^2 costh2 = mu^2 sinth2 = 1 - mu^2 rho2 = r^2 + a^2*mu^2 + b^2*sinth2 # Explicit expressions: Delta_r_expr = (r^2+a^2)*(r^2+b^2)*sig - 2*m Delta_th_expr = 1 - a^2*l^2*costh2 - b^2*l^2*sinth2 Xi_a_expr = 1 - a^2*l^2 Xi_b_expr = 1 - b^2*l^2 if keep_Delta: Delta_r = var('Delta_r', latex_name=r'\Delta_r', domain='real') Delta_th = var('Delta_th', latex_name=r'\Delta_\theta', domain='real') if a == b: Xi_a = var('Xi', latex_name=r'\Xi', domain='real') Xi_b = Xi_a else: Xi_a = var('Xi_a', latex_name=r'\Xi_a', domain='real') Xi_b = var('Xi_b', latex_name=r'\Xi_b', domain='real') else: Delta_r = Delta_r_expr Delta_th = Delta_th_expr Xi_a = Xi_a_expr Xi_b = Xi_b_expr

The metric is set by its components in the coordinate frame associated with the Boyer-Lindquist-type coordinates, which is the current manifold's default frame. These components can be deduced from Eq. (5.22) of the article S.W. Hawking, C.J. Hunter & M.M. Taylor-Robinson, Phys. Rev. D 59, 064005 (1999) (the check of agreement with this equation is performed below):

G = M.metric() tmp = 1/rho2*( -Delta_r + Delta_th*(a^2*sinth2 + b^2*mu^2) + a^2*b^2*sig ) G[0,0] = tmp.simplify_full() tmp = a*sinth2/(rho2*Xi_a)*( Delta_r - (r^2+a^2)*(Delta_th + b^2*sig) ) G[0,3] = tmp.simplify_full() tmp = b*mu^2/(rho2*Xi_b)*( Delta_r - (r^2+b^2)*(Delta_th + a^2*sig) ) G[0,4] = tmp.simplify_full() G[1,1] = (rho2/Delta_r).simplify_full() G[2,2] = (rho2/Delta_th/(1-mu^2)).simplify_full() tmp = sinth2/(rho2*Xi_a^2)*( - Delta_r*a^2*sinth2 + (r^2+a^2)^2*(Delta_th + sig*b^2*sinth2) ) G[3,3] = tmp.simplify_full() tmp = a*b*sinth2*mu^2/(rho2*Xi_a*Xi_b)*( - Delta_r + sig*(r^2+a^2)*(r^2+b^2) ) G[3,4] = tmp.simplify_full() tmp = mu^2/(rho2*Xi_b^2)*( - Delta_r*b^2*mu^2 + (r^2+b^2)^2*(Delta_th + sig*a^2*mu^2) ) G[4,4] = tmp.simplify_full()
G.display_comp(only_nonredundant=True)
Gtttt=a4n2(Δθa2μ2Δθa2(a42+Δθa2μ2)n2+Δr)r2r4+(a2μ2(a2μ2a2)n2)r2Gtϕtϕ=(Δθaμ2+(a32μ2a32)n2Δθa)r4+(a5μ2a5)n2(Δθa3(Δθa3Δra)μ2+(a52+a3(a52+a3)μ2)n2Δra)r2Ξar4+(Ξaa2μ2(Ξaa2μ2Ξaa2)n2)r2Gtψtψ=a5μ2n3+(a32+Δθa)μ2nr4+((a52+Δθa3)μ2n3+(a3Δra)μ2n)r2Ξbr4+(Ξba2μ2(Ξba2μ2Ξba2)n2)r2Grrrr=a2μ2(a2μ2a2)n2+r2ΔrGμμμμ=a2μ2(a2μ2a2)n2+r2Δθμ2ΔθGϕϕϕϕ=(Δθμ2(a22μ42a22μ2+a22)n2Δθ)r6+(2Δθa2μ22Δθa2(2a42+(2a42+a2)μ42(2a42+a2)μ2+a2)n2)r4(a6μ42a6μ2+a6)n2+(Δra2μ4Δθa4+Δra2+(Δθa42Δra2)μ2(a62+(a62+2a4)μ4+2a42(a62+2a4)μ2)n2)r2Ξa2r4+(Ξa2a2μ2(Ξa2a2μ2Ξa2a2)n2)r2Gϕψϕψ=(a22μ4a22μ2)nr6+((a42μ4a42μ2)n3+((a42+a2)μ4(a42+a2)μ2)n)r4+(a6μ4a6μ2)n3+(((a62+a4)μ4(a62+a4)μ2)n3+((a4Δra2)μ4(a4Δra2)μ2)n)r2ΞaΞbr4+(ΞaΞba2μ2(ΞaΞba2μ2ΞaΞba2)n2)r2Gψψψψ=a6μ4n4+(a22μ4+Δθμ2)r6+(a2μ4+2(a42μ4+Δθa2μ2)n2)r4+((2a4Δra2)μ4n2+(a62μ4+Δθa4μ2)n4)r2Ξb2r4+(Ξb2a2μ2(Ξb2a2μ2Ξb2a2)n2)r2\renewcommand{\Bold}[1]{\mathbf{#1}}\begin{array}{lcl} G_{ \, t \, t }^{ \phantom{\, t}\phantom{\, t} } & = & \frac{a^{4} n^{2} - {\left({\Delta_\theta} a^{2} {\mu}^{2} - {\Delta_\theta} a^{2} - {\left(a^{4} {\ell}^{2} + {\Delta_\theta} a^{2} {\mu}^{2}\right)} n^{2} + {\Delta_r}\right)} r^{2}}{r^{4} + {\left(a^{2} {\mu}^{2} - {\left(a^{2} {\mu}^{2} - a^{2}\right)} n^{2}\right)} r^{2}} \\ G_{ \, t \, {\phi} }^{ \phantom{\, t}\phantom{\, {\phi}} } & = & \frac{{\left({\Delta_\theta} a {\mu}^{2} + {\left(a^{3} {\ell}^{2} {\mu}^{2} - a^{3} {\ell}^{2}\right)} n^{2} - {\Delta_\theta} a\right)} r^{4} + {\left(a^{5} {\mu}^{2} - a^{5}\right)} n^{2} - {\left({\Delta_\theta} a^{3} - {\left({\Delta_\theta} a^{3} - {\Delta_r} a\right)} {\mu}^{2} + {\left(a^{5} {\ell}^{2} + a^{3} - {\left(a^{5} {\ell}^{2} + a^{3}\right)} {\mu}^{2}\right)} n^{2} - {\Delta_r} a\right)} r^{2}}{{\Xi_a} r^{4} + {\left({\Xi_a} a^{2} {\mu}^{2} - {\left({\Xi_a} a^{2} {\mu}^{2} - {\Xi_a} a^{2}\right)} n^{2}\right)} r^{2}} \\ G_{ \, t \, {\psi} }^{ \phantom{\, t}\phantom{\, {\psi}} } & = & -\frac{a^{5} {\mu}^{2} n^{3} + {\left(a^{3} {\ell}^{2} + {\Delta_\theta} a\right)} {\mu}^{2} n r^{4} + {\left({\left(a^{5} {\ell}^{2} + {\Delta_\theta} a^{3}\right)} {\mu}^{2} n^{3} + {\left(a^{3} - {\Delta_r} a\right)} {\mu}^{2} n\right)} r^{2}}{{\Xi_b} r^{4} + {\left({\Xi_b} a^{2} {\mu}^{2} - {\left({\Xi_b} a^{2} {\mu}^{2} - {\Xi_b} a^{2}\right)} n^{2}\right)} r^{2}} \\ G_{ \, r \, r }^{ \phantom{\, r}\phantom{\, r} } & = & \frac{a^{2} {\mu}^{2} - {\left(a^{2} {\mu}^{2} - a^{2}\right)} n^{2} + r^{2}}{{\Delta_r}} \\ G_{ \, {\mu} \, {\mu} }^{ \phantom{\, {\mu}}\phantom{\, {\mu}} } & = & -\frac{a^{2} {\mu}^{2} - {\left(a^{2} {\mu}^{2} - a^{2}\right)} n^{2} + r^{2}}{{\Delta_\theta} {\mu}^{2} - {\Delta_\theta}} \\ G_{ \, {\phi} \, {\phi} }^{ \phantom{\, {\phi}}\phantom{\, {\phi}} } & = & -\frac{{\left({\Delta_\theta} {\mu}^{2} - {\left(a^{2} {\ell}^{2} {\mu}^{4} - 2 \, a^{2} {\ell}^{2} {\mu}^{2} + a^{2} {\ell}^{2}\right)} n^{2} - {\Delta_\theta}\right)} r^{6} + {\left(2 \, {\Delta_\theta} a^{2} {\mu}^{2} - 2 \, {\Delta_\theta} a^{2} - {\left(2 \, a^{4} {\ell}^{2} + {\left(2 \, a^{4} {\ell}^{2} + a^{2}\right)} {\mu}^{4} - 2 \, {\left(2 \, a^{4} {\ell}^{2} + a^{2}\right)} {\mu}^{2} + a^{2}\right)} n^{2}\right)} r^{4} - {\left(a^{6} {\mu}^{4} - 2 \, a^{6} {\mu}^{2} + a^{6}\right)} n^{2} + {\left({\Delta_r} a^{2} {\mu}^{4} - {\Delta_\theta} a^{4} + {\Delta_r} a^{2} + {\left({\Delta_\theta} a^{4} - 2 \, {\Delta_r} a^{2}\right)} {\mu}^{2} - {\left(a^{6} {\ell}^{2} + {\left(a^{6} {\ell}^{2} + 2 \, a^{4}\right)} {\mu}^{4} + 2 \, a^{4} - 2 \, {\left(a^{6} {\ell}^{2} + 2 \, a^{4}\right)} {\mu}^{2}\right)} n^{2}\right)} r^{2}}{{\Xi_a}^{2} r^{4} + {\left({\Xi_a}^{2} a^{2} {\mu}^{2} - {\left({\Xi_a}^{2} a^{2} {\mu}^{2} - {\Xi_a}^{2} a^{2}\right)} n^{2}\right)} r^{2}} \\ G_{ \, {\phi} \, {\psi} }^{ \phantom{\, {\phi}}\phantom{\, {\psi}} } & = & -\frac{{\left(a^{2} {\ell}^{2} {\mu}^{4} - a^{2} {\ell}^{2} {\mu}^{2}\right)} n r^{6} + {\left({\left(a^{4} {\ell}^{2} {\mu}^{4} - a^{4} {\ell}^{2} {\mu}^{2}\right)} n^{3} + {\left({\left(a^{4} {\ell}^{2} + a^{2}\right)} {\mu}^{4} - {\left(a^{4} {\ell}^{2} + a^{2}\right)} {\mu}^{2}\right)} n\right)} r^{4} + {\left(a^{6} {\mu}^{4} - a^{6} {\mu}^{2}\right)} n^{3} + {\left({\left({\left(a^{6} {\ell}^{2} + a^{4}\right)} {\mu}^{4} - {\left(a^{6} {\ell}^{2} + a^{4}\right)} {\mu}^{2}\right)} n^{3} + {\left({\left(a^{4} - {\Delta_r} a^{2}\right)} {\mu}^{4} - {\left(a^{4} - {\Delta_r} a^{2}\right)} {\mu}^{2}\right)} n\right)} r^{2}}{{\Xi_a} {\Xi_b} r^{4} + {\left({\Xi_a} {\Xi_b} a^{2} {\mu}^{2} - {\left({\Xi_a} {\Xi_b} a^{2} {\mu}^{2} - {\Xi_a} {\Xi_b} a^{2}\right)} n^{2}\right)} r^{2}} \\ G_{ \, {\psi} \, {\psi} }^{ \phantom{\, {\psi}}\phantom{\, {\psi}} } & = & \frac{a^{6} {\mu}^{4} n^{4} + {\left(a^{2} {\ell}^{2} {\mu}^{4} + {\Delta_\theta} {\mu}^{2}\right)} r^{6} + {\left(a^{2} {\mu}^{4} + 2 \, {\left(a^{4} {\ell}^{2} {\mu}^{4} + {\Delta_\theta} a^{2} {\mu}^{2}\right)} n^{2}\right)} r^{4} + {\left({\left(2 \, a^{4} - {\Delta_r} a^{2}\right)} {\mu}^{4} n^{2} + {\left(a^{6} {\ell}^{2} {\mu}^{4} + {\Delta_\theta} a^{4} {\mu}^{2}\right)} n^{4}\right)} r^{2}}{{\Xi_b}^{2} r^{4} + {\left({\Xi_b}^{2} a^{2} {\mu}^{2} - {\left({\Xi_b}^{2} a^{2} {\mu}^{2} - {\Xi_b}^{2} a^{2}\right)} n^{2}\right)} r^{2}} \end{array}

Check of Eq. (2.9)

We need the 1-forms dt\mathrm{d}t, dr\mathrm{d}r, dμ\mathrm{d}\mu, dϕ\mathrm{d}\phi and dψ\mathrm{d}\psi:

dt, dr, dmu, dph, dps = (BL.coframe()[i] for i in M.irange()) dt, dr, dmu, dph, dps
(dt,dr,dμ,dϕ,dψ)\renewcommand{\Bold}[1]{\mathbf{#1}}\left(\mathrm{d} t, \mathrm{d} r, \mathrm{d} {\mu}, \mathrm{d} {\phi}, \mathrm{d} {\psi}\right)
print(dt)
1-form dt on the 5-dimensional Lorentzian manifold M

In agreement with μ=cosθ\mu = \cos\theta, we introduce the 1-form dθ=dμ/sinθ,\mathrm{d}\theta = - \mathrm{d}\mu /\sin\theta , with sinθ=1μ2\sin\theta = \sqrt{1-\mu^2} since θ(0,π2)\theta\in\left(0, \frac{\pi}{2}\right) :

dth = - 1/sqrt(1 - mu^2)*dmu
s1 = dt - a*sinth2/Xi_a*dph - b*costh2/Xi_b*dps s1.display()
dt+(aμ2aΞa)dϕaμ2nΞbdψ\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{d} t + \left( \frac{a {\mu}^{2} - a}{{\Xi_a}} \right) \mathrm{d} {\phi} -\frac{a {\mu}^{2} n}{{\Xi_b}} \mathrm{d} {\psi}
s2 = a*dt - (r^2 + a^2)/Xi_a*dph s2.display()
adt+(a2+r2Ξa)dϕ\renewcommand{\Bold}[1]{\mathbf{#1}}a \mathrm{d} t + \left( -\frac{a^{2} + r^{2}}{{\Xi_a}} \right) \mathrm{d} {\phi}
s3 = b*dt - (r^2 + b^2)/Xi_b*dps s3.display()
andt+(a2n2+r2Ξb)dψ\renewcommand{\Bold}[1]{\mathbf{#1}}a n \mathrm{d} t + \left( -\frac{a^{2} n^{2} + r^{2}}{{\Xi_b}} \right) \mathrm{d} {\psi}
s4 = a*b*dt - b*(r^2 + a^2)*sinth2/Xi_a * dph - a*(r^2 + b^2)*costh2/Xi_b * dps s4.display()
a2ndt+((aμ2a)nr2+(a3μ2a3)nΞa)dϕ+(a3μ2n2+aμ2r2Ξb)dψ\renewcommand{\Bold}[1]{\mathbf{#1}}a^{2} n \mathrm{d} t + \left( \frac{{\left(a {\mu}^{2} - a\right)} n r^{2} + {\left(a^{3} {\mu}^{2} - a^{3}\right)} n}{{\Xi_a}} \right) \mathrm{d} {\phi} + \left( -\frac{a^{3} {\mu}^{2} n^{2} + a {\mu}^{2} r^{2}}{{\Xi_b}} \right) \mathrm{d} {\psi}
G0 = - Delta_r/rho2 * s1*s1 + Delta_th*sinth2/rho2 * s2*s2 \ + Delta_th*costh2/rho2 * s3*s3 + rho2/Delta_r * dr*dr \ + rho2/Delta_th * dth*dth + sig/rho2 * s4*s4

Check of Eq. (2.9):

G0 == G
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}

Einstein equation

The Ricci tensor of gg is

if not keep_Delta: # Ric = G.ricci() # print(Ric) pass
if not keep_Delta: # show(Ric.display_comp(only_nonredundant=True)) pass

Let us check that gg is a solution of the vacuum Einstein equation with the cosmological constant Λ=62\Lambda = - 6 \ell^2:

Lambda = -6*l^2 if not keep_Delta: # print(Ric == 2/3*Lambda*G) pass

String worldsheet

The string worldsheet as a 2-dimensional Lorentzian submanifold of M\mathcal{M}:

W = Manifold(2, 'W', ambient=M, structure='Lorentzian', latex_name=r'\mathcal{W}') print(W)
2-dimensional Lorentzian submanifold W immersed in the 5-dimensional Lorentzian manifold M

Let us assume that the string worldsheet is parametrized by (t,r)(t,r):

XW.<t,r> = W.chart(r't r:(0,+oo)') XW
(W,(t,r))\renewcommand{\Bold}[1]{\mathbf{#1}}\left(\mathcal{W},(t, r)\right)

The string embedding in Kerr-AdS spacetime, as an expansion about a straight string solution in AdS (Eq. (4.5) of the paper):

Mu0 = var('Mu0', latex_name=r'\mu_0', domain='real') Phi0 = var('Phi0', latex_name=r'\Phi_0', domain='real') Psi0 = var('Psi0', latex_name=r'\Psi_0', domain='real') cosTh0 = Mu0 sinTh0 = sqrt(1 - Mu0^2) mu_s = Mu0 + (a+b)^2*function('mu_1')(r) ph_s = Phi0 + a*l^2*t + a*function('phi_1')(r) ps_s = Psi0 + b*l^2*t + b*function('psi_1')(r) F = W.diff_map(M, {(XW, BL): [t, r, mu_s, ph_s, ps_s]}, name='F') W.set_embedding(F) F.display()
F:WM(t,r)(t,r,μ,ϕ,ψ)=(t,r,(an+a)2μ1(r)+μ0,a2t+aϕ1(r)+Φ0,a2nt+anψ1(r)+Ψ0)\renewcommand{\Bold}[1]{\mathbf{#1}}\begin{array}{llcl} F:& \mathcal{W} & \longrightarrow & \mathcal{M} \\ & \left(t, r\right) & \longmapsto & \left(t, r, {\mu}, {\phi}, {\psi}\right) = \left(t, r, {\left(a n + a\right)}^{2} \mu_{1}\left(r\right) + {\mu_0}, a {\ell}^{2} t + a \phi_{1}\left(r\right) + {\Phi_0}, a {\ell}^{2} n t + a n \psi_{1}\left(r\right) + {\Psi_0}\right) \end{array}
F.jacobian_matrix()
(10010(a2n2+2a2n+a2)rμ1(r)a2arϕ1(r)a2nanrψ1(r))\renewcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{rr} 1 & 0 \\ 0 & 1 \\ 0 & {\left(a^{2} n^{2} + 2 \, a^{2} n + a^{2}\right)} \frac{\partial}{\partial r}\mu_{1}\left(r\right) \\ a {\ell}^{2} & a \frac{\partial}{\partial r}\phi_{1}\left(r\right) \\ a {\ell}^{2} n & a n \frac{\partial}{\partial r}\psi_{1}\left(r\right) \end{array}\right)

Induced metric on the string worldsheet

Because of the bug #27492, which impedes parallel computations involving symbolic functions, such as μ1\mu_1, we switch back to serial computations:

Parallelism().set(nproc=1)

The metric on the string worldsheet W\mathcal{W} is the metric gg induced by the spacetime metric GG, i.e. the pullback of GG by the embedding FF:

g = W.induced_metric()

Nambu-Goto action

The determinant of gg is

detg = g.determinant().expr()

Expanding at fourth order in aa (will be required latter):

detg_a4 = detg.series(a, 5).truncate().simplify_full()
detg_a40 = detg_a4
detg_a4 = detg_a40
detg_a4 = detg_a4.subs({Xi_a: Xi_a_expr, Xi_b: Xi_b_expr, Delta_r: Delta_r_expr, Delta_th: Delta_th_expr}) detg_a4 = detg_a4.subs({mu: mu_s})
detg_a4 = detg_a4.series(a, 5).truncate().simplify_full()
detg_a4.variables()
(μ0,a,,m,n,r)\renewcommand{\Bold}[1]{\mathbf{#1}}\left({\mu_0}, a, {\ell}, m, n, r\right)

For the time being, only the expansion at second order in aa is required:

detg_a2 = detg_a4.series(a, 3).truncate().simplify_full()

The Nambu-Goto Lagrangian at second order in aa:

L_a2 = (sqrt(-detg_a2)).series(a, 3).truncate().simplify_full() L_a2
(3μ02a24n23(μ021)a2422)r4(μ021)a2(4μ02a22mμ02a2)n2+((μ021)a24r8+2(μ021)a22r64(μ021)a2mr2+4(μ021)a2m2(4(μ021)a22m(μ021)a2)r4)rϕ1(r)2(μ02a24n2r8+2μ02a22n2r64μ02a2mn2r2+4μ02a2m2n2(4μ02a22mμ02a2)n2r4)rψ1(r)2+4((μ021)a22+1)m2r22(2r4+r22m)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{{\left(3 \, {\mu_0}^{2} a^{2} {\ell}^{4} n^{2} - 3 \, {\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{4} - 2 \, {\ell}^{2}\right)} r^{4} - {\left({\mu_0}^{2} - 1\right)} a^{2} - {\left(4 \, {\mu_0}^{2} a^{2} {\ell}^{2} m - {\mu_0}^{2} a^{2}\right)} n^{2} + {\left({\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{4} r^{8} + 2 \, {\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} r^{6} - 4 \, {\left({\mu_0}^{2} - 1\right)} a^{2} m r^{2} + 4 \, {\left({\mu_0}^{2} - 1\right)} a^{2} m^{2} - {\left(4 \, {\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} m - {\left({\mu_0}^{2} - 1\right)} a^{2}\right)} r^{4}\right)} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} - {\left({\mu_0}^{2} a^{2} {\ell}^{4} n^{2} r^{8} + 2 \, {\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{6} - 4 \, {\mu_0}^{2} a^{2} m n^{2} r^{2} + 4 \, {\mu_0}^{2} a^{2} m^{2} n^{2} - {\left(4 \, {\mu_0}^{2} a^{2} {\ell}^{2} m - {\mu_0}^{2} a^{2}\right)} n^{2} r^{4}\right)} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} + 4 \, {\left({\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} + 1\right)} m - 2 \, r^{2}}{2 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)}}
L_a2.numerator()
μ02a24n2r8rψ1(r)2μ02a24r8rϕ1(r)2+a24r8rϕ1(r)2+2μ02a22n2r6rψ1(r)24μ02a22mn2r4rψ1(r)23μ02a24n2r42μ02a22r6rϕ1(r)2+4μ02a22mr4rϕ1(r)2+3μ02a24r4+2a22r6rϕ1(r)2+μ02a2n2r4rψ1(r)24a22mr4rϕ1(r)24μ02a2mn2r2rψ1(r)23a24r4μ02a2r4rϕ1(r)2+4μ02a2m2n2rψ1(r)2+4μ02a22mn2+4μ02a2mr2rϕ1(r)24μ02a2m2rϕ1(r)2+a2r4rϕ1(r)24μ02a22m4a2mr2rϕ1(r)2μ02a2n2+22r4+4a2m2rϕ1(r)2+4a22m+μ02a2a2+2r24m\renewcommand{\Bold}[1]{\mathbf{#1}}{\mu_0}^{2} a^{2} {\ell}^{4} n^{2} r^{8} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} - {\mu_0}^{2} a^{2} {\ell}^{4} r^{8} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + a^{2} {\ell}^{4} r^{8} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + 2 \, {\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{6} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} - 4 \, {\mu_0}^{2} a^{2} {\ell}^{2} m n^{2} r^{4} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} - 3 \, {\mu_0}^{2} a^{2} {\ell}^{4} n^{2} r^{4} - 2 \, {\mu_0}^{2} a^{2} {\ell}^{2} r^{6} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + 4 \, {\mu_0}^{2} a^{2} {\ell}^{2} m r^{4} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + 3 \, {\mu_0}^{2} a^{2} {\ell}^{4} r^{4} + 2 \, a^{2} {\ell}^{2} r^{6} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + {\mu_0}^{2} a^{2} n^{2} r^{4} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} - 4 \, a^{2} {\ell}^{2} m r^{4} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} - 4 \, {\mu_0}^{2} a^{2} m n^{2} r^{2} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} - 3 \, a^{2} {\ell}^{4} r^{4} - {\mu_0}^{2} a^{2} r^{4} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + 4 \, {\mu_0}^{2} a^{2} m^{2} n^{2} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} + 4 \, {\mu_0}^{2} a^{2} {\ell}^{2} m n^{2} + 4 \, {\mu_0}^{2} a^{2} m r^{2} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} - 4 \, {\mu_0}^{2} a^{2} m^{2} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + a^{2} r^{4} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} - 4 \, {\mu_0}^{2} a^{2} {\ell}^{2} m - 4 \, a^{2} m r^{2} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} - {\mu_0}^{2} a^{2} n^{2} + 2 \, {\ell}^{2} r^{4} + 4 \, a^{2} m^{2} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + 4 \, a^{2} {\ell}^{2} m + {\mu_0}^{2} a^{2} - a^{2} + 2 \, r^{2} - 4 \, m
L_a2.denominator()
22r4+2r24m\renewcommand{\Bold}[1]{\mathbf{#1}}2 \, {\ell}^{2} r^{4} + 2 \, r^{2} - 4 \, m

Euler-Lagrange equations

def euler_lagrange(lagr, qs, var): r""" Derive the Euler-Lagrange equations from a given Lagrangian. INPUT: - ``lagr`` -- symbolic expression representing the Lagrangian density - ``qs`` -- either a single symbolic function or a list/tuple of symbolic functions, representing the `q`'s; these functions must appear in ``lagr`` up to at most their first derivatives - ``var`` -- either a single variable, typically `t` (1-dimensional problem) or a list/tuple of symbolic variables OUTPUT: - list of Euler-Lagrange equations; if only one function is involved, the single Euler-Lagrannge equation is returned instead. """ if not isinstance(qs, (list, tuple)): qs = [qs] if not isinstance(var, (list, tuple)): var = [var] n = len(qs) d = len(var) qv = [SR.var('qxxxx{}'.format(q)) for q in qs] dqv = [[SR.var('qxxxx{}_{}'.format(q, v)) for v in var] for q in qs] subs = {qs[i](*var): qv[i] for i in range(n)} subs_inv = {qv[i]: qs[i](*var) for i in range(n)} for i in range(n): subs.update({diff(qs[i](*var), var[j]): dqv[i][j] for j in range(d)}) subs_inv.update({dqv[i][j]: diff(qs[i](*var), var[j]) for j in range(d)}) lg = lagr.substitute(subs) eqs = [] for i in range(n): dLdq = diff(lg, qv[i]).simplify_full() dLdq = dLdq.substitute(subs_inv) ddL = 0 for j in range(d): h = diff(lg, dqv[i][j]).simplify_full() h = h.substitute(subs_inv) ddL += diff(h, var[j]) eqs.append((dLdq - ddL).simplify_full() == 0) if n == 1: return eqs[0] return eqs

We compute the Euler-Lagrange equations at order a2a^2 for ϕ1\phi_1 and ψ1\psi_1:

eqs = euler_lagrange(L_a2, [phi_1, psi_1], r) eqs
[2(2(μ021)a22r3+(μ021)a2r)rϕ1(r)+((μ021)a22r4+(μ021)a2r22(μ021)a2m)2(r)2ϕ1(r)=0,2(2μ02a22n2r3+μ02a2n2r)rψ1(r)(μ02a22n2r4+μ02a2n2r22μ02a2mn2)2(r)2ψ1(r)=0]\renewcommand{\Bold}[1]{\mathbf{#1}}\left[2 \, {\left(2 \, {\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} r^{3} + {\left({\mu_0}^{2} - 1\right)} a^{2} r\right)} \frac{\partial}{\partial r}\phi_{1}\left(r\right) + {\left({\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} r^{4} + {\left({\mu_0}^{2} - 1\right)} a^{2} r^{2} - 2 \, {\left({\mu_0}^{2} - 1\right)} a^{2} m\right)} \frac{\partial^{2}}{(\partial r)^{2}}\phi_{1}\left(r\right) = 0, -2 \, {\left(2 \, {\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{3} + {\mu_0}^{2} a^{2} n^{2} r\right)} \frac{\partial}{\partial r}\psi_{1}\left(r\right) - {\left({\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{4} + {\mu_0}^{2} a^{2} n^{2} r^{2} - 2 \, {\mu_0}^{2} a^{2} m n^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\psi_{1}\left(r\right) = 0\right]

Solving the equation for ϕ1\phi_1 (Eq. (4.8))

eq_phi1 = eqs[0] eq_phi1
2(2(μ021)a22r3+(μ021)a2r)rϕ1(r)+((μ021)a22r4+(μ021)a2r22(μ021)a2m)2(r)2ϕ1(r)=0\renewcommand{\Bold}[1]{\mathbf{#1}}2 \, {\left(2 \, {\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} r^{3} + {\left({\mu_0}^{2} - 1\right)} a^{2} r\right)} \frac{\partial}{\partial r}\phi_{1}\left(r\right) + {\left({\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} r^{4} + {\left({\mu_0}^{2} - 1\right)} a^{2} r^{2} - 2 \, {\left({\mu_0}^{2} - 1\right)} a^{2} m\right)} \frac{\partial^{2}}{(\partial r)^{2}}\phi_{1}\left(r\right) = 0
eq_phi1 = (eq_phi1/(a^2*(Mu0^2-1))).simplify_full() eq_phi1
2(22r3+r)rϕ1(r)+(2r4+r22m)2(r)2ϕ1(r)=0\renewcommand{\Bold}[1]{\mathbf{#1}}2 \, {\left(2 \, {\ell}^{2} r^{3} + r\right)} \frac{\partial}{\partial r}\phi_{1}\left(r\right) + {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} \frac{\partial^{2}}{(\partial r)^{2}}\phi_{1}\left(r\right) = 0
phi1_sol(r) = desolve(eq_phi1, phi_1(r), ivar=r) phi1_sol(r)
K112r4+r22mdr+K2\renewcommand{\Bold}[1]{\mathbf{#1}}K_{1} \int \frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}\,{d r} + K_{2}

We recover Eqs. (4.8) with K1=pK_1 = \mathfrak{p} and K2=0K_2=0.

The symbolic constants K1K_1 and K2K_2 are actually denoted _K1 and _K2 by SageMath, as the print reveals:

print(phi1_sol(r))
_K1*integrate(1/(l^2*r^4 + r^2 - 2*m), r) + _K2

Hence we perform the substitutions with SR.var('_K1') and SR.var('_K2'):

pf = var("pf", latex_name=r"\mathfrak{p}") phi1_sol(r) = phi1_sol(r).subs({SR.var('_K1'): pf, SR.var('_K2'): 0}) phi1_sol(r)
p12r4+r22mdr\renewcommand{\Bold}[1]{\mathbf{#1}}{\mathfrak{p}} \int \frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}\,{d r}

Solving the equation for ψ1\psi_1 (Eq. (4.8))

eq_psi1 = eqs[1] eq_psi1
2(2μ02a22n2r3+μ02a2n2r)rψ1(r)(μ02a22n2r4+μ02a2n2r22μ02a2mn2)2(r)2ψ1(r)=0\renewcommand{\Bold}[1]{\mathbf{#1}}-2 \, {\left(2 \, {\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{3} + {\mu_0}^{2} a^{2} n^{2} r\right)} \frac{\partial}{\partial r}\psi_{1}\left(r\right) - {\left({\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{4} + {\mu_0}^{2} a^{2} n^{2} r^{2} - 2 \, {\mu_0}^{2} a^{2} m n^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\psi_{1}\left(r\right) = 0
eq_psi1 = (eq_psi1/(a^2*Mu0^2)).simplify_full() eq_psi1
2(22n2r3+n2r)rψ1(r)(2n2r4+n2r22mn2)2(r)2ψ1(r)=0\renewcommand{\Bold}[1]{\mathbf{#1}}-2 \, {\left(2 \, {\ell}^{2} n^{2} r^{3} + n^{2} r\right)} \frac{\partial}{\partial r}\psi_{1}\left(r\right) - {\left({\ell}^{2} n^{2} r^{4} + n^{2} r^{2} - 2 \, m n^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\psi_{1}\left(r\right) = 0
psi1_sol(r) = desolve(eq_psi1, psi_1(r), ivar=r) psi1_sol(r)
K112r4+r22mdr+K2\renewcommand{\Bold}[1]{\mathbf{#1}}K_{1} \int \frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}\,{d r} + K_{2}

We recover Eq. (4.8) with K1=qK_1 = \mathfrak{q} and K2=0K_2=0.

qf = var('qf', latex_name=r"\mathfrak{q}") psi1_sol(r) = psi1_sol(r).subs({SR.var('_K1'): qf, SR.var('_K2'): 0}) psi1_sol(r)
q12r4+r22mdr\renewcommand{\Bold}[1]{\mathbf{#1}}{\mathfrak{q}} \int \frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}\,{d r}

Nambu-Goto Lagrangian at fourth order in aa

L_a4 = (sqrt(-detg_a4)).series(a, 5).truncate().simplify_full()
eqs = euler_lagrange(L_a4, [phi_1, psi_1, mu_1], r)

The equation for μ1\mu_1

eq_mu1 = eqs[2] eq_mu1
4(μ03μ0)a42m(μ03μ0)a4(4(μ03μ0)a42m(μ03μ0)a4)n4+3((μ03μ0)a44n4+2(μ03μ0)a44n32(μ03μ0)a44n(μ03μ0)a44)r42(4(μ03μ0)a42m(μ03μ0)a4)n3+(4(μ03μ0)a4m2n2+((μ03μ0)a44n2+2(μ03μ0)a44n+(μ03μ0)a44)r8+8(μ03μ0)a4m2n+4(μ03μ0)a4m2+2((μ03μ0)a42n2+2(μ03μ0)a42n+(μ03μ0)a42)r6(4(μ03μ0)a42m(μ03μ0)a4+(4(μ03μ0)a42m(μ03μ0)a4)n2+2(4(μ03μ0)a42m(μ03μ0)a4)n)r44((μ03μ0)a4mn2+2(μ03μ0)a4mn+(μ03μ0)a4m)r2)rϕ1(r)2(4(μ03μ0)a4m2n4+8(μ03μ0)a4m2n3+4(μ03μ0)a4m2n2+((μ03μ0)a44n4+2(μ03μ0)a44n3+(μ03μ0)a44n2)r8+2((μ03μ0)a42n4+2(μ03μ0)a42n3+(μ03μ0)a42n2)r6((4(μ03μ0)a42m(μ03μ0)a4)n4+2(4(μ03μ0)a42m(μ03μ0)a4)n3+(4(μ03μ0)a42m(μ03μ0)a4)n2)r44((μ03μ0)a4mn4+2(μ03μ0)a4mn3+(μ03μ0)a4mn2)r2)rψ1(r)2+2(4(μ03μ0)a42m(μ03μ0)a4)n2(2(a44n4+4a44n3+6a44n2+4a44n+a44)r7+3(a42n4+4a42n3+6a42n2+4a42n+a42)r5(4a42m+(4a42ma4)n4a4+4(4a42ma4)n3+6(4a42ma4)n2+4(4a42ma4)n)r32(a4mn4+4a4mn3+6a4mn2+4a4mn+a4m)r)rμ1(r)(4a4m2n4+16a4m2n3+(a44n4+4a44n3+6a44n2+4a44n+a44)r8+24a4m2n2+16a4m2n+2(a42n4+4a42n3+6a42n2+4a42n+a42)r6+4a4m2(4a42m+(4a42ma4)n4a4+4(4a42ma4)n3+6(4a42ma4)n2+4(4a42ma4)n)r44(a4mn4+4a4mn3+6a4mn2+4a4mn+a4m)r2)2(r)2μ1(r)(μ021)2r4+(μ021)r22(μ021)m=0\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} - {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n^{4} + 3 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n^{3} - 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4}\right)} r^{4} - 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n^{3} + {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m^{2} n^{2} + {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n + {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4}\right)} r^{8} + 8 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m^{2} n + 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m^{2} + 2 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} n + {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2}\right)} r^{6} - {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} + {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n^{2} + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n\right)} r^{4} - 4 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m n + {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m\right)} r^{2}\right)} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} - {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m^{2} n^{4} + 8 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m^{2} n^{3} + 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m^{2} n^{2} + {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{4} n^{2}\right)} r^{8} + 2 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} n^{2}\right)} r^{6} - {\left({\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n^{4} + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n^{3} + {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n^{2}\right)} r^{4} - 4 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} m n^{2}\right)} r^{2}\right)} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4} {\ell}^{2} m - {\left({\mu_0}^{3} - {\mu_0}\right)} a^{4}\right)} n - 2 \, {\left(2 \, {\left(a^{4} {\ell}^{4} n^{4} + 4 \, a^{4} {\ell}^{4} n^{3} + 6 \, a^{4} {\ell}^{4} n^{2} + 4 \, a^{4} {\ell}^{4} n + a^{4} {\ell}^{4}\right)} r^{7} + 3 \, {\left(a^{4} {\ell}^{2} n^{4} + 4 \, a^{4} {\ell}^{2} n^{3} + 6 \, a^{4} {\ell}^{2} n^{2} + 4 \, a^{4} {\ell}^{2} n + a^{4} {\ell}^{2}\right)} r^{5} - {\left(4 \, a^{4} {\ell}^{2} m + {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n^{4} - a^{4} + 4 \, {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n^{3} + 6 \, {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n^{2} + 4 \, {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n\right)} r^{3} - 2 \, {\left(a^{4} m n^{4} + 4 \, a^{4} m n^{3} + 6 \, a^{4} m n^{2} + 4 \, a^{4} m n + a^{4} m\right)} r\right)} \frac{\partial}{\partial r}\mu_{1}\left(r\right) - {\left(4 \, a^{4} m^{2} n^{4} + 16 \, a^{4} m^{2} n^{3} + {\left(a^{4} {\ell}^{4} n^{4} + 4 \, a^{4} {\ell}^{4} n^{3} + 6 \, a^{4} {\ell}^{4} n^{2} + 4 \, a^{4} {\ell}^{4} n + a^{4} {\ell}^{4}\right)} r^{8} + 24 \, a^{4} m^{2} n^{2} + 16 \, a^{4} m^{2} n + 2 \, {\left(a^{4} {\ell}^{2} n^{4} + 4 \, a^{4} {\ell}^{2} n^{3} + 6 \, a^{4} {\ell}^{2} n^{2} + 4 \, a^{4} {\ell}^{2} n + a^{4} {\ell}^{2}\right)} r^{6} + 4 \, a^{4} m^{2} - {\left(4 \, a^{4} {\ell}^{2} m + {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n^{4} - a^{4} + 4 \, {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n^{3} + 6 \, {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n^{2} + 4 \, {\left(4 \, a^{4} {\ell}^{2} m - a^{4}\right)} n\right)} r^{4} - 4 \, {\left(a^{4} m n^{4} + 4 \, a^{4} m n^{3} + 6 \, a^{4} m n^{2} + 4 \, a^{4} m n + a^{4} m\right)} r^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\mu_{1}\left(r\right)}{{\left({\mu_0}^{2} - 1\right)} {\ell}^{2} r^{4} + {\left({\mu_0}^{2} - 1\right)} r^{2} - 2 \, {\left({\mu_0}^{2} - 1\right)} m} = 0
eq_mu1.lhs().denominator().simplify_full()
(μ021)2r4+(μ021)r22(μ021)m\renewcommand{\Bold}[1]{\mathbf{#1}}{\left({\mu_0}^{2} - 1\right)} {\ell}^{2} r^{4} + {\left({\mu_0}^{2} - 1\right)} r^{2} - 2 \, {\left({\mu_0}^{2} - 1\right)} m
eq_mu1 = eq_mu1.lhs().numerator().simplify_full() == 0 #eq_mu1
eq_mu1 = (eq_mu1/a^4).simplify_full() eq_mu1
(4(μ03μ0)2mμ03+μ0)n43((μ03μ0)4n4+2(μ03μ0)4n32(μ03μ0)4n(μ03μ0)4)r44(μ03μ0)2m+2(4(μ03μ0)2mμ03+μ0)n3+μ03(((μ03μ0)4n2+2(μ03μ0)4n+(μ03μ0)4)r8+2((μ03μ0)2n2+2(μ03μ0)2n+(μ03μ0)2)r6+4(μ03μ0)m2n2(4(μ03μ0)2mμ03+(4(μ03μ0)2mμ03+μ0)n2+2(4(μ03μ0)2mμ03+μ0)n+μ0)r4+8(μ03μ0)m2n+4(μ03μ0)m24((μ03μ0)mn2+2(μ03μ0)mn+(μ03μ0)m)r2)rϕ1(r)2+(((μ03μ0)4n4+2(μ03μ0)4n3+(μ03μ0)4n2)r8+4(μ03μ0)m2n4+2((μ03μ0)2n4+2(μ03μ0)2n3+(μ03μ0)2n2)r6+8(μ03μ0)m2n3+4(μ03μ0)m2n2((4(μ03μ0)2mμ03+μ0)n4+2(4(μ03μ0)2mμ03+μ0)n3+(4(μ03μ0)2mμ03+μ0)n2)r44((μ03μ0)mn4+2(μ03μ0)mn3+(μ03μ0)mn2)r2)rψ1(r)22(4(μ03μ0)2mμ03+μ0)n+2(2(4n4+44n3+64n2+44n+4)r7+3(2n4+42n3+62n2+42n+2)r5((42m1)n4+4(42m1)n3+42m+6(42m1)n2+4(42m1)n1)r32(mn4+4mn3+6mn2+4mn+m)r)rμ1(r)+((4n4+44n3+64n2+44n+4)r8+2(2n4+42n3+62n2+42n+2)r6+4m2n4+16m2n3((42m1)n4+4(42m1)n3+42m+6(42m1)n2+4(42m1)n1)r4+24m2n2+16m2n4(mn4+4mn3+6mn2+4mn+m)r2+4m2)2(r)2μ1(r)μ0=0\renewcommand{\Bold}[1]{\mathbf{#1}}{\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{4} - 3 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{3} - 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n - {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4}\right)} r^{4} - 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{3} + {\mu_0}^{3} - {\left({\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n + {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4}\right)} r^{8} + 2 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} n + {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2}\right)} r^{6} + 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m^{2} n^{2} - {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{2} + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n + {\mu_0}\right)} r^{4} + 8 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m^{2} n + 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m^{2} - 4 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} m n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m n + {\left({\mu_0}^{3} - {\mu_0}\right)} m\right)} r^{2}\right)} \frac{\partial}{\partial r}\phi_{1}\left(r\right)^{2} + {\left({\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{2}\right)} r^{8} + 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m^{2} n^{4} + 2 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} n^{2}\right)} r^{6} + 8 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m^{2} n^{3} + 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m^{2} n^{2} - {\left({\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{4} + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{3} + {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{2}\right)} r^{4} - 4 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} m n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} m n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} m n^{2}\right)} r^{2}\right)} \frac{\partial}{\partial r}\psi_{1}\left(r\right)^{2} - 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n + 2 \, {\left(2 \, {\left({\ell}^{4} n^{4} + 4 \, {\ell}^{4} n^{3} + 6 \, {\ell}^{4} n^{2} + 4 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{7} + 3 \, {\left({\ell}^{2} n^{4} + 4 \, {\ell}^{2} n^{3} + 6 \, {\ell}^{2} n^{2} + 4 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{5} - {\left({\left(4 \, {\ell}^{2} m - 1\right)} n^{4} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{3} + 4 \, {\ell}^{2} m + 6 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{3} - 2 \, {\left(m n^{4} + 4 \, m n^{3} + 6 \, m n^{2} + 4 \, m n + m\right)} r\right)} \frac{\partial}{\partial r}\mu_{1}\left(r\right) + {\left({\left({\ell}^{4} n^{4} + 4 \, {\ell}^{4} n^{3} + 6 \, {\ell}^{4} n^{2} + 4 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{8} + 2 \, {\left({\ell}^{2} n^{4} + 4 \, {\ell}^{2} n^{3} + 6 \, {\ell}^{2} n^{2} + 4 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{6} + 4 \, m^{2} n^{4} + 16 \, m^{2} n^{3} - {\left({\left(4 \, {\ell}^{2} m - 1\right)} n^{4} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{3} + 4 \, {\ell}^{2} m + 6 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{4} + 24 \, m^{2} n^{2} + 16 \, m^{2} n - 4 \, {\left(m n^{4} + 4 \, m n^{3} + 6 \, m n^{2} + 4 \, m n + m\right)} r^{2} + 4 \, m^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\mu_{1}\left(r\right) - {\mu_0} = 0

We plug the solutions obtained previously for ϕ1(r)\phi_1(r) and ψ1(r)\psi_1(r) in this equation:

eq_mu1 = eq_mu1.substitute_function(phi_1, phi1_sol).substitute_function(psi_1, psi1_sol) eq_mu1 = eq_mu1.simplify_full() eq_mu1
(4(μ03μ0)2mμ03+μ0)n43((μ03μ0)4n4+2(μ03μ0)4n32(μ03μ0)4n(μ03μ0)4)r44(μ03μ0)2m+2(4(μ03μ0)2mμ03+μ0)n3+μ03(μ03+(μ03μ0)n2+2(μ03μ0)nμ0)p2+((μ03μ0)n4+2(μ03μ0)n3+(μ03μ0)n2)q22(4(μ03μ0)2mμ03+μ0)n+2(2(4n4+44n3+64n2+44n+4)r7+3(2n4+42n3+62n2+42n+2)r5((42m1)n4+4(42m1)n3+42m+6(42m1)n2+4(42m1)n1)r32(mn4+4mn3+6mn2+4mn+m)r)rμ1(r)+((4n4+44n3+64n2+44n+4)r8+2(2n4+42n3+62n2+42n+2)r6+4m2n4+16m2n3((42m1)n4+4(42m1)n3+42m+6(42m1)n2+4(42m1)n1)r4+24m2n2+16m2n4(mn4+4mn3+6mn2+4mn+m)r2+4m2)2(r)2μ1(r)μ0=0\renewcommand{\Bold}[1]{\mathbf{#1}}{\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{4} - 3 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{3} - 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n - {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4}\right)} r^{4} - 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{3} + {\mu_0}^{3} - {\left({\mu_0}^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} n - {\mu_0}\right)} {\mathfrak{p}}^{2} + {\left({\left({\mu_0}^{3} - {\mu_0}\right)} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} n^{2}\right)} {\mathfrak{q}}^{2} - 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n + 2 \, {\left(2 \, {\left({\ell}^{4} n^{4} + 4 \, {\ell}^{4} n^{3} + 6 \, {\ell}^{4} n^{2} + 4 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{7} + 3 \, {\left({\ell}^{2} n^{4} + 4 \, {\ell}^{2} n^{3} + 6 \, {\ell}^{2} n^{2} + 4 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{5} - {\left({\left(4 \, {\ell}^{2} m - 1\right)} n^{4} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{3} + 4 \, {\ell}^{2} m + 6 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{3} - 2 \, {\left(m n^{4} + 4 \, m n^{3} + 6 \, m n^{2} + 4 \, m n + m\right)} r\right)} \frac{\partial}{\partial r}\mu_{1}\left(r\right) + {\left({\left({\ell}^{4} n^{4} + 4 \, {\ell}^{4} n^{3} + 6 \, {\ell}^{4} n^{2} + 4 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{8} + 2 \, {\left({\ell}^{2} n^{4} + 4 \, {\ell}^{2} n^{3} + 6 \, {\ell}^{2} n^{2} + 4 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{6} + 4 \, m^{2} n^{4} + 16 \, m^{2} n^{3} - {\left({\left(4 \, {\ell}^{2} m - 1\right)} n^{4} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{3} + 4 \, {\ell}^{2} m + 6 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{4} + 24 \, m^{2} n^{2} + 16 \, m^{2} n - 4 \, {\left(m n^{4} + 4 \, m n^{3} + 6 \, m n^{2} + 4 \, m n + m\right)} r^{2} + 4 \, m^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\mu_{1}\left(r\right) - {\mu_0} = 0

Check of Eq. (4.9)

lhs = eq_mu1.lhs() lhs = lhs.simplify_full() lhs
(4(μ03μ0)2mμ03+μ0)n43((μ03μ0)4n4+2(μ03μ0)4n32(μ03μ0)4n(μ03μ0)4)r44(μ03μ0)2m+2(4(μ03μ0)2mμ03+μ0)n3+μ03(μ03+(μ03μ0)n2+2(μ03μ0)nμ0)p2+((μ03μ0)n4+2(μ03μ0)n3+(μ03μ0)n2)q22(4(μ03μ0)2mμ03+μ0)n+2(2(4n4+44n3+64n2+44n+4)r7+3(2n4+42n3+62n2+42n+2)r5((42m1)n4+4(42m1)n3+42m+6(42m1)n2+4(42m1)n1)r32(mn4+4mn3+6mn2+4mn+m)r)rμ1(r)+((4n4+44n3+64n2+44n+4)r8+2(2n4+42n3+62n2+42n+2)r6+4m2n4+16m2n3((42m1)n4+4(42m1)n3+42m+6(42m1)n2+4(42m1)n1)r4+24m2n2+16m2n4(mn4+4mn3+6mn2+4mn+m)r2+4m2)2(r)2μ1(r)μ0\renewcommand{\Bold}[1]{\mathbf{#1}}{\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{4} - 3 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{3} - 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n - {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4}\right)} r^{4} - 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m + 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{3} + {\mu_0}^{3} - {\left({\mu_0}^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} n^{2} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} n - {\mu_0}\right)} {\mathfrak{p}}^{2} + {\left({\left({\mu_0}^{3} - {\mu_0}\right)} n^{4} + 2 \, {\left({\mu_0}^{3} - {\mu_0}\right)} n^{3} + {\left({\mu_0}^{3} - {\mu_0}\right)} n^{2}\right)} {\mathfrak{q}}^{2} - 2 \, {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n + 2 \, {\left(2 \, {\left({\ell}^{4} n^{4} + 4 \, {\ell}^{4} n^{3} + 6 \, {\ell}^{4} n^{2} + 4 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{7} + 3 \, {\left({\ell}^{2} n^{4} + 4 \, {\ell}^{2} n^{3} + 6 \, {\ell}^{2} n^{2} + 4 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{5} - {\left({\left(4 \, {\ell}^{2} m - 1\right)} n^{4} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{3} + 4 \, {\ell}^{2} m + 6 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{3} - 2 \, {\left(m n^{4} + 4 \, m n^{3} + 6 \, m n^{2} + 4 \, m n + m\right)} r\right)} \frac{\partial}{\partial r}\mu_{1}\left(r\right) + {\left({\left({\ell}^{4} n^{4} + 4 \, {\ell}^{4} n^{3} + 6 \, {\ell}^{4} n^{2} + 4 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{8} + 2 \, {\left({\ell}^{2} n^{4} + 4 \, {\ell}^{2} n^{3} + 6 \, {\ell}^{2} n^{2} + 4 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{6} + 4 \, m^{2} n^{4} + 16 \, m^{2} n^{3} - {\left({\left(4 \, {\ell}^{2} m - 1\right)} n^{4} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{3} + 4 \, {\ell}^{2} m + 6 \, {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 4 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{4} + 24 \, m^{2} n^{2} + 16 \, m^{2} n - 4 \, {\left(m n^{4} + 4 \, m n^{3} + 6 \, m n^{2} + 4 \, m n + m\right)} r^{2} + 4 \, m^{2}\right)} \frac{\partial^{2}}{(\partial r)^{2}}\mu_{1}\left(r\right) - {\mu_0}
s = lhs.coefficient(diff(mu_1(r), r, 2)) # coefficient of mu_1'' s.factor()
(2r4+r22m)2(n+1)4\renewcommand{\Bold}[1]{\mathbf{#1}}{\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)}^{2} {\left(n + 1\right)}^{4}
s1 = (lhs/s - diff(mu_1(r), r, 2)).simplify_full() s1
(μ03μ0)n2q23((μ03μ0)4n2(μ03μ0)4)r44(μ03μ0)2m+μ03+(4(μ03μ0)2mμ03+μ0)n2(μ03μ0)p2+2(2(4n2+24n+4)r7+3(2n2+22n+2)r5(42m+(42m1)n2+2(42m1)n1)r32(mn2+2mn+m)r)rμ1(r)μ0(4n2+24n+4)r8+2(2n2+22n+2)r6(42m+(42m1)n2+2(42m1)n1)r4+4m2n2+8m2n4(mn2+2mn+m)r2+4m2\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{{\left({\mu_0}^{3} - {\mu_0}\right)} n^{2} {\mathfrak{q}}^{2} - 3 \, {\left({\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4} n^{2} - {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{4}\right)} r^{4} - 4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m + {\mu_0}^{3} + {\left(4 \, {\left({\mu_0}^{3} - {\mu_0}\right)} {\ell}^{2} m - {\mu_0}^{3} + {\mu_0}\right)} n^{2} - {\left({\mu_0}^{3} - {\mu_0}\right)} {\mathfrak{p}}^{2} + 2 \, {\left(2 \, {\left({\ell}^{4} n^{2} + 2 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{7} + 3 \, {\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{5} - {\left(4 \, {\ell}^{2} m + {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 2 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{3} - 2 \, {\left(m n^{2} + 2 \, m n + m\right)} r\right)} \frac{\partial}{\partial r}\mu_{1}\left(r\right) - {\mu_0}}{{\left({\ell}^{4} n^{2} + 2 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{8} + 2 \, {\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{6} - {\left(4 \, {\ell}^{2} m + {\left(4 \, {\ell}^{2} m - 1\right)} n^{2} + 2 \, {\left(4 \, {\ell}^{2} m - 1\right)} n - 1\right)} r^{4} + 4 \, m^{2} n^{2} + 8 \, m^{2} n - 4 \, {\left(m n^{2} + 2 \, m n + m\right)} r^{2} + 4 \, m^{2}}
b1 = s1.coefficient(diff(mu_1(r), r)).factor() b1
2(22r2+1)r2r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, {\left(2 \, {\ell}^{2} r^{2} + 1\right)} r}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}
b2 = (s1 - b1*diff(mu_1(r), r)).simplify_full().factor() b2
(34n2r434r442mn2n2q2+42m+n2+p21)(μ0+1)(μ01)μ0(2r4+r22m)2(n+1)2\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{{\left(3 \, {\ell}^{4} n^{2} r^{4} - 3 \, {\ell}^{4} r^{4} - 4 \, {\ell}^{2} m n^{2} - n^{2} {\mathfrak{q}}^{2} + 4 \, {\ell}^{2} m + n^{2} + {\mathfrak{p}}^{2} - 1\right)} {\left({\mu_0} + 1\right)} {\left({\mu_0} - 1\right)} {\mu_0}}{{\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)}^{2} {\left(n + 1\right)}^{2}}

The equation for μ1\mu_1 is thus:

eq_mu1 = diff(mu_1(r), r, 2) + b1*diff(mu_1(r), r) + b2 == 0 eq_mu1
2(22r2+1)rrμ1(r)2r4+r22m(34n2r434r442mn2n2q2+42m+n2+p21)(μ0+1)(μ01)μ0(2r4+r22m)2(n+1)2+2(r)2μ1(r)=0\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, {\left(2 \, {\ell}^{2} r^{2} + 1\right)} r \frac{\partial}{\partial r}\mu_{1}\left(r\right)}{{\ell}^{2} r^{4} + r^{2} - 2 \, m} - \frac{{\left(3 \, {\ell}^{4} n^{2} r^{4} - 3 \, {\ell}^{4} r^{4} - 4 \, {\ell}^{2} m n^{2} - n^{2} {\mathfrak{q}}^{2} + 4 \, {\ell}^{2} m + n^{2} + {\mathfrak{p}}^{2} - 1\right)} {\left({\mu_0} + 1\right)} {\left({\mu_0} - 1\right)} {\mu_0}}{{\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)}^{2} {\left(n + 1\right)}^{2}} + \frac{\partial^{2}}{(\partial r)^{2}}\mu_{1}\left(r\right) = 0

Let us define Θ2:=2Θ0 \Theta_2 := 2 \Theta_0

Th2 = var('Th2', latex_name=r'\Theta_2', domain='real')

Given that μ1(r)=sinΘ0  θ1(r)=1μ02  θ1(r) \mu_1(r) = - \sin\Theta_0 \; \theta_1(r) = - \sqrt{1-\mu_0^2} \; \theta_1(r) and sin2Θ0=2μ01μ02,\sin2\Theta_0 = 2\mu_0\sqrt{1-\mu_0^2}, we get the equation for Υ:=θ1=μ11μ0\Upsilon := \theta_1' = - \frac{\mu_1'}{\sqrt{1 - \mu_0}}:

Y = function('Y', latex_name=r'\Upsilon') eq_Y = diff(Y(r), r) + b1*Y(r) \ - (b2/(2*(1-Mu0^2)*Mu0)*sin(Th2)).factor() == 0 eq_Y
2(22r2+1)rΥ(r)2r4+r22m(34n2r434r442mn2n2q2+42m+n2+p21)sin(Θ2)2(2r4+r22m)2(n+1)2+rΥ(r)=0\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, {\left(2 \, {\ell}^{2} r^{2} + 1\right)} r \Upsilon\left(r\right)}{{\ell}^{2} r^{4} + r^{2} - 2 \, m} - \frac{{\left(3 \, {\ell}^{4} n^{2} r^{4} - 3 \, {\ell}^{4} r^{4} - 4 \, {\ell}^{2} m n^{2} - n^{2} {\mathfrak{q}}^{2} + 4 \, {\ell}^{2} m + n^{2} + {\mathfrak{p}}^{2} - 1\right)} \sin\left({\Theta_2}\right)}{2 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)}^{2} {\left(n + 1\right)}^{2}} + \frac{\partial}{\partial r}\Upsilon\left(r\right) = 0

This agrees with Eq. (4.9) of the paper.

Solving the equation for Υ:=θ1\Upsilon := \theta_1'

We use the function desolve to solve the differential equation for Υ\Upsilon:

Y_sol(r) = desolve(eq_Y, Y(r), ivar=r)

The solution involves an integral that SageMath is not capable to evaluate with the default integrator. Trying to display Y_sol would make SageMath hang. Instead, we print Y_sol to get the unvaluated form of the integral, in order to compute it by means of FriCAS:

print(Y_sol(r))
1/2*(2*_C + 3*(l^2*n*sin(Th2) - l^2*sin(Th2))*r/(n + 1) - integrate((n^2*qf^2 + 2*l^2*m - (2*l^2*m + 1)*n^2 + 3*(l^2*n^2 - l^2)*r^2 - pf^2 + 1)/(l^2*r^4 + r^2 - 2*m), r)*sin(Th2)/(n^2 + 2*n + 1))/(l^2*r^4 + r^2 - 2*m)

The solution involves some constant, denoted _C by SageMath. We rename it C_1 and rewrite the above solution as

C_1 = var('C_1') Integ(r) = function('Integ')(r) Y_sol0(r) = 1/2*(2*C_1 + 3*(l^2*n*sin(Th2) - l^2*sin(Th2))*r/(n + 1) \ - Integ(r)*sin(Th2)/(n^2 + 2*n + 1))/(l^2*r^4 + r^2 - 2*m) Y_sol0(r)
2C1+3(2nsin(Θ2)2sin(Θ2))rn+1Integ(r)sin(Θ2)n2+2n+12(2r4+r22m)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, C_{1} + \frac{3 \, {\left({\ell}^{2} n \sin\left({\Theta_2}\right) - {\ell}^{2} \sin\left({\Theta_2}\right)\right)} r}{n + 1} - \frac{{\rm Integ}\left(r\right) \sin\left({\Theta_2}\right)}{n^{2} + 2 \, n + 1}}{2 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)}}

Integ(r) represents the integral I(r)I(r), whose integrand, F(r)F(r) say, is read from the output of print(Y_sol(r)):

F(r) = (n^2*qf^2 + 2*l^2*m - (2*l^2*m + 1)*n^2 + 3*(l^2*n^2 - l^2)*r^2 - pf^2 + 1) \ /(l^2*r^4 + r^2 - 2*m) F(r)
n2q2+22m(22m+1)n2+3(2n22)r2p2+12r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{2} m - {\left(2 \, {\ell}^{2} m + 1\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} r^{2} - {\mathfrak{p}}^{2} + 1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}

We split the integral in two parts: I(r)=F1  s1(r)+F2  s2(r) I(r) = F_1 \; s_1(r) + F_2 \; s_2(r) with s1(r):=rrˉ22rˉ4+rˉ22mdrˉ,s2(r):=rdrˉ2rˉ4+rˉ22m s_1(r) := \int^r \frac{\bar{r}^2}{\ell^2 \bar{r}^4 + \bar{r}^2 - 2m} \, \mathrm{d}\bar{r}, \qquad s_2(r) := \int^r \frac{\mathrm{d}\bar{r}}{\ell^2 \bar{r}^4 + \bar{r}^2 - 2m} and

F1 = 3*(l^2*n^2 - l^2) F1
32n232\renewcommand{\Bold}[1]{\mathbf{#1}}3 \, {\ell}^{2} n^{2} - 3 \, {\ell}^{2}
F2 = n^2*qf^2 + 2*l^2*m - (2*l^2*m + 1)*n^2 - pf^2 + 1 F2
n2q2+22m(22m+1)n2p2+1\renewcommand{\Bold}[1]{\mathbf{#1}}n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{2} m - {\left(2 \, {\ell}^{2} m + 1\right)} n^{2} - {\mathfrak{p}}^{2} + 1

Check:

bool(F(r) == F1*r^2/(l^2*r^4 + r^2 - 2*m) + F2/(l^2*r^4 + r^2 - 2*m))
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}

Let us evaluate s1(r)s_1(r) by means of FriCAS:

s1 = integrate(r^2/(l^2*r^4 + r^2 - 2*m), r, algorithm='fricas') s1
121284m+286m+4+184m+2log(12(84m+2)84m+286m+4+184m+286m+4+r)121284m+286m+4+184m+2log(12(84m+2)84m+286m+4+184m+286m+4+r)121284m+286m+4184m+2log(12(84m+2)84m+286m+4184m+286m+4+r)+121284m+286m+4184m+2log(12(84m+2)84m+286m+4184m+286m+4+r)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, \sqrt{\frac{1}{2}} \sqrt{-\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + 1}{8 \, {\ell}^{4} m + {\ell}^{2}}} \log\left(\frac{\sqrt{\frac{1}{2}} {\left(8 \, {\ell}^{4} m + {\ell}^{2}\right)} \sqrt{-\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + 1}{8 \, {\ell}^{4} m + {\ell}^{2}}}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + r\right) - \frac{1}{2} \, \sqrt{\frac{1}{2}} \sqrt{-\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + 1}{8 \, {\ell}^{4} m + {\ell}^{2}}} \log\left(-\frac{\sqrt{\frac{1}{2}} {\left(8 \, {\ell}^{4} m + {\ell}^{2}\right)} \sqrt{-\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + 1}{8 \, {\ell}^{4} m + {\ell}^{2}}}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + r\right) - \frac{1}{2} \, \sqrt{\frac{1}{2}} \sqrt{\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} - 1}{8 \, {\ell}^{4} m + {\ell}^{2}}} \log\left(\frac{\sqrt{\frac{1}{2}} {\left(8 \, {\ell}^{4} m + {\ell}^{2}\right)} \sqrt{\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} - 1}{8 \, {\ell}^{4} m + {\ell}^{2}}}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + r\right) + \frac{1}{2} \, \sqrt{\frac{1}{2}} \sqrt{\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} - 1}{8 \, {\ell}^{4} m + {\ell}^{2}}} \log\left(-\frac{\sqrt{\frac{1}{2}} {\left(8 \, {\ell}^{4} m + {\ell}^{2}\right)} \sqrt{\frac{\frac{8 \, {\ell}^{4} m + {\ell}^{2}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} - 1}{8 \, {\ell}^{4} m + {\ell}^{2}}}}{\sqrt{8 \, {\ell}^{6} m + {\ell}^{4}}} + r\right)
s1 = s1.canonicalize_radical().simplify_log() s1
282m82m+1+1log(2(82m+1)14r82m82m+1+12(82m+1)14r+82m82m+1+1)+282m82m+11log(2(82m+1)14r+82m82m+112(82m+1)14r82m82m+11)4(82m+1)34\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{\sqrt{2} \sqrt{8 \, {\ell}^{2} m - \sqrt{8 \, {\ell}^{2} m + 1} + 1} \log\left(\frac{\sqrt{2} {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell} r - \sqrt{8 \, {\ell}^{2} m - \sqrt{8 \, {\ell}^{2} m + 1} + 1}}{\sqrt{2} {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell} r + \sqrt{8 \, {\ell}^{2} m - \sqrt{8 \, {\ell}^{2} m + 1} + 1}}\right) + \sqrt{2} \sqrt{-8 \, {\ell}^{2} m - \sqrt{8 \, {\ell}^{2} m + 1} - 1} \log\left(\frac{\sqrt{2} {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell} r + \sqrt{-8 \, {\ell}^{2} m - \sqrt{8 \, {\ell}^{2} m + 1} - 1}}{\sqrt{2} {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell} r - \sqrt{-8 \, {\ell}^{2} m - \sqrt{8 \, {\ell}^{2} m + 1} - 1}}\right)}{4 \, {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{3}{4}} {\ell}}

Check:

diff(s1, r).simplify_full()
r22r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{r^{2}}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}

Similarly, we evaluate s2(r)s_2(r) by means of FriCAS:

s2 = integrate(1/(l^2*r^4 + r^2 - 2*m), r, algorithm='fricas') s2
1482m2+m82m3+m2+182m2+mlog(22r+12(82m82m2+m82m3+m2+1)82m2+m82m3+m2+182m2+m)+1482m2+m82m3+m2+182m2+mlog(22r12(82m82m2+m82m3+m2+1)82m2+m82m3+m2+182m2+m)1482m2+m82m3+m2182m2+mlog(22r+12(82m+82m2+m82m3+m2+1)82m2+m82m3+m2182m2+m)+1482m2+m82m3+m2182m2+mlog(22r12(82m+82m2+m82m3+m2+1)82m2+m82m3+m2182m2+m)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{1}{4} \, \sqrt{\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1}{8 \, {\ell}^{2} m^{2} + m}} \log\left(2 \, {\ell}^{2} r + \frac{1}{2} \, {\left(8 \, {\ell}^{2} m - \frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1\right)} \sqrt{\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1}{8 \, {\ell}^{2} m^{2} + m}}\right) + \frac{1}{4} \, \sqrt{\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1}{8 \, {\ell}^{2} m^{2} + m}} \log\left(2 \, {\ell}^{2} r - \frac{1}{2} \, {\left(8 \, {\ell}^{2} m - \frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1\right)} \sqrt{\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1}{8 \, {\ell}^{2} m^{2} + m}}\right) - \frac{1}{4} \, \sqrt{-\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} - 1}{8 \, {\ell}^{2} m^{2} + m}} \log\left(2 \, {\ell}^{2} r + \frac{1}{2} \, {\left(8 \, {\ell}^{2} m + \frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1\right)} \sqrt{-\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} - 1}{8 \, {\ell}^{2} m^{2} + m}}\right) + \frac{1}{4} \, \sqrt{-\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} - 1}{8 \, {\ell}^{2} m^{2} + m}} \log\left(2 \, {\ell}^{2} r - \frac{1}{2} \, {\left(8 \, {\ell}^{2} m + \frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} + 1\right)} \sqrt{-\frac{\frac{8 \, {\ell}^{2} m^{2} + m}{\sqrt{8 \, {\ell}^{2} m^{3} + m^{2}}} - 1}{8 \, {\ell}^{2} m^{2} + m}}\right)
s2 = s2.canonicalize_radical().simplify_log() s2
82m+82m+11log(4(82m+1)142mr82m+82m+11(82m+1+1)4(82m+1)142mr+82m+82m+11(82m+1+1))+82m+82m+1+1log(4(82m+1)142mr82m+82m+1+1(82m+11)4(82m+1)142mr+82m+82m+1+1(82m+11))4(82m+1)34m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{\sqrt{-8 \, {\ell}^{2} m + \sqrt{8 \, {\ell}^{2} m + 1} - 1} \log\left(\frac{4 \, {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell}^{2} \sqrt{m} r - \sqrt{-8 \, {\ell}^{2} m + \sqrt{8 \, {\ell}^{2} m + 1} - 1} {\left(\sqrt{8 \, {\ell}^{2} m + 1} + 1\right)}}{4 \, {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell}^{2} \sqrt{m} r + \sqrt{-8 \, {\ell}^{2} m + \sqrt{8 \, {\ell}^{2} m + 1} - 1} {\left(\sqrt{8 \, {\ell}^{2} m + 1} + 1\right)}}\right) + \sqrt{8 \, {\ell}^{2} m + \sqrt{8 \, {\ell}^{2} m + 1} + 1} \log\left(\frac{4 \, {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell}^{2} \sqrt{m} r - \sqrt{8 \, {\ell}^{2} m + \sqrt{8 \, {\ell}^{2} m + 1} + 1} {\left(\sqrt{8 \, {\ell}^{2} m + 1} - 1\right)}}{4 \, {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{1}{4}} {\ell}^{2} \sqrt{m} r + \sqrt{8 \, {\ell}^{2} m + \sqrt{8 \, {\ell}^{2} m + 1} + 1} {\left(\sqrt{8 \, {\ell}^{2} m + 1} - 1\right)}}\right)}{4 \, {\left(8 \, {\ell}^{2} m + 1\right)}^{\frac{3}{4}} \sqrt{m}}

Check:

diff(s2, r).simplify_full()
12r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}

In the above expressions for s1(r)s_1(r) and s2(r)s_2(r) there appears 1+82m\sqrt{1 + 8 \ell^2 m}, which can be rewritten 1+82m=22rH2+1 \sqrt{1 + 8 \ell^2 m} = 2 \ell^2 r_H^2 + 1 where rHr_H is the positive root of 2rH4+rH22m=0\ell^2 r_H^4 + r_H^2 - 2m = 0. More precisely, we perform the following substitution: m=12rH2(2rH2+1) m = \frac{1}{2} r_H^2 (\ell^2 r_H^2 + 1)

rH = var('rH', latex_name=r'r_H', domain='real') assume(rH > 0) m_rH = rH^2*(l^2*rH^2 + 1)/2 s1 = s1.subs({m: m_rH}).canonicalize_radical().simplify_log() s1
rHlog(r+rHrrH)+i2rH2+1log(ri2rH2+1r+i2rH2+1)2(23rH2+)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{{\ell} {r_H} \log\left(\frac{r + {r_H}}{r - {r_H}}\right) + i \, \sqrt{{\ell}^{2} {r_H}^{2} + 1} \log\left(\frac{{\ell} r - i \, \sqrt{{\ell}^{2} {r_H}^{2} + 1}}{{\ell} r + i \, \sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right)}{2 \, {\left(2 \, {\ell}^{3} {r_H}^{2} + {\ell}\right)}}

In the second log\log, we recognize the arccot\mathrm{arccot} function, via the identity arccotx=i2ln(xix+i). \mathrm{arccot}\, x = \frac{i}{2} \ln\left( \frac{x - i}{x + i} \right) . Given that arccotx=π/2arctanx\mathrm{arccot}\, x = \pi/2 - \mathrm{arctan}\, x, we use this identity as iln(xix+i)=2arctan(x)+π i \ln\left( \frac{x - i}{x + i} \right) = - 2 \, \mathrm{arctan}(x) + \pi

Thus, we perform the following substitution, disregarding the additive constant π\pi:

s1 = s1.subs({I*sqrt(l^2*rH^2 + 1)*log((l*r - I*sqrt(l^2*rH^2 + 1))/(l*r + I*sqrt(l^2*rH^2 + 1))): -2*sqrt(l^2*rH^2 + 1)*atan(l*r/sqrt(l^2*rH^2 + 1))}) s1
rHlog(r+rHrrH)22rH2+1arctan(r2rH2+1)2(23rH2+)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{{\ell} {r_H} \log\left(\frac{r + {r_H}}{r - {r_H}}\right) - 2 \, \sqrt{{\ell}^{2} {r_H}^{2} + 1} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right)}{2 \, {\left(2 \, {\ell}^{3} {r_H}^{2} + {\ell}\right)}}

Let us check that we have indeed a primitive of rr22r4+r22mr\mapsto \frac{r^2}{\ell^2 r^4 + r^2 - 2m}:

Ds1 = diff(s1, r).simplify_full() Ds1
r22r42rH4+r2rH2\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{r^{2}}{{\ell}^{2} r^{4} - {\ell}^{2} {r_H}^{4} + r^{2} - {r_H}^{2}}
rH_m = sqrt(sqrt(1 + 8*l^2*m) - 1)/(sqrt(2)*l) Ds1.subs({rH: rH_m}).simplify_full()
r22r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{r^{2}}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}

Similarly, let us express s2s_2 in terms of rHr_H:

s2 = s2.subs({m: m_rH}).canonicalize_radical().simplify_log() s2
irHlog(i2rH2+2rH2+1rii2rH2+2rH2+1r+i)+2rH2+1log(rrHr+rH)2(22rH3+rH)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{i \, {\ell} {r_H} \log\left(\frac{-i \, {\ell}^{2} {r_H}^{2} + \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} r - i}{i \, {\ell}^{2} {r_H}^{2} + \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} r + i}\right) + \sqrt{{\ell}^{2} {r_H}^{2} + 1} \log\left(\frac{r - {r_H}}{r + {r_H}}\right)}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{3} + {r_H}\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

Again, we use the identity iln(xix+i)=2arctan(x)+π i \ln\left( \frac{x - i}{x + i} \right) = - 2 \, \mathrm{arctan}(x) + \pi to rewrite s2s_2 as

s2 = s2.subs({I*l*rH*log((-I*l^2*rH^2 + sqrt(l^2*rH^2 + 1)*l*r - I)/(I*l^2*rH^2 + sqrt(l^2*rH^2 + 1)*l*r + I)): -2*l*rH*atan(l*r/sqrt(l^2*rH^2 + 1))}) s2
2rHarctan(r2rH2+1)2rH2+1log(rrHr+rH)2(22rH3+rH)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{2 \, {\ell} {r_H} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right) - \sqrt{{\ell}^{2} {r_H}^{2} + 1} \log\left(\frac{r - {r_H}}{r + {r_H}}\right)}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{3} + {r_H}\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

Let us also replace ln(rrHr+rH)\ln\left(\frac{r - r_H}{r + r_H}\right) by ln(r+rHrrH)-\ln\left(\frac{r + r_H}{r - r_H}\right) in order to have the same log term as in s1(r)s_1(r):

s2 = s2.subs({log((r - rH)/(r + rH)): - log((r + rH)/(r - rH))}) s2
2rHarctan(r2rH2+1)+2rH2+1log(r+rHrrH)2(22rH3+rH)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{2 \, {\ell} {r_H} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right) + \sqrt{{\ell}^{2} {r_H}^{2} + 1} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{3} + {r_H}\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

Let us check that we have indeed a primitive of r12r4+r22mr\mapsto \frac{1}{\ell^2 r^4 + r^2 - 2m}:

Ds2 = diff(s2, r).simplify_full() Ds2
12r42rH4+r2rH2\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{{\ell}^{2} r^{4} - {\ell}^{2} {r_H}^{4} + r^{2} - {r_H}^{2}}
Ds2.subs({rH: rH_m}).simplify_full()
12r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}

The full integral is thus

Integ0 = (F1*s1 + F2*s2).simplify_full() Integ0
(n2q2+22m(22m+1)n2+3(2n22)rH2p2+1)2rH2+1log(r+rHrrH)2(3(3n23)rH3(n2q2+23m2(3m+2)n2p2+4)rH)arctan(r2rH2+1)2(22rH3+rH)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{{\left(n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{2} m - {\left(2 \, {\ell}^{2} m + 1\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} {r_H}^{2} - {\mathfrak{p}}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1} \log\left(\frac{r + {r_H}}{r - {r_H}}\right) - 2 \, {\left(3 \, {\left({\ell}^{3} n^{2} - {\ell}^{3}\right)} {r_H}^{3} - {\left({\ell} n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{3} m - 2 \, {\left({\ell}^{3} m + 2 \, {\ell}\right)} n^{2} - {\ell} {\mathfrak{p}}^{2} + 4 \, {\ell}\right)} {r_H}\right)} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right)}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{3} + {r_H}\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

so that the solution is

Y_sol(r) = Y_sol0(r).subs({Integ(r): Integ0}).simplify_full() Y_sol(r)
2(3(3n2sin(Θ2)3sin(Θ2))rH3(n2q2sin(Θ2)+23msin(Θ2)p2sin(Θ2)2(3msin(Θ2)+2sin(Θ2))n2+4sin(Θ2))rH)arctan(r2rH2+1)2rH2+1(4(2C12n2+4C12n+2C12+3(4n2sin(Θ2)4sin(Θ2))r)rH3+2(2C1n2+4C1n+3(2n2sin(Θ2)2sin(Θ2))r+2C1)rH+(n2q2sin(Θ2)+22msin(Θ2)(22msin(Θ2)+sin(Θ2))n2+3(2n2sin(Θ2)2sin(Θ2))rH2p2sin(Θ2)+sin(Θ2))log(r+rHrrH))42rH2+1(2(22mn2(4n2+24n+4)r4+42mn+22m(2n2+22n+2)r2)rH3((2n2+22n+2)r42mn2+(n2+2n+1)r24mn2m)rH)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, {\left(3 \, {\left({\ell}^{3} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{3} \sin\left({\Theta_2}\right)\right)} {r_H}^{3} - {\left({\ell} n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, {\ell}^{3} m \sin\left({\Theta_2}\right) - {\ell} {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) - 2 \, {\left({\ell}^{3} m \sin\left({\Theta_2}\right) + 2 \, {\ell} \sin\left({\Theta_2}\right)\right)} n^{2} + 4 \, {\ell} \sin\left({\Theta_2}\right)\right)} {r_H}\right)} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right) - \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(4 \, {\left(2 \, C_{1} {\ell}^{2} n^{2} + 4 \, C_{1} {\ell}^{2} n + 2 \, C_{1} {\ell}^{2} + 3 \, {\left({\ell}^{4} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{4} \sin\left({\Theta_2}\right)\right)} r\right)} {r_H}^{3} + 2 \, {\left(2 \, C_{1} n^{2} + 4 \, C_{1} n + 3 \, {\left({\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{2} \sin\left({\Theta_2}\right)\right)} r + 2 \, C_{1}\right)} {r_H} + {\left(n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, {\ell}^{2} m \sin\left({\Theta_2}\right) - {\left(2 \, {\ell}^{2} m \sin\left({\Theta_2}\right) + \sin\left({\Theta_2}\right)\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{2} \sin\left({\Theta_2}\right)\right)} {r_H}^{2} - {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) + \sin\left({\Theta_2}\right)\right)} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)\right)}}{4 \, \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(2 \, {\left(2 \, {\ell}^{2} m n^{2} - {\left({\ell}^{4} n^{2} + 2 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{4} + 4 \, {\ell}^{2} m n + 2 \, {\ell}^{2} m - {\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{2}\right)} {r_H}^{3} - {\left({\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{4} - 2 \, m n^{2} + {\left(n^{2} + 2 \, n + 1\right)} r^{2} - 4 \, m n - 2 \, m\right)} {r_H}\right)}}
Y_sol(r).numerator().simplify_full()
2(3(3n2sin(Θ2)3sin(Θ2))rH3(n2q2sin(Θ2)+23msin(Θ2)p2sin(Θ2)2(3msin(Θ2)+2sin(Θ2))n2+4sin(Θ2))rH)arctan(r2rH2+1)+2rH2+1(4(2C12n2+4C12n+2C12+3(4n2sin(Θ2)4sin(Θ2))r)rH3+2(2C1n2+4C1n+3(2n2sin(Θ2)2sin(Θ2))r+2C1)rH+(n2q2sin(Θ2)+22msin(Θ2)(22msin(Θ2)+sin(Θ2))n2+3(2n2sin(Θ2)2sin(Θ2))rH2p2sin(Θ2)+sin(Θ2))log(r+rHrrH))\renewcommand{\Bold}[1]{\mathbf{#1}}-2 \, {\left(3 \, {\left({\ell}^{3} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{3} \sin\left({\Theta_2}\right)\right)} {r_H}^{3} - {\left({\ell} n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, {\ell}^{3} m \sin\left({\Theta_2}\right) - {\ell} {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) - 2 \, {\left({\ell}^{3} m \sin\left({\Theta_2}\right) + 2 \, {\ell} \sin\left({\Theta_2}\right)\right)} n^{2} + 4 \, {\ell} \sin\left({\Theta_2}\right)\right)} {r_H}\right)} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right) + \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(4 \, {\left(2 \, C_{1} {\ell}^{2} n^{2} + 4 \, C_{1} {\ell}^{2} n + 2 \, C_{1} {\ell}^{2} + 3 \, {\left({\ell}^{4} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{4} \sin\left({\Theta_2}\right)\right)} r\right)} {r_H}^{3} + 2 \, {\left(2 \, C_{1} n^{2} + 4 \, C_{1} n + 3 \, {\left({\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{2} \sin\left({\Theta_2}\right)\right)} r + 2 \, C_{1}\right)} {r_H} + {\left(n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, {\ell}^{2} m \sin\left({\Theta_2}\right) - {\left(2 \, {\ell}^{2} m \sin\left({\Theta_2}\right) + \sin\left({\Theta_2}\right)\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - {\ell}^{2} \sin\left({\Theta_2}\right)\right)} {r_H}^{2} - {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) + \sin\left({\Theta_2}\right)\right)} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)\right)}
Y_sol(r).denominator().factor()
4(2r4+r22m)(22rH2+1)2rH2+1(n+1)2rH\renewcommand{\Bold}[1]{\mathbf{#1}}4 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(n + 1\right)}^{2} {r_H}

Let us check that Y_sol is indeed a solution of the differential equation for Υ\Upsilon:

eq_Y.substitute_function(Y, Y_sol).subs({rH: rH_m}).simplify_full()
0=0\renewcommand{\Bold}[1]{\mathbf{#1}}0 = 0
print(Y_sol(r))
1/4*(2*(3*(l^3*n^2*sin(Th2) - l^3*sin(Th2))*rH^3 - (l*n^2*qf^2*sin(Th2) + 2*l^3*m*sin(Th2) - l*pf^2*sin(Th2) - 2*(l^3*m*sin(Th2) + 2*l*sin(Th2))*n^2 + 4*l*sin(Th2))*rH)*arctan(l*r/sqrt(l^2*rH^2 + 1)) - sqrt(l^2*rH^2 + 1)*(4*(2*C_1*l^2*n^2 + 4*C_1*l^2*n + 2*C_1*l^2 + 3*(l^4*n^2*sin(Th2) - l^4*sin(Th2))*r)*rH^3 + 2*(2*C_1*n^2 + 4*C_1*n + 3*(l^2*n^2*sin(Th2) - l^2*sin(Th2))*r + 2*C_1)*rH + (n^2*qf^2*sin(Th2) + 2*l^2*m*sin(Th2) - (2*l^2*m*sin(Th2) + sin(Th2))*n^2 + 3*(l^2*n^2*sin(Th2) - l^2*sin(Th2))*rH^2 - pf^2*sin(Th2) + sin(Th2))*log((r + rH)/(r - rH))))/(sqrt(l^2*rH^2 + 1)*(2*(2*l^2*m*n^2 - (l^4*n^2 + 2*l^4*n + l^4)*r^4 + 4*l^2*m*n + 2*l^2*m - (l^2*n^2 + 2*l^2*n + l^2)*r^2)*rH^3 - ((l^2*n^2 + 2*l^2*n + l^2)*r^4 - 2*m*n^2 + (n^2 + 2*n + 1)*r^2 - 4*m*n - 2*m)*rH))

Check of Eq. (4.10) (expression of θ1=Υ\theta'_1 = \Upsilon)

The term involving the constant C1C_1 agrees with that of Eq. (4.10):

s = Y_sol(r).coefficient(C_1).simplify_full() s
12r4+r22m\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{{\ell}^{2} r^{4} + r^{2} - 2 \, m}
fr4 = l^2*r^4 + r^2 - 2*m bool(s == 1/fr4)
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}

Let us remove it from Υ\Upsilon and divide the result by sin(2Θ0)\sin(2\Theta_0):

Y1 = ((Y_sol(r) - s*C_1)/sin(Th2)).simplify_full() Y1
2(3(3n23)rH3(n2q2+23m2(3m+2)n2p2+4)rH)arctan(r2rH2+1)(12(4n24)rrH3+6(2n22)rrH+(n2q2+22m(22m+1)n2+3(2n22)rH2p2+1)log(r+rHrrH))2rH2+142rH2+1(2(22mn2(4n2+24n+4)r4+42mn+22m(2n2+22n+2)r2)rH3((2n2+22n+2)r42mn2+(n2+2n+1)r24mn2m)rH)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, {\left(3 \, {\left({\ell}^{3} n^{2} - {\ell}^{3}\right)} {r_H}^{3} - {\left({\ell} n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{3} m - 2 \, {\left({\ell}^{3} m + 2 \, {\ell}\right)} n^{2} - {\ell} {\mathfrak{p}}^{2} + 4 \, {\ell}\right)} {r_H}\right)} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right) - {\left(12 \, {\left({\ell}^{4} n^{2} - {\ell}^{4}\right)} r {r_H}^{3} + 6 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} r {r_H} + {\left(n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{2} m - {\left(2 \, {\ell}^{2} m + 1\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} {r_H}^{2} - {\mathfrak{p}}^{2} + 1\right)} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}{4 \, \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(2 \, {\left(2 \, {\ell}^{2} m n^{2} - {\left({\ell}^{4} n^{2} + 2 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{4} + 4 \, {\ell}^{2} m n + 2 \, {\ell}^{2} m - {\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{2}\right)} {r_H}^{3} - {\left({\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{4} - 2 \, m n^{2} + {\left(n^{2} + 2 \, n + 1\right)} r^{2} - 4 \, m n - 2 \, m\right)} {r_H}\right)}}

The coefficient of the arctan term is

s = Y1.coefficient(arctan(l*r/sqrt(l^2*rH^2 + 1))).simplify_full().factor() s
(32n2rH2+22mn2n2q232rH222m+4n2+p24)2(2r4+r22m)(22rH2+1)2rH2+1(n+1)2\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{{\left(3 \, {\ell}^{2} n^{2} {r_H}^{2} + 2 \, {\ell}^{2} m n^{2} - n^{2} {\mathfrak{q}}^{2} - 3 \, {\ell}^{2} {r_H}^{2} - 2 \, {\ell}^{2} m + 4 \, n^{2} + {\mathfrak{p}}^{2} - 4\right)} {\ell}}{2 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(n + 1\right)}^{2}}

The numerator of this term agrees with Eq. (4.10):

s.numerator().subs({m: m_rH}).factor()
(4n2rH44rH4+42n2rH2n2q242rH2+4n2+p24)\renewcommand{\Bold}[1]{\mathbf{#1}}-{\left({\ell}^{4} n^{2} {r_H}^{4} - {\ell}^{4} {r_H}^{4} + 4 \, {\ell}^{2} n^{2} {r_H}^{2} - n^{2} {\mathfrak{q}}^{2} - 4 \, {\ell}^{2} {r_H}^{2} + 4 \, n^{2} + {\mathfrak{p}}^{2} - 4\right)} {\ell}
bool(s.numerator().subs({m: m_rH}) == l*((1 - n^2)*(l^2*rH^2 + 2)^2 - pf^2 + n^2*qf^2))
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}

The denominator agrees with Eq. (4.10) as well:

s.denominator()
2(2r4+r22m)(22rH2+1)2rH2+1(n+1)2\renewcommand{\Bold}[1]{\mathbf{#1}}2 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(n + 1\right)}^{2}

Let us remove the arctan term from Υ\Upsilon:

Y2 = (Y1 - s*arctan(l*r/sqrt(l^2*rH^2 + 1))).simplify_full() Y2
12(4n24)rrH3+6(2n22)rrH+(n2q2+22m(22m+1)n2+3(2n22)rH2p2+1)log(r+rHrrH)4(2(22mn2(4n2+24n+4)r4+42mn+22m(2n2+22n+2)r2)rH3((2n2+22n+2)r42mn2+(n2+2n+1)r24mn2m)rH)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{12 \, {\left({\ell}^{4} n^{2} - {\ell}^{4}\right)} r {r_H}^{3} + 6 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} r {r_H} + {\left(n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{2} m - {\left(2 \, {\ell}^{2} m + 1\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} {r_H}^{2} - {\mathfrak{p}}^{2} + 1\right)} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)}{4 \, {\left(2 \, {\left(2 \, {\ell}^{2} m n^{2} - {\left({\ell}^{4} n^{2} + 2 \, {\ell}^{4} n + {\ell}^{4}\right)} r^{4} + 4 \, {\ell}^{2} m n + 2 \, {\ell}^{2} m - {\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{2}\right)} {r_H}^{3} - {\left({\left({\ell}^{2} n^{2} + 2 \, {\ell}^{2} n + {\ell}^{2}\right)} r^{4} - 2 \, m n^{2} + {\left(n^{2} + 2 \, n + 1\right)} r^{2} - 4 \, m n - 2 \, m\right)} {r_H}\right)}}
Y2.numerator().simplify_full()
12(4n24)rrH3+6(2n22)rrH+(n2q2+22m(22m+1)n2+3(2n22)rH2p2+1)log(r+rHrrH)\renewcommand{\Bold}[1]{\mathbf{#1}}12 \, {\left({\ell}^{4} n^{2} - {\ell}^{4}\right)} r {r_H}^{3} + 6 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} r {r_H} + {\left(n^{2} {\mathfrak{q}}^{2} + 2 \, {\ell}^{2} m - {\left(2 \, {\ell}^{2} m + 1\right)} n^{2} + 3 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} {r_H}^{2} - {\mathfrak{p}}^{2} + 1\right)} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)
Y2.denominator().factor()
4(2r4+r22m)(22rH2+1)(n+1)2rH\renewcommand{\Bold}[1]{\mathbf{#1}}4 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} {\left(n + 1\right)}^{2} {r_H}

The coefficient of the log term is

s = Y2.coefficient(log((r + rH)/(r - rH))).simplify_full().factor() s
32n2rH222mn2+n2q232rH2+22mn2p2+14(2r4+r22m)(22rH2+1)(n+1)2rH\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{3 \, {\ell}^{2} n^{2} {r_H}^{2} - 2 \, {\ell}^{2} m n^{2} + n^{2} {\mathfrak{q}}^{2} - 3 \, {\ell}^{2} {r_H}^{2} + 2 \, {\ell}^{2} m - n^{2} - {\mathfrak{p}}^{2} + 1}{4 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} {\left(n + 1\right)}^{2} {r_H}}
s.numerator().subs({m: m_rH}).factor()
4n2rH4+4rH4+22n2rH2+n2q222rH2n2p2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-{\ell}^{4} n^{2} {r_H}^{4} + {\ell}^{4} {r_H}^{4} + 2 \, {\ell}^{2} n^{2} {r_H}^{2} + n^{2} {\mathfrak{q}}^{2} - 2 \, {\ell}^{2} {r_H}^{2} - n^{2} - {\mathfrak{p}}^{2} + 1

Check against Eq. (4.10):

bool(s.numerator().subs({m: m_rH}) == (1 - n^2)*(l^2*rH^2 - 1)^2 - pf^2 + n^2*qf^2)
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}
s.denominator()
4(2r4+r22m)(22rH2+1)(n+1)2rH\renewcommand{\Bold}[1]{\mathbf{#1}}4 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} {\left(n + 1\right)}^{2} {r_H}

Given that artanhx=12ln(1+x1x) \mathrm{artanh}\, x = \frac{1}{2} \ln\left( \frac{1 + x}{1 - x} \right) we have ln(x+1x1)=2artanh(1x) \ln \left( \frac{x + 1}{x - 1} \right) = 2\, \mathrm{artanh}\left(\frac{1}{x}\right)

Hence the term in ln(r+rHrrH)\ln\left(\frac{r + r_H}{r - r_H}\right) agrees with the corresponding term in Eq. (4.10).

Finally, the last term in Υ\Upsilon is

Y3 = (Y2 - s*log((r + rH)/(r - rH))).simplify_full() Y3.factor()
32(n1)r2(2r4+r22m)(n+1)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{3 \, {\ell}^{2} {\left(n - 1\right)} r}{2 \, {\left({\ell}^{2} r^{4} + r^{2} - 2 \, m\right)} {\left(n + 1\right)}}

This term agrees with Eq. (4.10), given the simplification 1n2(1+n)2=n1n+1\frac{1 - n^2}{(1 + n)^2} = -\frac{n - 1}{n + 1}.

Conclusion: we have full agreement with Eq. (4.10).

Conjugate momenta

def conjugate_momenta(lagr, qs, var): r""" Compute the conjugate momenta from a given Lagrangian. INPUT: - ``lagr`` -- symbolic expression representing the Lagrangian density - ``qs`` -- either a single symbolic function or a list/tuple of symbolic functions, representing the `q`'s; these functions must appear in ``lagr`` up to at most their first derivatives - ``var`` -- either a single variable, typically `t` (1-dimensional problem) or a list/tuple of symbolic variables; in the latter case the time coordinate must the first one OUTPUT: - list of conjugate momenta; if only one function is involved, the single conjugate momentum is returned instead. """ if not isinstance(qs, (list, tuple)): qs = [qs] if not isinstance(var, (list, tuple)): var = [var] n = len(qs) d = len(var) dqvt = [SR.var('qxxxx{}_t'.format(q)) for q in qs] subs = {diff(qs[i](*var), var[0]): dqvt[i] for i in range(n)} subs_inv = {dqvt[i]: diff(qs[i](*var), var[0]) for i in range(n)} lg = lagr.substitute(subs) ps = [diff(lg, dotq).simplify_full().substitute(subs_inv) for dotq in dqvt] if n == 1: return ps[0] return ps
pis = conjugate_momenta(L_a2, [phi_1, psi_1], r) pis
[(μ021)a22r4rϕ1(r)(μ021)a2r2rϕ1(r)+2(μ021)a2mrϕ1(r),μ02a22n2r4rψ1(r)+μ02a2n2r2rψ1(r)2μ02a2mn2rψ1(r)]\renewcommand{\Bold}[1]{\mathbf{#1}}\left[-{\left({\mu_0}^{2} - 1\right)} a^{2} {\ell}^{2} r^{4} \frac{\partial}{\partial r}\phi_{1}\left(r\right) - {\left({\mu_0}^{2} - 1\right)} a^{2} r^{2} \frac{\partial}{\partial r}\phi_{1}\left(r\right) + 2 \, {\left({\mu_0}^{2} - 1\right)} a^{2} m \frac{\partial}{\partial r}\phi_{1}\left(r\right), {\mu_0}^{2} a^{2} {\ell}^{2} n^{2} r^{4} \frac{\partial}{\partial r}\psi_{1}\left(r\right) + {\mu_0}^{2} a^{2} n^{2} r^{2} \frac{\partial}{\partial r}\psi_{1}\left(r\right) - 2 \, {\mu_0}^{2} a^{2} m n^{2} \frac{\partial}{\partial r}\psi_{1}\left(r\right)\right]

Check of Eq. (4.15):

pi_phi_r = (pis[0]/a).substitute_function(phi_1, phi1_sol).simplify_full() pi_phi_r
(μ021)ap\renewcommand{\Bold}[1]{\mathbf{#1}}-{\left({\mu_0}^{2} - 1\right)} a {\mathfrak{p}}
pi_psi_r = (pis[1]/a).substitute_function(psi_1, psi1_sol).simplify_full() pi_psi_r
μ02an2q\renewcommand{\Bold}[1]{\mathbf{#1}}{\mu_0}^{2} a n^{2} {\mathfrak{q}}

Check of Eq. (4.14):

We start from πθr\pi^r_\theta as given by Eq. (4.13):

pi_theta = (fr4*(a + b)^2*Y_sol(r)).simplify_full() pi_theta
2(3(a23n2sin(Θ2)a23sin(Θ2))rH3(a2n2q2sin(Θ2)+2a23msin(Θ2)a2p2sin(Θ2)+4a2sin(Θ2)2(a23msin(Θ2)+2a2sin(Θ2))n2)rH)arctan(r2rH2+1)2rH2+1(4(2C1a22n2+4C1a22n+2C1a22+3(a24n2sin(Θ2)a24sin(Θ2))r)rH3+2(2C1a2n2+4C1a2n+2C1a2+3(a22n2sin(Θ2)a22sin(Θ2))r)rH+(a2n2q2sin(Θ2)+2a22msin(Θ2)a2p2sin(Θ2)(2a22msin(Θ2)+a2sin(Θ2))n2+3(a22n2sin(Θ2)a22sin(Θ2))rH2+a2sin(Θ2))log(r+rHrrH))4(22rH3+rH)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{2 \, {\left(3 \, {\left(a^{2} {\ell}^{3} n^{2} \sin\left({\Theta_2}\right) - a^{2} {\ell}^{3} \sin\left({\Theta_2}\right)\right)} {r_H}^{3} - {\left(a^{2} {\ell} n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, a^{2} {\ell}^{3} m \sin\left({\Theta_2}\right) - a^{2} {\ell} {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) + 4 \, a^{2} {\ell} \sin\left({\Theta_2}\right) - 2 \, {\left(a^{2} {\ell}^{3} m \sin\left({\Theta_2}\right) + 2 \, a^{2} {\ell} \sin\left({\Theta_2}\right)\right)} n^{2}\right)} {r_H}\right)} \arctan\left(\frac{{\ell} r}{\sqrt{{\ell}^{2} {r_H}^{2} + 1}}\right) - \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(4 \, {\left(2 \, C_{1} a^{2} {\ell}^{2} n^{2} + 4 \, C_{1} a^{2} {\ell}^{2} n + 2 \, C_{1} a^{2} {\ell}^{2} + 3 \, {\left(a^{2} {\ell}^{4} n^{2} \sin\left({\Theta_2}\right) - a^{2} {\ell}^{4} \sin\left({\Theta_2}\right)\right)} r\right)} {r_H}^{3} + 2 \, {\left(2 \, C_{1} a^{2} n^{2} + 4 \, C_{1} a^{2} n + 2 \, C_{1} a^{2} + 3 \, {\left(a^{2} {\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - a^{2} {\ell}^{2} \sin\left({\Theta_2}\right)\right)} r\right)} {r_H} + {\left(a^{2} n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, a^{2} {\ell}^{2} m \sin\left({\Theta_2}\right) - a^{2} {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) - {\left(2 \, a^{2} {\ell}^{2} m \sin\left({\Theta_2}\right) + a^{2} \sin\left({\Theta_2}\right)\right)} n^{2} + 3 \, {\left(a^{2} {\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - a^{2} {\ell}^{2} \sin\left({\Theta_2}\right)\right)} {r_H}^{2} + a^{2} \sin\left({\Theta_2}\right)\right)} \log\left(\frac{r + {r_H}}{r - {r_H}}\right)\right)}}{4 \, {\left(2 \, {\ell}^{2} {r_H}^{3} + {r_H}\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

Let us perform an expansion in 1/r1/r for r+r\rightarrow +\infty:

assume(l > 0) u = var('u') assume(u > 0) s = pi_theta.subs({r: 1/u}).simplify_log() s = s.taylor(u, 0, 2) s = s.subs({u: 1/r}) s
32(a22n2sin(Θ2)a22sin(Θ2))r+3(a2n2sin(Θ2)a2sin(Θ2))2r+πa2n2q2sin(Θ2)+2πa23msin(Θ2)πa2p2sin(Θ2)+4πa2sin(Θ2)2(πa23msin(Θ2)+2πa2sin(Θ2))n23(πa23n2sin(Θ2)πa23sin(Θ2))rH2+4(C1a2n2+2C1a2n+C1a2+2(C1a22n2+2C1a22n+C1a22)rH2)2rH2+14(22rH2+1)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{3}{2} \, {\left(a^{2} {\ell}^{2} n^{2} \sin\left({\Theta_2}\right) - a^{2} {\ell}^{2} \sin\left({\Theta_2}\right)\right)} r + \frac{3 \, {\left(a^{2} n^{2} \sin\left({\Theta_2}\right) - a^{2} \sin\left({\Theta_2}\right)\right)}}{2 \, r} + \frac{\pi a^{2} {\ell} n^{2} {\mathfrak{q}}^{2} \sin\left({\Theta_2}\right) + 2 \, \pi a^{2} {\ell}^{3} m \sin\left({\Theta_2}\right) - \pi a^{2} {\ell} {\mathfrak{p}}^{2} \sin\left({\Theta_2}\right) + 4 \, \pi a^{2} {\ell} \sin\left({\Theta_2}\right) - 2 \, {\left(\pi a^{2} {\ell}^{3} m \sin\left({\Theta_2}\right) + 2 \, \pi a^{2} {\ell} \sin\left({\Theta_2}\right)\right)} n^{2} - 3 \, {\left(\pi a^{2} {\ell}^{3} n^{2} \sin\left({\Theta_2}\right) - \pi a^{2} {\ell}^{3} \sin\left({\Theta_2}\right)\right)} {r_H}^{2} + 4 \, {\left(C_{1} a^{2} n^{2} + 2 \, C_{1} a^{2} n + C_{1} a^{2} + 2 \, {\left(C_{1} a^{2} {\ell}^{2} n^{2} + 2 \, C_{1} a^{2} {\ell}^{2} n + C_{1} a^{2} {\ell}^{2}\right)} {r_H}^{2}\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}{4 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

We consider πθr(a2/2)sin(2Θ0)\frac{\pi^r_\theta}{(a^2/2)\sin(2\Theta_0)}:

s1 = (s/(a^2/2*sin(Th2))).expand() s1
3π3n2rH22(22rH2+1)2rH2+1π3mn2(22rH2+1)2rH2+1+32n2r+4C12n2rH2(22rH2+1)sin(Θ2)+πn2q22(22rH2+1)2rH2+1+3π3rH22(22rH2+1)2rH2+1+8C12nrH2(22rH2+1)sin(Θ2)+π3m(22rH2+1)2rH2+132r+4C12rH2(22rH2+1)sin(Θ2)2πn2(22rH2+1)2rH2+1πp22(22rH2+1)2rH2+1+3n2r+2C1n2(22rH2+1)sin(Θ2)+2π(22rH2+1)2rH2+1+4C1n(22rH2+1)sin(Θ2)3r+2C1(22rH2+1)sin(Θ2)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{3 \, \pi {\ell}^{3} n^{2} {r_H}^{2}}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} - \frac{\pi {\ell}^{3} m n^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} + 3 \, {\ell}^{2} n^{2} r + \frac{4 \, C_{1} {\ell}^{2} n^{2} {r_H}^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{\pi {\ell} n^{2} {\mathfrak{q}}^{2}}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} + \frac{3 \, \pi {\ell}^{3} {r_H}^{2}}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} + \frac{8 \, C_{1} {\ell}^{2} n {r_H}^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{\pi {\ell}^{3} m}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} - 3 \, {\ell}^{2} r + \frac{4 \, C_{1} {\ell}^{2} {r_H}^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} - \frac{2 \, \pi {\ell} n^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} - \frac{\pi {\ell} {\mathfrak{p}}^{2}}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} + \frac{3 \, n^{2}}{r} + \frac{2 \, C_{1} n^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{2 \, \pi {\ell}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}} + \frac{4 \, C_{1} n}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} - \frac{3}{r} + \frac{2 \, C_{1}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)}

The term in factor of C1C_1 is

s1.coefficient(C_1)
42n2rH2(22rH2+1)sin(Θ2)+82nrH2(22rH2+1)sin(Θ2)+42rH2(22rH2+1)sin(Θ2)+2n2(22rH2+1)sin(Θ2)+4n(22rH2+1)sin(Θ2)+2(22rH2+1)sin(Θ2)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{4 \, {\ell}^{2} n^{2} {r_H}^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{8 \, {\ell}^{2} n {r_H}^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{4 \, {\ell}^{2} {r_H}^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{2 \, n^{2}}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{4 \, n}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)} + \frac{2}{{\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sin\left({\Theta_2}\right)}

It is a constant term, of the type C~1/sin(2Θ0)\tilde{C}_1/\sin(2\Theta_0), in agreement with Eq. (4.14). We remove it from the main term:

s2 = (s1 - s1.coefficient(C_1)*C_1).simplify_full() s2
12(4n24+(6n26)r2)rH4+6(2n22)r2+18(2n2+(4n24)r22)rH2+6n22rH2+1(3(π3n2π3)rrH2(πn2q2+2π3mπp22(π3m+2π)n2+4π)r)62(24rrH4+32rrH2+r)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{12 \, {\left({\ell}^{4} n^{2} - {\ell}^{4} + {\left({\ell}^{6} n^{2} - {\ell}^{6}\right)} r^{2}\right)} {r_H}^{4} + 6 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} r^{2} + 18 \, {\left({\ell}^{2} n^{2} + {\left({\ell}^{4} n^{2} - {\ell}^{4}\right)} r^{2} - {\ell}^{2}\right)} {r_H}^{2} + 6 \, n^{2} - \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(3 \, {\left(\pi {\ell}^{3} n^{2} - \pi {\ell}^{3}\right)} r {r_H}^{2} - {\left(\pi {\ell} n^{2} {\mathfrak{q}}^{2} + 2 \, \pi {\ell}^{3} m - \pi {\ell} {\mathfrak{p}}^{2} - 2 \, {\left(\pi {\ell}^{3} m + 2 \, \pi {\ell}\right)} n^{2} + 4 \, \pi {\ell}\right)} r\right)} - 6}{2 \, {\left(2 \, {\ell}^{4} r {r_H}^{4} + 3 \, {\ell}^{2} r {r_H}^{2} + r\right)}}
s2.numerator().simplify_full()
12(4n24+(6n26)r2)rH4+6(2n22)r2+18(2n2+(4n24)r22)rH2+6n22rH2+1(3(π3n2π3)rrH2(πn2q2+2π3mπp22(π3m+2π)n2+4π)r)6\renewcommand{\Bold}[1]{\mathbf{#1}}12 \, {\left({\ell}^{4} n^{2} - {\ell}^{4} + {\left({\ell}^{6} n^{2} - {\ell}^{6}\right)} r^{2}\right)} {r_H}^{4} + 6 \, {\left({\ell}^{2} n^{2} - {\ell}^{2}\right)} r^{2} + 18 \, {\left({\ell}^{2} n^{2} + {\left({\ell}^{4} n^{2} - {\ell}^{4}\right)} r^{2} - {\ell}^{2}\right)} {r_H}^{2} + 6 \, n^{2} - \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\left(3 \, {\left(\pi {\ell}^{3} n^{2} - \pi {\ell}^{3}\right)} r {r_H}^{2} - {\left(\pi {\ell} n^{2} {\mathfrak{q}}^{2} + 2 \, \pi {\ell}^{3} m - \pi {\ell} {\mathfrak{p}}^{2} - 2 \, {\left(\pi {\ell}^{3} m + 2 \, \pi {\ell}\right)} n^{2} + 4 \, \pi {\ell}\right)} r\right)} - 6
s2.denominator().factor()
2(22rH2+1)(2rH2+1)r\renewcommand{\Bold}[1]{\mathbf{#1}}2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} {\left({\ell}^{2} {r_H}^{2} + 1\right)} r

Let divide both the numerator and denominator by rr:

s2n = (s2.numerator()/r).expand() s2n
126n2rrH4126rrH4+184n2rrH2+124n2rH4r3π2rH2+13n2rH22π2rH2+13mn2184rrH2124rH4r+π2rH2+1n2q2+3π2rH2+13rH2+2π2rH2+13m+62n2r+182n2rH2r4π2rH2+1n2π2rH2+1p262r182rH2r+4π2rH2+1+6n2r6r\renewcommand{\Bold}[1]{\mathbf{#1}}12 \, {\ell}^{6} n^{2} r {r_H}^{4} - 12 \, {\ell}^{6} r {r_H}^{4} + 18 \, {\ell}^{4} n^{2} r {r_H}^{2} + \frac{12 \, {\ell}^{4} n^{2} {r_H}^{4}}{r} - 3 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} n^{2} {r_H}^{2} - 2 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} m n^{2} - 18 \, {\ell}^{4} r {r_H}^{2} - \frac{12 \, {\ell}^{4} {r_H}^{4}}{r} + \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} n^{2} {\mathfrak{q}}^{2} + 3 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} {r_H}^{2} + 2 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} m + 6 \, {\ell}^{2} n^{2} r + \frac{18 \, {\ell}^{2} n^{2} {r_H}^{2}}{r} - 4 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} n^{2} - \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} {\mathfrak{p}}^{2} - 6 \, {\ell}^{2} r - \frac{18 \, {\ell}^{2} {r_H}^{2}}{r} + 4 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} + \frac{6 \, n^{2}}{r} - \frac{6}{r}
s2d = (s2.denominator()/r).factor() s2d
2(22rH2+1)(2rH2+1)\renewcommand{\Bold}[1]{\mathbf{#1}}2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} {\left({\ell}^{2} {r_H}^{2} + 1\right)}

The coefficient of the term in rr is

s = s2n.coefficient(r).factor() s/s2d
32(n+1)(n1)\renewcommand{\Bold}[1]{\mathbf{#1}}3 \, {\ell}^{2} {\left(n + 1\right)} {\left(n - 1\right)}

This is in agreement with Eq. (4.14).

We remove it:

s3n = (s2n - s*r).simplify_full().expand() s3n
124n2rH4r3π2rH2+13n2rH22π2rH2+13mn2124rH4r+π2rH2+1n2q2+3π2rH2+13rH2+2π2rH2+13m+182n2rH2r4π2rH2+1n2π2rH2+1p2182rH2r+4π2rH2+1+6n2r6r\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{12 \, {\ell}^{4} n^{2} {r_H}^{4}}{r} - 3 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} n^{2} {r_H}^{2} - 2 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} m n^{2} - \frac{12 \, {\ell}^{4} {r_H}^{4}}{r} + \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} n^{2} {\mathfrak{q}}^{2} + 3 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} {r_H}^{2} + 2 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} m + \frac{18 \, {\ell}^{2} n^{2} {r_H}^{2}}{r} - 4 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} n^{2} - \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} {\mathfrak{p}}^{2} - \frac{18 \, {\ell}^{2} {r_H}^{2}}{r} + 4 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} + \frac{6 \, n^{2}}{r} - \frac{6}{r}

The coefficient of the term in 1/r1/r is

s = s3n.coefficient(r^(-1)).factor() s/s2d
3(n+1)(n1)\renewcommand{\Bold}[1]{\mathbf{#1}}3 \, {\left(n + 1\right)} {\left(n - 1\right)}

This is in agreement with Eq. (4.14).

Finally the remaining term is

s4n = (s3n - s/r).simplify_full().expand() s4n
3π2rH2+13n2rH22π2rH2+13mn2+π2rH2+1n2q2+3π2rH2+13rH2+2π2rH2+13m4π2rH2+1n2π2rH2+1p2+4π2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-3 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} n^{2} {r_H}^{2} - 2 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} m n^{2} + \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} n^{2} {\mathfrak{q}}^{2} + 3 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} {r_H}^{2} + 2 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}^{3} m - 4 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} n^{2} - \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell} {\mathfrak{p}}^{2} + 4 \, \pi \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}
s4n.factor()
π(32n2rH2+22mn2n2q232rH222m+4n2+p24)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-\pi {\left(3 \, {\ell}^{2} n^{2} {r_H}^{2} + 2 \, {\ell}^{2} m n^{2} - n^{2} {\mathfrak{q}}^{2} - 3 \, {\ell}^{2} {r_H}^{2} - 2 \, {\ell}^{2} m + 4 \, n^{2} + {\mathfrak{p}}^{2} - 4\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1} {\ell}
s = (s4n.factor()/s2d).subs({m: m_rH}).factor() s
π(4n2rH44rH4+42n2rH2n2q242rH2+4n2+p24)2(22rH2+1)2rH2+1\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{\pi {\left({\ell}^{4} n^{2} {r_H}^{4} - {\ell}^{4} {r_H}^{4} + 4 \, {\ell}^{2} n^{2} {r_H}^{2} - n^{2} {\mathfrak{q}}^{2} - 4 \, {\ell}^{2} {r_H}^{2} + 4 \, n^{2} + {\mathfrak{p}}^{2} - 4\right)} {\ell}}{2 \, {\left(2 \, {\ell}^{2} {r_H}^{2} + 1\right)} \sqrt{{\ell}^{2} {r_H}^{2} + 1}}

The denominator clearly agrees with Eq. (4.14); the numerator agrees as well:

bool(s.numerator()/(pi*l) == (1 - n^2)*(l^2*rH^2 + 2)^2 - pf^2 + n^2*qf^2)
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}

Conclusion: we have full agreement with Eq. (4.14).