Metric Perturbations with xAct

As you have surely noticed by now, computing metric perturbations by hand is painful. Luckily, people have had this problem before, and automatised it. This notebook gives you a small illustration of one way to compute them with Mathematica, using the free package xAct.
To get started, first you have to install the package. Go to xact.es and follow the instructions there. Once it is properly installed, this notebook should work. On the webpage, you also find the documentation, and many introductory talks and examples.

Reading the package

To be able to work with xAct, we first have to load the package:

Block[{Print},
<<xAct`xTras`;
];

The Block environment suppresses some output that the package produces. It is not important for us.

Defining a manifold and metric perturbations

As a next step, we define the spacetime manifold. Our spacetime has the name “spacetime” (how fitting, right?), is four-dimensional, and has indices {μ,ν,ρ,σ,α,β,γ,δ,κ,λ,τ,ω}:

DefManifold[spacetime,4,{μ,ν,ρ,σ,α,β,γ,δ,κ,λ,τ,ω}];

You can add more indices later if you need them. Also, xAct creates more indices automatically if it needs them!
Next, let us define a metric with name “metric” on spacetime, with a Lorentzian signature (that’s the -1 in the first argument). The metric comes with lower indices, these are denoted with a minus sign in xAct. xAct automatically defines the metric-compatible covariant derivative that we will call “CD”, and a notation for it.

DefMetric[-1,metric[-μ,-ν],CD,{";","𝒟"}];

As you can see from the wall of text, xAct automatically defines all curvature tensors that one can derive from the metric. Finally, let us define metric perturbations that we call “pert” (that will be printed as h), and that come with a perturbation parameter ϵ:

DefMetricPerturbation[metric,pert,ϵ,PrintAs->"h"];

For reasons that you will see very soon, let us also define a symmetric tensor that we call hpert:

DefTensor[hpert[-μ,-ν],spacetime,Symmetric[{-μ,-ν}],PrintAs->ColorString["𝒽",Lighter@Blue]];

Now we are all set up to do some computations!

Basic perturbations

Let us first understand how xAct defines metric perturbations. For this, we use “Perturbed”, which takes perturbations of its first argument up to a given order which is the second argument:

Perturbed[metric[-μ,-ν],3]

You can see now that xAct defines very general perturbations! We will restrict ourselves to a linear parameterisation (this is the reason why we introduced hpert):

metricpertrule={pert[LI[1],X__]:>hpert[X],pert[LI[n_/;n>1],__]:>0};

Perturbed[metric[-μ,-ν],3]/.metricpertrule

Using our rule, we can throw away higher order perturbations.

We can use Perturbed to compute perturbations of curvature tensors. However, Perturbed first produces rather abstract perturbations when it acts on anything that is not the metric:

Perturbed[RicciScalarCD[],1]

To translate this into metric perturbations, we use ExpandPerturbation:

(Perturbed[RicciScalarCD[],1]//ExpandPerturbation)/.metricpertrule

Neat! But somehow, things are not yet simplified. For this, we can use the command CollectTensors, which does two things: it brings all tensors into a canonical form, and then simplifies the prefactors:

(Perturbed[RicciScalarCD[],1]//ExpandPerturbation)/.metricpertrule//CollectTensors

Often, we only want the perturbation and not the zero-th order term. Then we can use Perturbation. This picks out the n-th perturbation (without the factor of 1/n! and without ϵ):

(Perturbation[RicciScalarCD[],1]//ExpandPerturbation)/.metricpertrule//CollectTensors

Computing Einstein’s equations

We can now compute Einstein’s equations! For this, let us first define Newton’s constant and the cosmological constant. We declare them as constant symbols so that xAct knows that they are indeed constants:

DefConstantSymbol[GNewton,PrintAs->"(Subscript[G, N])"];
DefConstantSymbol[CC,PrintAs->"Λ"];

Now let us compute the linear perturbation of the Einstein-Hilbert action (for linear perturbations, we can drop the second argument in Perturbation):

SEHpert=(Perturbation[(Sqrt[-Detmetric[]])(((RicciScalarCD[]-2CC)/(16π GNewton)))]//ExpandPerturbation)/.metricpertrule//CollectTensors

We can see that there is still the determinant of the metric, and in the following, we will divide by it (this is part of the definition of a functional derivative). To compute the equations of motion, we have to take the variational derivative with respect to h. The first brackets indicate the field with respect to which we take the derivative of, and the derivative which we can partially integrate. The second bracket has the expression that we take the derivative of. For convenience, we also multiply by (-16π G) and get:

-16π GNewton VarD[hpert[-μ,-ν],CD][(SEHpert/(Sqrt[-Detmetric[]]))]//CollectTensors

Einstein’s equations! Fantastico. Now, of course there is a direct command for this in xAct:

-16π GNewton VarL[metric[-μ,-ν]][((RicciScalarCD[]-2CC)/(16π GNewton))]//CollectTensors

Computing equations of motion for general theories

Let us finally look at the equations of motion for more general theories. First, we define couplings as constant symbols:

DefConstantSymbol[gR2,PrintAs->"(Subscript[g, (R^2)])"];
DefConstantSymbol[gRic2,PrintAs->"(Subscript[g, (Ric^2)])"];
DefConstantSymbol[gRiem2,PrintAs->"(Subscript[g, (Riem^2)])"];

Let’s compute equations of motion!

eom1=VarL[metric[-μ,-ν]][((RicciScalarCD[]-2CC)/(16π GNewton))+gR2 ((RicciScalarCD[])^2)+gRic2 RicciCD[α,β]RicciCD[-α,-β]+gRiem2 RiemannCD[α,β,γ,δ]RiemannCD[-α,-β,-γ,-δ]]//CollectTensors

These can still be simplified. Let us sort covariant derivatives (SortCovDsToDiv), first so that we can use the differential Bianchi identity, and then actually use the Bianchi identity (CurvatureRelationsBianchi):

eom1//SortCovDsToDiv[RiemannCD]//SortCovDsToDiv[RicciCD]//CollectTensors
eom2=%/.CurvatureRelationsBianchi[CD]//CollectTensors

Due to using the Bianchi identity for the Riemann tensor, we have more terms that we can reduce this way, so let’s repeat this step:

eom2//SortCovDsToDiv[RicciCD]//CollectTensors
eom3=%/.CurvatureRelationsBianchi[CD]//CollectTensors

Finally, let us sort all remaining covariant derivatives into a canonical order (SortCovDs):

eom4=eom3//SortCovDs//CollectTensors

Fun fact: in four dimensions only, there is an additional simplification. Namely, if two Riemann tensors are contracted over three pairs of indices, the expression can be written as their full contraction plus terms that involve Ricci tensors and Ricci scalars. This can be derived by an “over-antisymmetrisation”, that is, an antisymmetrisation over five indices which vanishes in four dimensions. Specifically, we consider (Subscript[g, μ[ν])(Subscript[(R^(αβ)), αβ])(Subscript[(R^(γδ)), γδ]]) which is zero in four dimensions:

antisymeqzero=Antisymmetrize[metric[μ,ν]RiemannCD[κ,λ,-α,-β]RiemannCD[τ,ω,-γ,-δ],{ν,-α,-β,-γ,-δ}]metric[-κ,α]metric[-λ,β]metric[-τ,γ]metric[-ω,δ]//CollectTensors

We can solve the first identity (don’t worry about the options):

riemsq4did=Quiet[SolveTensors[antisymeqzero==0,UseSymmetries->False,MetricOn->None][[1]]]

Using this identity for the equations of motion:

eom5=eom4/.riemsq4did//CollectTensors

We can compute the linear perturbations from this equation and then go to Minkowski space with ToFlat:

Perturbation[eom5]/.CC->0//ExpandPerturbation//ToFlat//CollectTensors
%/.metricpertrule//SortCovDs//SortCovDsToBox[hpert]//SortCovDsToDiv[hpert]//ToFlat//CollectTensors