Wolfram Function Repository

Instant-use add-on functions for the Wolfram Language

Function Repository Resource:

 Geodesic

Source Notebook

Compute the geodesics for a parametrized surface

Contributed by: Wolfram Staff (original content by Alfred Gray)

ResourceFunction["Geodesic"][s,{u,v},t,{u0,v0},θ0]

computes the geodesics for surface s with parameters u and v, emanating from the point parametrized by u0,v0 and proceeding in the direction θ0.

Details and Options

A geodesic is a curve that locally minimizes length traversed.
The result returned by ResourceFunction["Geodesic"] is a set of differential equations of u and v in the variable t.
The system of equations generated by ResourceFunction["Geodesic"] has the form , where are the Christoffel symbols of the second kind that can be effectively computed by the resource function ChristoffelSymbol.
Initial conditions are of the form u(0)=u0, v(0)=v0 and u'(0)=cos(θ0), v'(0)=sin(θ0).

Examples

Basic Examples (7) 

A sphere:

In[1]:=
sphere[a_][u_, v_] := {a Cos[v] Cos[u], a Cos[v] Sin[u], a Sin[v]}

Equations for geodesics on a sphere:

In[2]:=
ge = ResourceFunction["Geodesic"][sphere[1][u, v], {u, v}, t, {0, 0}, \[Theta]]
Out[2]=

A set of numerical solutions of equations for geodesics of a sphere:

