- #1
- 2,355
- 10
For those who use GRTensorII, here is a simple and effective method for
This method requires no particular knowledge of GRTensorII.
I'll explain using a simple example, the Melvin electrovacuum, which models a cylindrically symmetric magnetostatic field, in which the EM field energy is the sole source of the gravitational field. Here is an (extensively commented) GRTensorII file you can use to define the static coframe in this spacetime:
Now fire up Maple and in your session type the following commands:
Here is what we did:
[tex]
\begin{array}{rcl}
\dot{t} & = & \frac{E}{(1+q^2 \, r^2)^2} \\
\dot{z} & = & \frac{A}{(1+q^2 \, r^2)^2} \\
\dot{r} & = & \pm \sqrt{\frac{E^2-A^2}{(1+q^2 \, r^2)^3} - \frac{L^2 \, (1+q^2 \, r^2)}{r^2}} \\
\dot{\phi} & = & \frac{L \, (1+q^2 \,r^2)^2}{r^2}
\end{array}
[/tex]
where [itex]E, \, A, \, L[/itex] are the invariants of null geodesic motion, and where the choice of sign determines an ingoing or outgoing motion.
We can immediately write down a null geodesic vector field whose integral curves form a null geodesic congruence:
[tex]
\vec{k} =
\frac{1}{(1+q^2 \, r^2)^2} \; \left( E \, \partial_t + A \, \partial_z \right)
\pm \sqrt{ \left( \frac{E^2-A^2}{(1+q^2 \, r^2)^3} - \frac{L^2 \, (1+q^2 \, r^2)}{r^2} \right)} \; \partial_r
+ \frac{L \, (1+q^2 \,r^2)^2}{r^2} \; \partial_\phi
[/tex]
Far from the axis of cylindrical symmetry, this reduces to
[tex]
\vec{k}_{\rm far} = \left( E \, \partial_t + A \, \partial_z \right)
\pm \sqrt{E^2-A^2} \; \partial_r + \frac{L}{r^2} \; \partial_\phi
[/tex]
which corresponds to a plane wave in cylindrical coordinates in flat spacetime, with a constant wave vector pointing in an arbitrary direction (determined by E,A,L). This also shows that for our null geodesics, we can loosely interpret E as "the energy of a photon measured at spatial infinity", and likewise L as "linear momentum parallel to the axis of cylindrical symmetry" and "angular momentum about the axis of cylindrical symmetry". This "photon" doesn't remain "at spatial infinity", but because E, A, L are constants of geodesic motion, these parameters still have meaning as it approaches its point of closest approach to the axis of cylindrical symmetry, even though we interpret them in terms of measurements made when the "photon" was "at spatial infinity".
By using [itex]-1 = ds^2[/itex] instead, you can find first integrals for an arbitrary timelike geodesic in the Melvin electrovacuum. Using this you can study such interesting phenomena as the light bending and geodesic precession due to the gravitational attraction of the magnetic field, which is strongly concentrated near the symmetry axis.
In spacetimes with smaller Lie groups of self-isometries, it may not be possibly to obtain in such trivial fashion a complete list of first integrals for null, timelike, or spacelike geodesics, but one almost always obtains some valuable information. The case of the Kerr vacuum is remarkable because, although we have only a two dimensional Lie group of self-isometries, and thus we can expect only three first integrals for (null, timelike or spacelike) geodesicsl. In this case, it turns out that there also exists a more subtle symmetry, given by a nontrivial Killing tensor field. Brandon Carter was able to exploit this to solve the geodesic equations for the Kerr vacuum. Trivial Killing tensor fields are common and useless; nontrivial Killing tensor fields seem to be quite rare, in fact, IIRC no example is known other than the Kerr vacuum!
At this point, if you want to examine the Killing equations with which we started, you can try this
The first command instructs Maple to abbreviate partials by subscripts. The second prints out each equation, one line per equation. From this you can see that as I said, this is a first order system of PDEs. Compare with the "triangularized" version output by casesplit!
I should say that casesplit uses nontrivial Groebner-basis like methods in differential rings, which we can think of as a sophisticated analog of using Gaussian reduction to "triangularize" a system of linear equations. In both cases, the point is that we can systematically solve a triangularized system by solving for the variable for which we have obtained equations involving only that variable, then plugging the result into the equations involving that variable and one other, solving, and so forth until we are done.
There is a general principle in computational algebra which says that the fastest methods tend to be non-deterministic. Thus, casesplit may sometimes present results in a different order, so the directions I gave above must be interpreted with this in mind. Also, casesplit is not perfect and sometimes the above method will fail to find all the solutions it should, so be on the lookout for such problems.
- determining the Lie algebra of Killing vector fields
- using the result to find first integrals of the geodesic equations via Noether's theorem on symmetries
This method requires no particular knowledge of GRTensorII.
I'll explain using a simple example, the Melvin electrovacuum, which models a cylindrically symmetric magnetostatic field, in which the EM field energy is the sole source of the gravitational field. Here is an (extensively commented) GRTensorII file you can use to define the static coframe in this spacetime:
Code:
Ndim_ := 4:
x1_ := t:
x2_ := z:
x3_ := r:
x4_ := phi:
eta11_ := -1:
eta22_ := 1:
eta33_ := 1:
eta44_ := 1:
bd11_ := -(1+q^2*r^2):
bd22_ := (1+q^2*r^2):
bd33_ := (1+q^2*r^2):
bd44_ := r/(1+q^2*r^2):
Info_ := `Melvin non-null electrovacuum (1964); cyl chart; static obsvrs`:
# Chart -infty < t,z < infty, 0 < r < infty, -Pi < u < Pi
# (Discovered by Bonnor 1954 and Witten 1962 and then Melvin 1964)
#
# Four dimensional Lie algebra of Killing vector fields:
# @_t time translation
# @_z spatial translation in z direction
# z @_t + t @_z boost in z direction
# @_phi axial rotation
# So this is boost/rotation symmetric as well as cyl sym static.
# Geodesic equations:
# From Lagrangian we obtain quadratic invariant
# epsilon = (1+q^2*r^2)^2*(-t1^2+z1^2+r1^2)+phi1^2*r^2/(1+q^2*r^2)^2
# From Killing vectors we obtain four linear quadratic invariants:
# B = t1*(1+q^2*r^2)^2
# P = z1*(1+q^2*r^2)^2
# K = (z*t1+t*z1)*(1+q^2*r^2)^2
# J = phi1*r^2/(1+q^2*r^2)^2
# Solving { epsilon=0,B,P,J} for {t1,z1,r1,phi1} we obtain
# t1 = B/(1+q^2*r^2)^2
# z1 = P/(1+q^2*r^2)^2
# r1 = sqrt((B^2-P^2)*r^2-(1+q^2*r^2)^4*J^2)/(1+q^2*r^2)^2/r
# phi1 = J/r^2*(1+q^2*r^2)^2
# Helical null congruence
# J/sqrt(B^2-P^2) = r/(1+q^2*r^2)^2
# 1/(1+q*r^2) (B e_1 + P e_2) + sqrt(B^2-P^2)/r e_4
#
# Physical experience of static observers:
# Acceleration vector
# -2*q^2*r/(1+q^2*r^2)^2 e_3
# Expansion and vorticity tensor vanish
# Vector Potential
# -q/2*(1+q^2*r^2) @_phi
# Axial regular magnetic field concentrated near z = 0
# q/(1+q^2*r^2)^2 e_2
# Einstein tensor
# G^(ab) = 4*q^2/(1+q^2*r^2) diag(1,-1,1,1)
# Principal gravitational field invariants
# RR = 64*q^4*(3*r^4*q^4-6*r^2*q^2+5)/(1+r^2*q^2)^8
# Rstar = 0
# diRiem = 256*r^2*q^8*(133-126*r^2*q^2+45*r^4*q^4)/(1+r^2*q^2)^12
# CC = 192*q^4*(r*q-1)^2*(r*q+1)^2/(1+r^2*q^2)^8
# Cstar = 0
# diWeyl = 768*q^8*r^2*(15*r^4*q^4-42*r^2*q^2+31)/(1+r^2*q^2)^12
# Curvature everywhere finite and concentrated near r = 0
# Tidal tensor or electro-Riemann shows less axial tension than electro-Weyl tensor
# E_(ab) = EW_(ab) + 4*q^2/(1+q^2*r^2)^4 diag(1,0,0)
# This describes tidal accelerations of clouds of electrically neutral test particles.
# Weyl tensor is "pure electric"
# EW_(ab) = 2*q^2*(1-r^2*q^2)/(1+q^2*r^2)^4 diag(-2,1,1)
# Orthogonal hyperslices have three-dimensional Riemann tensor
# rv^2_(323) = rv^2_(424) = 2*q^2*(q^2*r^2-1)/(q*2*r^2+1)^4
# rv^3_(434) = 4*q^2*(q^2*r^2-2)/(q^2*r^2+1)^4
#
# For cylindrically symmetric functions, the Laplace-Beltrami equation
# becomes the usual wave equation
# h_(tt) = h_(rr) + (h_r)/r
# A simple exact solution is
# h = cos(a*t)*BesselJ(0,a*r)
# The corresponding vector field is
# sin(a*t)*BesselJ(0,a*r)/(1+q^2*r^2) e_1
# -a*cos(a*t)*BesselJ(1,a*r)/(1+q^2*r^2) e_3
#
# References:
# Bonnor, prsl A67 (1954): 225.
# Melvin, Phys. Lett. 8 (1964): 65.
Code:
qload(Melvin_nnevac_mg_stataxs_reg_cyl);
grcalcalter(x(up),basisv(dn),basisv(up),g(dn,dn),g(up,up),1):grdisplay(_);
grcalcalter(LieD[kv0,g(dn,dn)],1):grdisplay(_);
keqs := [seq(seq(grcomponent(LieD[kv0,g(dn,dn)], [grcomponent(x(up),[j]),grcomponent(x(up),[k])]),j=1..4),k=1..4)]:
# Enter
# T(t,z,r,phi)
# Z(t,z,r,phi)
# R(t,z,r,phi)
# Phi(t,z,r,phi)
op(%)[1];
pdsolve(%,{T,Z,R,Phi});
subs(%,[T,Z,R,Phi](t,z,r,phi));
sol := subs(_C1=c[1],_C2=c[2],_C3=c[3],_C4=c[4],%);
# Output
# [c[1]*z+c[4], c[1]*t+c[2], 0, c[3]]
subs(c[1]=1,c[2]=0,c[3]=0,c[4]=0,sol);
subs(c[1]=0,c[2]=1,c[3]=0,c[4]=0,sol);
subs(c[1]=0,c[2]=0,c[3]=1,c[4]=0,sol);
subs(c[1]=0,c[2]=0,c[3]=0,c[4]=1,sol);
# Output shows four dimensional Lie algebra of Killing vectors
# z @_t + t @_z
# @_t, @_z, @_phi
# Transitive on cylinders r=r0
#
multiply(matrix([[1,0,0,0]]), matrix( [[-(1+q^2*r^2)^2,0,0,0],[0,(1+q^2*r^2)^2,0,0],[0,0,(1+q^2*r^2)^2,0],[0,0,0,r^2/(1+q^2*r^2)^2]]), matrix( [[t1],[z1],[r1],[phi1]]))[1,1] = -EE; t1=solve(%,t1)
# Output shows t1 = EE/(1+q^2*r^2)^2
multiply(matrix([[0,1,0,0]]), matrix( [[-(1+q^2*r^2)^2,0,0,0],[0,(1+q^2*r^2)^2,0,0],[0,0,(1+q^2*r^2)^2,0],[0,0,0,r^2/(1+q^2*r^2)^2]]), matrix( [[t1],[z1],[r1],[phi1]]))[1,1] = A; z1=solve(%,z1);
# Output shows z1 = A/(1+q^2*r^2)^2
multiply(matrix([[0,0,0,1]]), matrix( [[-(1+q^2*r^2)^2,0,0,0],[0,(1+q^2*r^2)^2,0,0],[0,0,(1+q^2*r^2)^2,0],[0,0,0,r^2/(1+q^2*r^2)^2]]), matrix( [[t1],[z1],[r1],[phi1]]))[1,1] = L; phi1=solve(%,phi1); lprint(%);
# Output shows phi1 = L*(1+q^2*r^2)^2/r^2
multiply(matrix([[z,t,0,0]]), matrix( [[-(1+q^2*r^2)^2,0,0,0],[0,(1+q^2*r^2)^2,0,0],[0,0,(1+q^2*r^2)^2,0],[0,0,0,r^2/(1+q^2*r^2)^2]]), matrix( [[t1],[z1],[r1],[phi1]]))[1,1] = B; factor(%); -t*z1+z*t1=solve(%,-z*t1+t*z1): simplify(%,factor);
# Output shows z*t1-t*z1 = B/(1+q^2*r^2)^2
[t1 = EE/(1+q^2*r^2)^2, z1 = A/(1+q^2*r^2)^2, phi1 = L*(1+q^2*r^2)^2/r^2];
grcomponent(g(dn,dn),[t,t])*t1^2 + grcomponent(g(dn,dn),[z,z])*z1^2 + grcomponent(g(dn,dn),[r,r])*r1^2 + grcomponent(g(dn,dn),[phi,phi])*phi1^2;
subs(%%,%): simplify(%,factor): collect(%,EE): collect(%,A): collect(%,L);
lprint(%);
# Output shows
# r1^2 = (EE^2-A^2)/(1+q^2*r^2)^3 - L^2*(1+q^2*r^2)/r^2
- We computed the Killing equations for a generic vector field
[tex]
T \, \partial_t + Z \, \partial_z + R \, \partial_r + \Phi \, \partial_\phi
[/tex]
You can call the undetermined functions [itex]T,Z,R,\Phi[/itex] anything you like, of course. - We put these PDEs in a list.
- The key step is feeding this list of first order PDEs to Maple's powerful "casesplit" commend, which "triangularizes" them. (The form of the Killing equations guarantees than in principle, this can be done.)
- We use Maple's "pdesolve" command solves these system.
- We extract a list of four vector fields which generate the Lie algebra,
- the boost [itex]z \, \partial_z + t \, \partial_z[/itex]
- time translation [itex]\partial_t[/itex]; this turns out to be a vorticity free timelike vector field, so the constant t slices form a family of spatial hyperslices everywhere orthogonal to the world lines of the static observers,
- translation along axis of cylindrical symmetry [itex]\partial_z[/itex], a spacelike Killing vector field,
- rotation about the axis of cylindrical symmetry [itex]\partial_\phi[/itex], a spacelike cyclic Killing vector field
- Courtesy of Noether's theorem on variational symmetries of a system of PDEs determined by a Lagrangian, we immediately obtain several first integrals (more than we need, actually). You can use any notation you like for the first derivatives of a geodesic wrt an arc length parameter you like, but I find t1, z1, r1, phi1 particularly convenient for working with Lagrangians (when at various moments we need to consider derivatives as independent variables, or to work purely algebraically, as here). I used EE for energy in the Maple session so that there is no possibility of confusion with the defined constant [itex]\exp(1)[/itex]. In particular, we obtain first integrals for [itex]\dot{t}, \; \dot{z}, \; \dot{\phi}[/itex], but we still need a first integral for [itex]\dot{r}[/itex].
- We plug these into the line element [itex] 0 = ds^2[/itex] to obtain the first integral for [itex]\dot{r}[/itex], in the case of null geodesics.
[tex]
\begin{array}{rcl}
\dot{t} & = & \frac{E}{(1+q^2 \, r^2)^2} \\
\dot{z} & = & \frac{A}{(1+q^2 \, r^2)^2} \\
\dot{r} & = & \pm \sqrt{\frac{E^2-A^2}{(1+q^2 \, r^2)^3} - \frac{L^2 \, (1+q^2 \, r^2)}{r^2}} \\
\dot{\phi} & = & \frac{L \, (1+q^2 \,r^2)^2}{r^2}
\end{array}
[/tex]
where [itex]E, \, A, \, L[/itex] are the invariants of null geodesic motion, and where the choice of sign determines an ingoing or outgoing motion.
We can immediately write down a null geodesic vector field whose integral curves form a null geodesic congruence:
[tex]
\vec{k} =
\frac{1}{(1+q^2 \, r^2)^2} \; \left( E \, \partial_t + A \, \partial_z \right)
\pm \sqrt{ \left( \frac{E^2-A^2}{(1+q^2 \, r^2)^3} - \frac{L^2 \, (1+q^2 \, r^2)}{r^2} \right)} \; \partial_r
+ \frac{L \, (1+q^2 \,r^2)^2}{r^2} \; \partial_\phi
[/tex]
Far from the axis of cylindrical symmetry, this reduces to
[tex]
\vec{k}_{\rm far} = \left( E \, \partial_t + A \, \partial_z \right)
\pm \sqrt{E^2-A^2} \; \partial_r + \frac{L}{r^2} \; \partial_\phi
[/tex]
which corresponds to a plane wave in cylindrical coordinates in flat spacetime, with a constant wave vector pointing in an arbitrary direction (determined by E,A,L). This also shows that for our null geodesics, we can loosely interpret E as "the energy of a photon measured at spatial infinity", and likewise L as "linear momentum parallel to the axis of cylindrical symmetry" and "angular momentum about the axis of cylindrical symmetry". This "photon" doesn't remain "at spatial infinity", but because E, A, L are constants of geodesic motion, these parameters still have meaning as it approaches its point of closest approach to the axis of cylindrical symmetry, even though we interpret them in terms of measurements made when the "photon" was "at spatial infinity".
By using [itex]-1 = ds^2[/itex] instead, you can find first integrals for an arbitrary timelike geodesic in the Melvin electrovacuum. Using this you can study such interesting phenomena as the light bending and geodesic precession due to the gravitational attraction of the magnetic field, which is strongly concentrated near the symmetry axis.
In spacetimes with smaller Lie groups of self-isometries, it may not be possibly to obtain in such trivial fashion a complete list of first integrals for null, timelike, or spacelike geodesics, but one almost always obtains some valuable information. The case of the Kerr vacuum is remarkable because, although we have only a two dimensional Lie group of self-isometries, and thus we can expect only three first integrals for (null, timelike or spacelike) geodesicsl. In this case, it turns out that there also exists a more subtle symmetry, given by a nontrivial Killing tensor field. Brandon Carter was able to exploit this to solve the geodesic equations for the Kerr vacuum. Trivial Killing tensor fields are common and useless; nontrivial Killing tensor fields seem to be quite rare, in fact, IIRC no example is known other than the Kerr vacuum!
At this point, if you want to examine the Killing equations with which we started, you can try this
Code:
declare(T(t,z,r,phi), Z(t,z,r,phi), R(t,z,r,phi), Phi(t,z,r,phi);
for guy in keqs do guy od;
I should say that casesplit uses nontrivial Groebner-basis like methods in differential rings, which we can think of as a sophisticated analog of using Gaussian reduction to "triangularize" a system of linear equations. In both cases, the point is that we can systematically solve a triangularized system by solving for the variable for which we have obtained equations involving only that variable, then plugging the result into the equations involving that variable and one other, solving, and so forth until we are done.
There is a general principle in computational algebra which says that the fastest methods tend to be non-deterministic. Thus, casesplit may sometimes present results in a different order, so the directions I gave above must be interpreted with this in mind. Also, casesplit is not perfect and sometimes the above method will fail to find all the solutions it should, so be on the lookout for such problems.
Last edited: