# Principal null directions in Kerr spacetime¶

This notebook demonstrates a few capabilities of SageMath in computations regarding Kerr spacetime. More precisely, it focuses on the principal null directions. This notebook makes use of SageMath tools developed through the SageManifolds project.

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

In [1]:
version()

Out[1]:
'SageMath version 8.8, Release Date: 2019-06-26'

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

In [2]:
%display latex


To speed up the computations, we ask for running them in parallel on 8 cores:

In [3]:
Parallelism().set(nproc=8)


## Spacetime manifold¶

We declare the Kerr spacetime (or more precisely the Boyer-Lindquist domain of Kerr spacetime) as a 4-dimensional Lorentzian manifold:

In [4]:
M = Manifold(4, 'M', latex_name=r'\mathcal{M}', structure='Lorentzian')
print(M)

4-dimensional Lorentzian manifold M


Let us declare the Boyer-Lindquist coordinates 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:

In [5]:
BL.<t,r,th,ph> = M.chart(r't r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi')
print(BL) ; BL

Chart (M, (t, r, th, ph))

Out[5]:
In [6]:
BL[0], BL[1]

Out[6]:

## Metric tensor

The 2 parameters $m$ and $a$ of the Kerr spacetime are declared as symbolic variables:

In [7]:
var('m, a', domain='real')

Out[7]:

We get the (yet undefined) spacetime metric by

In [8]:
g = M.metric()


The metric is set by its components in the coordinate frame associated with Boyer-Lindquist coordinates, which is the current manifold's default frame:

In [9]:
rho2 = r^2 + (a*cos(th))^2
Delta = r^2 -2*m*r + a^2
g[0,0] = -(1-2*m*r/rho2)
g[0,3] = -2*a*m*r*sin(th)^2/rho2
g[1,1], g[2,2] = rho2/Delta, rho2
g[3,3] = (r^2+a^2+2*m*r*(a*sin(th))^2/rho2)*sin(th)^2
g.display()

Out[9]:

A matrix view of the components with respect to the manifold's default vector frame:

In [10]:
g[:]

Out[10]:

The list of the non-vanishing components:

In [11]:
g.display_comp()

Out[11]:

## Levi-Civita Connection

The Levi-Civita connection $\nabla$ associated with $g$:

In [12]:
nabla = g.connection() ; print(nabla)

Levi-Civita connection nabla_g associated with the Lorentzian metric g on the 4-dimensional Lorentzian manifold M


Let us verify that the covariant derivative of $g$ with respect to $\nabla$ vanishes identically:

In [13]:
nabla(g).display()

Out[13]:

## Killing vectors

The default vector frame on the spacetime manifold is the coordinate basis associated with Boyer-Lindquist coordinates:

In [14]:
M.default_frame() is BL.frame()

Out[14]:
In [15]:
BL.frame()

Out[15]:

Let us consider the first vector field of this frame:

In [16]:
xi = BL.frame()[0] ; xi

Out[16]:
In [17]:
print(xi)

Vector field d/dt on the 4-dimensional Lorentzian manifold M


The 1-form associated to it by metric duality is

In [18]:
xi_form = xi.down(g) ; xi_form.display()

Out[18]:

Its covariant derivative is

In [19]:
nab_xi = nabla(xi_form) ; print(nab_xi) ; nab_xi.display()

Tensor field of type (0,2) on the 4-dimensional Lorentzian manifold M

Out[19]:

Let us check that the Killing equation is satisfied:

In [20]:
nab_xi.symmetrize() == 0

Out[20]:

Similarly, let us check that $\frac{\partial}{\partial\phi}$ is a Killing vector:

In [21]:
chi = BL.frame()[3] ; chi

Out[21]:
In [22]:
nabla(chi.down(g)).symmetrize() == 0

Out[22]:

## Principal null vectors¶

Let us consider the following vector fields $k$ and $\ell$, defined from their components with respect to Boyer-Lindquist coordinates:

In [23]:
k = M.vector_field(name='k')
k[:] = [(r^2+a^2)/(2*rho2), -Delta/(2*rho2), 0, a/(2*rho2)]
k.display()

Out[23]:
In [24]:
el = M.vector_field(name='el', latex_name=r'\ell')
el[:] = [(r^2+a^2)/Delta, 1, 0, a/Delta]
el.display()

Out[24]:

Let us check that $k$ and $\ell$ are null vectors:

In [25]:
g(k,k).expr()

Out[25]:
In [26]:
g(el,el).expr()

Out[26]:

Their scalar product is $-1$:

In [27]:
g(k,el).expr()

Out[27]:

Note that the scalar product (with respect to metric $g$) can also be computed by means of the method dot:

In [28]:
k.dot(el).expr()

Out[28]:

Let us evaluate the "acceleration" of $k$, i.e. $\nabla_k k$:

In [29]:
acc_k = nabla(k).contract(k)
acc_k.display()

Out[29]:

We check that $k$ is a pregeodesic vector, i.e. that $\nabla_k k = \kappa_k k$ for some scalar field $\kappa_k$:

In [30]:
for i in [0,1,3]:
show(acc_k[i] / k[i])

In [31]:
kappa_k = acc_k[[0]] / k[[0]]
kappa_k.display()

Out[31]:
In [32]:
acc_k == kappa_k * k

Out[32]:

Similarly let us evaluate the "acceleration" of $\ell$:

In [33]:
acc_l = nabla(el).contract(el)
acc_l.display()

Out[33]:

Hence $\ell$ obeys $\nabla_\ell\, \ell=0$, i.e. it is a geodesic vector.

### Check that $k$ and $\ell$ are doubly degenerate principal null vectors¶

$k$ defines a repeated null direction (one says that $k$ is a doubly degenerate principal null vector) iff it obeys the following identity: $$C^a_{\ \, mn[b} k_{c]} k^m k^n = 0 \qquad \mbox{(1)},$$ where $C$ is the Weyl tensor of the metric $g$.

The Weyl tensor is computed via

In [34]:
C = g.weyl()
print(C)

Tensor field C(g) of type (1,3) on the 4-dimensional Lorentzian manifold M


For instance, the component $C^0_{\ \, 010} = C^t_{\ \, trt}$ is

In [35]:
C[0, 1, 0, 1]

Out[35]:

We can check the identity (1) in a single line as follows:

In [36]:
(C.contract(1, 2, k*k, 0, 1) * k.down(g)).antisymmetrize(1, 2).display()

Out[36]:

Alternatively, we may check (1) by using index notations. First, we introduce the tensor $A^a_{\ \, b} = C^a_{\ \, mnb} k^m k^n$:

In [37]:
A = C['^a_{mnb}'] * (k*k)['^{mn}']
print(A)

Tensor field of type (1,1) on the 4-dimensional Lorentzian manifold M


and the 1-form $k_a = g_{am} k^m$:

In [38]:
kf = g['_{am}'] * k['^m']
print(kf)

1-form on the 4-dimensional Lorentzian manifold M


Then (1) is equivalent to $A^a_{\ \, [b} k_{c]} = 0$, which is readily checked as

In [39]:
(A*kf)['^a_{[bc]}'].display()

Out[39]:

Similarly, let us check that $\ell$ obeys (1) as well:

In [40]:
(C.contract(1, 2, el*el, 0, 1) * el.down(g)).antisymmetrize(1, 2).display()

Out[40]: