2
$\begingroup$

I'm using the Christoffel Symbols found on this link to generate a set of three coupled differential equations as solutions to the Geodesic Equation of a Sphere, I have: $$\frac{d^2r}{d\lambda^2}-r\frac{d\theta}{d\lambda}\frac{d\theta}{d\lambda}-r\space Sin(\theta)^2\frac{d\phi}{d\lambda}\frac{d\phi}{d\lambda}=0$$ $$\frac{d^2\theta}{d\lambda^2}+\frac{2}{r}\frac{dr}{d\lambda}\frac{d\theta}{d\lambda}-Cos(\theta)Sin(\theta)\frac{d\phi}{d\lambda}\frac{d\phi}{d\lambda}=0$$ $$\frac{d^2\phi}{d\lambda^2}+\frac{2}{r}\frac{d\theta}{d\lambda}\frac{dr}{d\lambda}+2\space Cot(\theta)\frac{d\phi}{d\lambda}\frac{d\theta}{d\lambda}=0$$ So far, so good. My eventual goal here is to find the geodesic on a sphere that expands with time (that is, a preliminary step to a full FLRW model). Now I want to program this formula in Mathematica to see what happens on the surface of a sphere. Here's the code:

X[r_, theta_, phi_] := r*Cos[phi]*Sin[theta];
Y[r_, theta_, phi_] := r*Sin[phi]*Sin[theta];
Z[r_, theta_, phi_] := r*Cos[theta];
R[r_, theta_, phi_] := {X[r, theta, phi], Y[r, theta, phi], Z[r, theta, phi]};

equations =
{
  r[0] == 1.,
  theta[0] == 1.,
  phi[0] == 1.,
  Derivative[1][r][0] == 0.,
  Derivative[1][theta][0] == 0.5,
  Derivative[1][phi][0] == 1.,
  Derivative[2][r][lambda] - r[lambda]*D[theta[lambda], lambda]*D[theta[lambda], lambda] - r[lambda]*Sin[theta[lambda]]^2*D[phi[lambda], lambda]*D[phi[lambda], lambda] == 0,
  Derivative[2][theta][lambda] + (2/r[lambda])*D[r[lambda], lambda]*D[theta[lambda], lambda] - Cos[theta[lambda]]*Sin[theta[lambda]]*D[phi[lambda], lambda]*D[phi[lambda], lambda] == 0,
  Derivative[2][phi][lambda] + (2/r[lambda])*D[phi[lambda], lambda]*D[r[lambda], lambda] + 2*Cot[theta[lambda]]*D[phi[lambda], lambda]*D[theta[lambda], lambda] == 0
};
S[lambda_] = {r[lambda], theta[lambda], phi[lambda]} /. First[NDSolve[equations, {r, theta, phi}, {lambda, 0, 5}]];
surface = ParametricPlot3D[{R[1., theta, phi]}, {theta, 0, Pi}, {phi, 0, 2*Pi}];
Show[{surface, ParametricPlot3D[R[S[lambda][[1]], S[lambda][[2]], S[lambda][[3]]], {lambda, 0, 5},
  PlotStyle -> Directive[Red, Thick]]}, ImageSize -> Medium]

Now, here's where I need some help. If I plot this, I get:

enter image description here

Which appears to be a tangent vector. I was expecting the $r$ value to remain constant around the sphere since I set the initial distance to 1. and the velocity to 0. If I change the second derivative formula to:

  Derivative[2][r][lambda] == 0,

Then I get the expected Geodesic: enter image description here Can someone explain to me what is going on? Did I make a mistake? Or should I not expect the $r$ coordinate to follow the Geodesic because it's normal to the surface?

| improve this question | | | | |
$\endgroup$
6
$\begingroup$

There is an obvious typo in the derivation of the equations of motion on a sphere. First, we write down the coordinates and the metric tensor for a sphere with a given radius r0=1 immersed in 4D:

n = 4; r0=1;
coord = {t, r, \[Theta], \[Phi]};
metric ={{-1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, r0^2, 0}, {0, 0, 0, 
  r0^2 Sin[\[Theta]]^2}};

We calculate the Christoffel symbols for this case

inversemetric = Simplify[Inverse[metric]]

affine := 
  affine = Simplify[
    Table[(1/2)*
      Sum[inversemetric[[i, 
         s]]*(D[metric[[s, j]], coord[[k]]] + 
          D[metric[[s, k]], coord[[j]]] - 
          D[metric[[j, k]], coord[[s]]]), {s, 1, n}], {i, 1, n}, {j, 
      1, n}, {k, 1, n}]];


listaffine := 
  Table[If[UnsameQ[affine[[i, j, k]], 
     0], {ToString[\[CapitalGamma][i, j, k]], affine[[i, j, k]]}], {i,
     1, n}, {j, 1, n}, {k, 1, j}];


TableForm[Partition[DeleteCases[Flatten[listaffine], Null], 2], 
 TableSpacing -> {2, 2}]

And so there are two non-zero components $$\Gamma[3,4,4]=-\cos\theta \sin\theta, \Gamma[4,4,3]=\cot \theta$$ Now we find the equations of motion $$\frac {d^2x^i}{ds^2}+\Gamma ^i_{kl}\frac{dx^k}{ds}\frac {dx^l}{ds}=0$$ Put

v = {t'[s], r'[s], \[Theta]'[s], \[Phi]'[s]};

And then calculate

Table[D[v[[i]], s] + ((affine.v).v)[[i]] == 0, {i, n}]

Finally, we have a system of equations

    eqs = {t''[s] == 0, 
   r''[s] == 
    0, -Cos[\[Theta][s]] Sin[\[Theta][s]] \[Phi]'[s]^2 + \[Theta]''[
      s] == 0, 
   2 Cot[\[Theta][s]] \[Theta]'[s] \[Phi]'[s] + \[Phi]''[s] == 0};

ic = {t[0] == 0, t'[0] == 1, 
   r[0] == 1., \[Theta][0] == 1., \[Phi][0] == 1., 
   r'[0] == 0., \[Theta]'[0] == 0.5, \[Phi]'[0] == 1.};

Now we can display the trajectory on the sphere

X[r_, theta_, phi_] := r*Cos[phi]*Sin[theta];
Y[r_, theta_, phi_] := r*Sin[phi]*Sin[theta];
Z[r_, theta_, phi_] := r*Cos[theta];
R[r_, theta_, phi_] := {X[r, theta, phi], Y[r, theta, phi], 
   Z[r, theta, phi]};

S1 = NDSolveValue[
   Flatten[{eqs, ic}], {r[s], \[Theta][s], \[Phi][s]}, {s, 0, 5}];
surface = 
  ParametricPlot3D[{R[1., theta, phi]}, {theta, 0, Pi}, {phi, 0, 
    2*Pi}];
Show[{surface, 
  ParametricPlot3D[R[S1[[1]], S1[[2]], S1[[3]]], {s, 0, 5}, 
   PlotStyle -> Directive[Red, Thick]]}, ImageSize -> Medium]

Figure 1

| improve this answer | | | | |
$\endgroup$
  • $\begingroup$ Could you direct me to my typo, please? Your solution is excellent, but it would be useful to me to see where my original solution broke. $\endgroup$ – Quarkly Dec 29 '19 at 18:29
  • $\begingroup$ You used Christoffel symbols for spherical coordinates in 3D, but not for a sphere. $\endgroup$ – Alex Trounev Dec 29 '19 at 18:38
  • $\begingroup$ Your equations describe a straight line in 3D (geodesic in empty space). $\endgroup$ – Alex Trounev Dec 29 '19 at 21:19
  • $\begingroup$ I believe your metric is incorrect. It should be: metric = {{1, 0, 0}, {0, r^2, 0}, {0, 0, r^2 Sin[\[Theta]]^2}}; $\endgroup$ – Quarkly Dec 29 '19 at 21:29
  • $\begingroup$ You wrote down the metric of empty space in spherical coordinates. I described the metric of a sphere immersed in 4D. $\endgroup$ – Alex Trounev Dec 29 '19 at 21:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.