T'b'(x') = Ta(x') d xa(x') / d x'b'As always, repeated indices are summed over, and the primed indices indicate the new coordinates. Note that the valence of the primed index is inverted: an upper index in the denominator is equivalent to a lower index in the numerator.
This transformation consists of two operations. The first involves a substitution of the original coordinates as functions of the new ones in the original tensor. The resulting tensor is then multiplied by the matrix whose elements are the derivatives of the old coordinates as functions of the new. For tensors with multiple covariant indices, the matrix multiplication is performed on each original index: the transformation is a multilinear operation, the tensor transforming linearly in each index.
As an example, we present the Mathematica code to perform a coordinate transformation
on a tensor with two covariant indices. Here, "xu" is the list of original coordinates,
and "xfu" is the list of those coordinates as functions of the new, primed coordinates
"xpu" (as mentioned in the Introduction, the "u"s indicate that these indices are
contravariant indices, and the "l"s denote covariant indices):
dxofxpdxplu = Table[ Table[
{j, dim}], {i, dim}];
{l, dim}], {k, dim}]],
The first statement performs the substitution, using the ReplaceAll (/.) and Rule (->) operators.
The use of Thread causes corresponding
elements of the lists to be substituted. For instance, xu[[1]] will be replaced by
xfu[[1]], xu[[2]] will be replaced by xfu[[2]], etc. The second statement creates the
matrix of derivatives, and the final statement performs the matrix multiplication
(and serves as the model for all index contraction in the sequel).
Note the use of the simpler function defined in the Introduction.
It will be useful to create a series of functions which perform the basic tensor
operations discussed in this section. In general, these operations can be performed on
tensors with any number of indices. As an example of how to program such a function,
we present "coordxform", which performs a coordinate transformation on a tensor with from
1 to 5 covariant indices:
coordxform[ T_List, x_List, xf_List, xp_List] := Block[ {dim, dimp, rank, Tp, Tofxp, dxofxpdxp}, (
If[ rank > 5,(Print[ "Error - tensor has too many indices for this function"]; Return[];),];
dim = Length[x];
If[ Length[T] == dim,,(Print[ "Error - tensor is incorrect dimension"]; Return[];)];
If[ Length[xf] == dim,,(Print[ "Error - function list x(xp) is incorrect dimension"]; Return[];)];
dimp = Length[xp];
Tofxp = T /. Thread[ x -> xf];
dxofxpdxp = Table[ Table[
{jj, dim}], {ii, dimp}];
{jj, dim}]],
{ll, dim}], {kk, dim}]],
{nn, dim}], {mm, dim}], {ll, dim}]],
{qq, dim}], {pp, dim}], {nn, dim}], {mm, dim}]],
{ww, dim}], {vv, dim}], {qq, dim}], {pp, dim}], {nn, dim}]],
A series of If statements are used to ensure that coordxform is called with
arguments that are appropriate to its definition and consistent with each other. Note
that Mathematica stores a matrix as a list of lists; its Length is therefore the
number of lists at the top level. These tests also contribute to making the code more
resilient in the face of accidental misuse.
The substitution operation and the creation of the derivative matrix are essentially identical
to those in the previous example. The remaining statements perform the number of
matrix multiplications appropriate for the rank of the tensor being transformed.
Note the index variables: the author has adopted the convention of using two letter index
variables inside functions so that it is less likely that a symbol occurring in an argument
will have the same name as an index variable.
Here is a Mathematica code example using the function coordxform define above. In this example,
we transform the Minkowski metric in spherical coordinates to eliminate the trigonometric functions
using the transformation
Manipulation of index valences are of course accomplished as matrix multiplications by the metric or its
inverse:
Tb = gb a Ta
The Christoffel Connection is defined as
Do[ Do[ Do[
D[ gll[[jj, ll]], xu[[kk]]] - D[ gll[[jj, kk]], xu[[ll]]]),
{ll, dim}] / 2],
{kk, jj + 1, dim}], {jj, dim - 1}], {ii, dim}];
The Christoffel Connection is of course used in the computation of the Riemann Tensor below, but it is
also used in the geodesic equation:
The Riemann Curvature Tensor is defined as
Ra b c d = - Ra b d c
Ra b c d = Rc d a b
R[a b c] d = 0
The following Mathematica code to compute the Riemann Curvature components
makes use of the first three symmetries:
Do[ Do[ Do[ Do[
{nn, dim}]) gll[[ll,mm]],
{ll, kk + 1, dim}], {kk, dim - 1}], {jj, ii + 1, dim}], {ii, dim - 1}];
{ll, kk + 1, dim}], {kk, dim - 1}], {jj, ii + 1, dim}], {ii, dim - 1}];
{ll, kk + 1, dim}], {kk, dim - 1}], {jj, ii + 1, dim}], {ii, dim - 1}];
We note that if Ra b a b is zero, then the a-b hypersurface is flat, but none of the
manifolds we will be examining are that trivial.
We will also be interested in the covariant derivative of the Riemann Curvature:
= Rb c e f ;a,
The following decomposition of the Riemann Curvature [Wald]
is sometimes useful:
Ga b c d = (ga d gb c - ga c gb d) R / (D - 1)(D - 2)
We have chosen to use a metric of positive signature. If one wished to work with a negative signature instead:
so Einstein's Equations become
We can use the following Mathematica code to see if any given metric is a solution:
{jj, dim}], {ii, dim}] === a * Tll]
We define the epsilon tensor as:
The next section discusses curvature invariants.
©2005, Kenneth R. Koehler. All Rights Reserved. This document may be freely reproduced provided that this copyright notice is included.
Please send comments or suggestions to the author.
Tofxpll = Tll /. Thread[ xu -> xfu];
"Tofxpll" might be read "T of x' with two lower indices", while "dxofxpdxplu" might be read
"the derivative of x as a function of x' with respect to x', whose first index is
lower and whose second index is upper". The use of meaningful variable names is always
encouraged as a way of making your code at least somewhat self-documenting.
D[ xfu[[j]], xpu[[i]]],
Tpll = Table[ Table[ simpler[
Sum[ Sum[
Tofxpll[[k, l]] dxofxpdxplu[[i, k]] dxofxpdxplu[[j, l]],
{j, dim}], {i, dim}];
coordxform::usage = "arguments are Tl...l(x), xu, xu(x'), x'u";
The first line documents the usage of coordxform, whose arguments are the original covariant tensor,
a list of the original coordinates, a list consisting of those coordinates as functions of the
new and a list of the new coordinates. Note that the dimension of the new coordinate system
need not be the same as that of the old. The function
definition includes patterns in the argument list which restrict the use of coordxform
to arguments which are lists. This is an effective control on the invocation of the function
which can help to prevent errors. Block is used to force the scope of the variables in
its first argument to be local to coordxform: another technique which can help make the
code less brittle.
rank = TensorRank[T];
D[ xf[[jj]], xp[[ii]]],
If[ rank == 1,
Tp = Table[ simpler[ Sum[
If[ rank == 2,
dxofxpdxp[[ii, jj]] Tofxp[[jj]],
{ii, dimp}],];
Tp = Table[ Table[ simpler[ Sum[ Sum[
If[ rank == 3,
dxofxpdxp[[ii, kk]] dxofxpdxp[[jj, ll]] Tofxp[[kk, ll]],
{jj, dimp}], {ii, dimp}],];
Tp = Table[ Table[ Table[ simpler[ Sum[ Sum[ Sum[
If[ rank == 4,
dxofxpdxp[[ii, ll]] dxofxpdxp[[jj, mm]] dxofxpdxp[[kk, nn]] Tofxp[[ll, mm, nn]],
{kk, dimp}], {jj, dimp}], {ii, dimp}],];
Tp = Table[ Table[ Table[ Table[ simpler[ Sum[ Sum[ Sum[ Sum[
If[ rank == 5,
dxofxpdxp[[ii, mm]] dxofxpdxp[[jj, nn]] dxofxpdxp[[kk, pp]] dxofxpdxp[[ll, qq]] Tofxp[[mm, nn, pp, qq]],
{ll, dimp}], {kk, dimp}], {jj, dimp}], {ii, dimp}],];
Tp = Table[ Table[ Table[ Table[ Table[ simpler[ Sum[ Sum[ Sum[ Sum[ Sum[
Return[Tp];)];
dxofxpdxp[[ii, nn]] dxofxpdxp[[jj, pp]] dxofxpdxp[[kk, qq]] dxofxpdxp[[ll, vv]] dxofxpdxp[[mm, ww]]
Tofxp[[nn, pp, qq, vv, ww]],
{mm, dimp}], {ll, dimp}], {kk, dimp}], {jj, dimp}], {ii, dimp}],];
c = cos (theta).
The third argument to coordxform is derived from the second argument and the inverse transformation.
TableForm is used to format the metrics conveniently:
?coordxform
This transformation is used extensively
throughout the text in order to improve Mathematica execution times.
arguments are Tl...l(x), xu, xu(x'), x'u
TableForm[ gll = {{1, 0, 0, 0}, {0, r^2, 0, 0}, {0, 0, r^2 Sin[theta]^2, 0}, {0, 0, 0, -1}}]
TableForm[ coordxform[ gll, {r, theta, phi, t}, {r, ArcCos[c], phi, t}, {r, c, phi, t}]]
1 0 0 0 0 r2 0 0 0 0 r2 Sin[theta]2 0 0 0 0 -1
1 0 0 0 0 - r2 / (-1+c)(1+c) 0 0 0 0 - (-1+c)(1+c) r2 0 0 0 0 -1
Tb = gb a Ta
The student is encouraged to write functions "raise" and "lower" analogous to coordxform. They will
each have two arguments: the original tensor and the index number affected, and can
be tested by performing raise on either index of the metric, and lower on either index of the
inverse metric. In either case, the result should be the Identity Matrix.
Γab c = ga e (db gc e +
dc gb e - de gb c) / 2,
and its lower indices are symmetric:
Γac b = Γab c
We use this symmetry to speed up the Mathematica code to compute the connection components:
Gull = Table[ Table[ Table[ 0, {kk, dim}], {jj, dim}], {ii, dim}];
Note that we must first create a zero tensor with the correct number of indices; since the remainder
of the code changes individual components, Mathematica requires that they exist before being modified.
This code assumes that the metric associated with this connection is stored in "gll", that the inverse metric
is "guu" and that xu is again the list of coordinates. We will assume that tensor names, once
adopted, are used consistently in all code examples unless otherwise stated.
Gull[[ii, jj, kk]] = simpler[ Sum[
Do[ Do[ Do[
guu[[ii, ll]] (D[ gll[[kk, ll]], xu[[jj]]] +
{kk, jj, dim}], {jj, dim}], {ii, dim}];
Gull[[ii, kk, jj]] = Gull[[ii, jj, kk]],
d2 xa / dτ2 +
Γab c (d xb / dτ)
(d xc / dτ) = 0
where τ is the affine parameter. Obviously, the nonzero components control
the geodesic deviation. Therefore if for some b,
Γab b is only nonzero
for a = b, then all curves tangent to xb are geodesics. In physical terms, their acceleration is
parallel to their velocity. We will call such directions geodesic directions. If b and c are members
of a subset
of the coordinate directions such that Γab c is only nonzero
when a is also a member of this subset, then these directions can be said to define a
geodesic hypersurface: any geodesic tangent to the hypersurface is wholly contained within it.
Ra b ce = db Γea c -
da Γeb c +
Γfa c Γef b -
Γfb c Γef a
and has the following symmetries:
Ra b c d = - Rb a c d
If we consider the indices in pairs, in D dimensions each antisymmetric pair has
NB: The definition of the Riemann Curvature Tensor is independent of metric signature.
However, lowering the final index will induce a relative negative sign if using negative
signature metrics. This effect occurs whenever the metric is involved (ie., the cosmological constant term in
Einstein's Equations, the pressure term in
perfect fluid Stress-Energy Tensors, etc.).
D (D - 1) / 2
independent components. The pairs are symmetric, so that the first three symmetries
imply that there are
(D (D - 1) / 2) (D (D - 1) / 2 + 1) / 2
independent components. But the last symmetry means that we need to subtract
D ! / (4 ! (D - 4) !)
components, leaving us with
D2 (D2 - 1) / 12
independent components. This implies that all D = 1 manifolds have zero curvature.
Rllll = Table[ Table[ Table[ Table[ 0, {ll, dim}], {kk, dim}], {jj, dim}], {ii, dim}];
computing the upper half-matrix for each pair of indicies and use the symmetries
to fill in the remaining components. The final symmetry is more trouble to program than it
is worth.
Rllll[[ii, jj, kk, mm]] = simpler[ Sum[
Do[ Do[ Do[ Do[
(D[ Gull[[ll, ii, kk]], xu[[jj]]] - D[ Gull[[ll, jj, kk]], xu[[ii]]] + Sum[
{mm, kk + 1, dim}], {kk, dim - 1}], {jj, ii + 1, dim}], {ii, dim - 1}];
Gull[[nn, ii, kk]] Gull[[ll, nn, jj]] - Gull[[nn, jj, kk]] Gull[[ll, nn, ii]],
{ll, dim}]],
Rllll[[jj, ii, kk, ll]] = - Rllll[[ii, jj, kk, ll]],
Do[ Do[ Do[ Do[
Rllll[[ii, jj, ll, kk]] = - Rllll[[ii, jj, kk, ll]],
Do[ Do[ Do[ Do[
Rllll[[jj, ii, ll, kk]] = Rllll[[ii, jj, kk, ll]],
Da Rb c e f = da Rb c e f -
Γha b Rh c e f -
the Ricci Tensor:
Γha c Rb h e f -
Γha e Rb c h f -
Γha f Rb c e h
Ra b = gc d Ra c b d
and the scalar curvature:
R = ga b Ra b
These tensors will be used repeatedly in the analyses below, and it makes sense for the
programmer to combine all of them into a single function so that for any metric, one function
call will compute them all. The inverse metric and Christoffel Symbols should of course be
included in the beginning of that function. It will also be useful for the student to write
a function analogous to coordxform which will return the covariant derivative of a tensor
with an arbitary number of covariant indices.
Ra b c d = Ca b c d + Ea b c d + Ga b c d
If we compute E (the Einstein Curvature Tensor) and G analogously to the Riemann Tensor,
the following Mathematica code will compute the Weyl Tensor in a single statement:
Ea b c d = (ga c Rb d - gb c Ra d -
ga d Rb c + gb d Ra c) / (D - 2)
Cllll = Rllll - Ellll - Gllll;
We will of course be interested in Einstein's Equations:
Ra b - R ga b / 2 + Λ ga b =
α Ta b
where the left hand side is called the Einstein Tensor, and the
stress-energy tensor T has dimensions of energy density.
α = 2 Area (SD - 2)
Provided that ga b is dimensionless, α can be multiplied by G / c4 to normalize the dimensions of all terms to
1/length2. However, if ga b is not dimensionless, as for instance in spherically symmetric metrics,
it is easier to work with Einstein's Equations with one index raised. Then gab = Iab and all is well
so long as Tab has dimensions of energy density.
Taking the trace of both sides with Λ = 0, we have
Ra b - R ga b / 2 - Λ ga b =
α Ta b
(1 - D / 2) R = α T
This fact implies that for vacuum solutions (where Ta b
is zero), Ra b (and all tensors and scalars derived from it) must be zero.
This means that the Weyl Tensor is equal to the Riemann Tensor for such solutions and so
we have the interpretation of the Weyl Tensor as that portion of the curvature which is not
due to local stress-energy.
It also implies that for vacuum solutions with nonzero cosmological constant, we have
R = D Λ / (D/2 - 1)
simpler[ Table[ Table[
In this example, L is the cosmological constant, and a replaces α as the
coefficient of the stress-energy tensor.
The "===" operator ensures that Mathematica will return a simple True or False. Of course this code
only tests for solutions on a specific coordinate chart; the process of solving Einstein's Equations
involves identifying the manifold or submanifold covered by a given chart on which a particular metric
solves the equations. Some things must still be done by human reasoning!
Rll[[ii, jj]] - R gll[[ii, jj]] / 2 + L gll[[ii, jj]],
εa b c d ... = | g |-1/2 Sign(a,b,c,d ...)
and can use the following Mathematica code to implement it in D = 4:
epsilonuuuu[ a_, b_, c_, d_] := Signature[ {a, b, c, d}] / simpler[ Sqrt[ +-Det[ gll]]];
Since Mathematica indices must start with 1, we assume that the timelike index corresponds to 4
(or in general, D). This means that in even dimensions this function will produce an extra minus
sign relative to the standard usage of 0 for the timelike index. The minus sign is selected for
metrics of Lorentz signature. We can of course define the symbol analogously for any fixed D.