In[3]:=
sg = Flatten[(NDSolve[#1, {u, v}, {t, 0, 4}] &) /@ Table[ge, {\[Theta], (2 \[Pi])/30, 2 \[Pi], (2 \[Pi])/30}], 1];

Evaluate the geodesic at a definite point:

In[4]:=
sphere[1][u[t], v[t]] /. sg [[1]] /. t -> .1
Out[4]=

Plot solutions in a plane:

In[5]:=
ParametricPlot[Evaluate[{u[t], v[t]} /. sg], {t, 0, \[Pi]}, PlotRange -> All, PlotStyle -> Hue[0]]
Out[5]=

Plot geodesics on a sphere:

In[6]:=
ParametricPlot3D[Evaluate[sphere[1][u[t], v[t]] /. sg], {t, 0, \[Pi]}]
Out[6]=

Plot the geodesic circles (locus surface points located at a given geodesic radius):

In[7]:=
Graphics3D[ Table[Line[Append[#, First[#]] &[sphere[1][u[t], v[t]] /. sg]], {t, 0, 3, 1/5}]]
Out[7]=

Scope (3) 

A torus:

In[8]:=
torus[a_, b_, c_][u_, v_] := {(a + b Cos[v]) Cos[u], (a + b Cos[v]) Sin[u], c Sin[v]}

The equations for geodesics:

In[9]:=
ResourceFunction["Geodesic"][ torus[3, 1, 1][u, v], {u, v}, t, {2, 3}, \[Theta]]
Out[9]=

Solve a geodesic for large t:

In[10]:=
Show[ParametricPlot3D[ Evaluate[ torus[3, 1, 1][u[t], v[t]] /. NDSolve[ResourceFunction["Geodesic"][torus[3, 1, 1][u, v], {u, v}, t, {\[Pi], \[Pi]}, .5], {u, v}, {t, 0, 150}]], {t, 0, 150}, Boxed -> False, Axes -> False], ParametricPlot3D[ Evaluate[torus[3, 1, 1][u, v]], {u, 0, 2 \[Pi]}, {v, 0, 2 \[Pi]}, PlotStyle -> Opacity[.5], Mesh -> False]]
Out[10]=

Solve for geodesics:

In[11]:=
sg = Flatten[ Table[NDSolve[ ResourceFunction["Geodesic"][torus[3, 1, 1][u, v], {u, v}, t, {\[Pi], \[Pi]}, \[Theta]], {u, v}, {t, 0, 4 \[Pi]}], {\[Theta], 2 \[Pi]/30, 2 \[Pi], 2 \[Pi]/30}], 1];

Plot solutions in a plane:

In[12]:=
ParametricPlot[Evaluate[{u[s], v[s]} /. sg], {s, 0, 4 \[Pi]}]
Out[12]=

Solve for a set of geodesics:

In[13]:=
Show[ParametricPlot3D[ Evaluate[torus[3, 1, 1][u, v]], {u, 0, 3 \[Pi]/2}, {v, \[Pi], 2 \[Pi]}, PlotStyle -> Opacity[.5], Mesh -> False, Boxed -> False, Axes -> False], ParametricPlot3D[ Evaluate[torus[3, 1, 1][u[t], v[t]] /. sg], {t, 0, 4 \[Pi]}]]
Out[13]=

Plot the geodesic circles:

In[14]:=
Show[Graphics[{Red, Table[Line[Append[#, First[#]] &[{u[t], v[t]} /. sg]], {t, 0, 2 \[Pi], 1/5}]}]]
Out[14]=

Plot the geodesic circles over the torus:

In[15]:=
sg = Flatten[ Table[NDSolve[ ResourceFunction["Geodesic"][torus[3, 1, 1][u, v], {u, v}, t, {\[Pi], \[Pi]}, \[Theta]], {u, v}, {t, 0, 2.5}], {\[Theta], 2 \[Pi]/200, 2 \[Pi], 2 \[Pi]/200}], 1];
In[16]:=
ta = Table[ Line[Append[#, First[#]] &[torus[3, 1, 1][u[t], v[t]] /. sg]], {t, 0, 2.5, .5}]; Show[Graphics3D[{Thickness[.005], Hue[.7], ta}], ParametricPlot3D[ Evaluate[torus[3, 1, 1][u, v]], {u, 0, 2 \[Pi]}, {v, 0, 2 \[Pi]}, PlotStyle -> Opacity[.5], Mesh -> False, Boxed -> False, Axes -> False], Boxed -> False]
Out[17]=

Define a monkey saddle:

In[18]:=
ms = Entity["Surface", "MonkeySaddle"]["ParametricEquations"][1][u, v]
Out[18]=

Get the equations for geodesics:

In[19]:=
ge = ResourceFunction["Geodesic"][ms, {u, v}, t, {0, 0}, \[Theta]] // FullSimplify
Out[19]=

Solve the equations for geodesics:

In[20]:=
sd = With[{u0 = 0, v0 = 0, a = 0, n = 100, tf = 2}, Flatten[((NDSolve[#1, {u, v}, {t, 0, 2}]) &) /@ Table[ge, {\[Theta], (2 \[Pi])/n, 2 \[Pi], (2 \[Pi])/n}], 1]];

Plot solutions for geodesics:

In[21]:=
ParametricPlot[Evaluate[{u[t], v[t]} /. sd], {t, 0, 2}]
Out[21]=

Plot the geodesics over the surface:

In[22]:=
Show[ParametricPlot3D[Evaluate[ms], {u, -2, 2}, {v, -\[Pi], \[Pi]}, PlotRange -> 1.5, PlotStyle -> Opacity[.5], Mesh -> False], ParametricPlot3D[ Evaluate[ms /. {u -> u[t], v -> v[t]} /. sd], {t, 0, 2}]]
Out[22]=

Plot the geodesic circles:

In[23]:=
Graphics[ Table[Line[Append[{u[t], v[t]} /. sd, {u[t], v[t]} /. sd[[1]]]], {t, 0, 1.7, .1}]]
Out[23]=

Plot the geodesics circles over the surface:

In[24]:=
Show[ParametricPlot3D[Evaluate[ms], {u, -2, 2}, {v, -\[Pi], \[Pi]}, PlotRange -> 1.5, PlotStyle -> Opacity[.5], Mesh -> False], Graphics3D[ Table[Line[ Append[ms /. {u -> u[t], v -> v[t]} /. sd, ms /. {u -> u[t], v -> v[t]} /. sd[[1]]]], {t, 0, 1.7, .1}]]]
Out[24]=

A pseudosphere:

In[25]:=
pseudosphere[a_][u_, v_] := a {Cos[u] Sin[v], Sin[u] Sin[v], Cos[v] + Log[.001 + Tan[v/2]]}

The geodesic equations:

In[26]:=
P = {5.3, 2.5};
In[27]:=
ResourceFunction["Geodesic"][pseudosphere[1][u, v], {u, v}, t, P, \[Theta]];

Plot the geodesic equations in 3D varying θ:

In[28]:=
plot = Show[ ParametricPlot3D[.99 pseudosphere[1][u, v], {u, 0, 2 \[Pi]}, {v, 0, \[Pi] - .01}, Mesh -> None, PlotStyle -> Opacity[.5]], Graphics3D[{PointSize[.03], Point[pseudosphere[1] @@ P]}]];
In[29]:=
Manipulate[ Show[ParametricPlot3D[ Evaluate[ pseudosphere[1][u[t], v[t]] /. NDSolve[ResourceFunction["Geodesic"][ pseudosphere[1][u, v], {u, v}, t, P, \[Theta]], {u, v}, {t, 0, 10}]], {t, 0, 10}, PlotRange -> {{-1, 1}, {-1, 1}, {0, 3.5`}}], plot], {{\[Theta], 1.3}, 0, \[Pi]}]
Out[29]=

A top view of a solution:

In[30]:=
Module[{\[Theta] = 1.3`}, ParametricPlot[ Evaluate[ Most[First[ pseudosphere[1][u[t], v[t]] /. NDSolve[ResourceFunction["Geodesic"][ pseudosphere[1][u, v], {u, v}, t, P, \[Theta]], {u, v}, {t, 0, 120}]]]], {t, 0, 120}, PlotRange -> 1]]
Out[30]=

Publisher

Enrique Zeleny

Version History

  • 1.0.0 – 08 April 2020

Source Metadata

  • Citation:
    • Abbena, E., Salamon, S., Gray, A., Modern Differential Geometry of Curves and Surfaces with Mathematica, 3rd ed. Boca Raton, FL: Chapman & Hall/CRC, 2006.

License Information