Jacobi Fields

Illustration of Jacobi Fields

This tutorial illustrates the usage of Jacobi Fields within Manopt.jl. For this tutorial you should be familiar with the basic terminology on a manifold like the exponential and logarithmic map as well as geodesics.

We first initialize the manifold

using Manopt

and we define some colors from Paul Tol

using Colors
black = RGBA{Float64}(colorant"#000000")
TolVibrantOrange = RGBA{Float64}(colorant"#EE7733")
TolVibrantCyan = RGBA{Float64}(colorant"#33BBEE")
TolVibrantTeal = RGBA{Float64}(colorant"#009988")

Assume we have two SnPoints on the equator of the Sphere(2) $\mathcal M = \mathbb S^2$

M = Sphere(2)
x,y = [ SnPoint([1.,0.,0.]), SnPoint([0.,1.,0.])]
2-element Array{SnPoint{Float64},1}:
 Sn([1.0, 0.0, 0.0])
 Sn([0.0, 1.0, 0.0])

their connecting geodesic (sampled at 100 points)

geodesicCurve = geodesic(M,x,y,100);
asyResolution = 2

looks as follows using renderAsymptote with the asyExportS2Signals export

renderAsymptote("jacobiGeodesic.asy",asyExportS2Signals;
    render = asyResolution,
    curves=[geodesicCurve], points = [ [x,y] ],
    colors=Dict(:curves => [black], :points => [TolVibrantOrange]),
    dotSize = 3.5, lineWidth = 0.75, cameraPosition = (1.,1.,.5)
)

A geodesic connecting two points on the equator

where $x$ is on the left. Then this tutorial solves the following task:

Given a direction $\xi_x\in T_x\mathcal M$, for example the SnTVector

ξx = SnTVector([0.,0.4,0.5])
SnT([0.0, 0.4, 0.5])

we move the start point $x$ into, how does any point on the geodesic move?

Or mathematically: Compute $D_x g(t; x,y)$ for some fixed $t\in[0,1]$ and a given direction $\xi_x$. Of course two cases are quite easy: For $t=0$ we are in $x$ and how $x$ “moves” is already known, so $D_x g(0;x,y) = \xi$. On the other side, for $t=1$, $g(1; x,y) = y$ which is fixed, so $D_x g(1; x,y)$ is the zero tangent vector (in $T_y\mathcal M$).

For all other cases we employ a jacobiField, which is a (tangent) vector field along the geodesic given as follows: The geodesic variation $\Gamma_{g,\xi}(s,t)$ is defined for some $\varepsilon > 0$ as

\[\Gamma_{g,\xi}(s,t):=\exp{\gamma_{x,\xi}(s)}[t\log_{g(s;x,\xi)}y],\qquad s\in(-\varepsilon,\varepsilon),\ t\in[0,1].\]

Intuitively we make a small step $s$ into direction $\xi$ using the geodesic $g(\cdot; x,\xi)$ and from $z=g(s; x,\xi)$ we follow (in $t$) the geodesic $g(\cdot; z,y)$. The corresponding Jacobi field~(J_{g,\xi}) along~(g(\cdot; x,y) is given

\[J_{g,\xi}(t):=\frac{D}{\partial s}\Gamma_{g,\xi}(s,t)\Bigl\rvert_{s=0}\]

which is an ODE and we know the boundary conditions $J_{g,\xi}(0)=\xi$ and $J_{g,\xi}(t) = 0$. In symmetric spaces we can compute the solution, since the system of ODEs decouples, see for example do Carmo, Chapter 4.2. Within Manopt.jl this is implemented as jacobiField(M,x,y,t,ξ[,β]), where the optional parameter (function) β specifies, which Jacobi field we want to evaluate and the one used here is the default.

We can hence evaluate that on the points on the geodesic at

T = [0:0.1:1.0...]

namely

Z = geodesic(M,x,y,T)

the geodesic moves as

ηx = jacobiField.(Ref(M), Ref(x), Ref(y), T, Ref(ξx) )
11-element Array{SnTVector{Float64},1}:
 SnT([0.0, 0.4, 0.5])                    
 SnT([-0.0563164, 0.355568, 0.493844])   
 SnT([-0.0988854, 0.304338, 0.475528])   
 SnT([-0.127117, 0.249482, 0.445503])    
 SnT([-0.141068, 0.194164, 0.404508])    
 SnT([-0.141421, 0.141421, 0.353553])    
 SnT([-0.129443, 0.0940456, 0.293893])   
 SnT([-0.106921, 0.0544789, 0.226995])   
 SnT([-0.0760845, 0.0247214, 0.154508])  
 SnT([-0.0395075, 0.00625738, 0.0782172])
 SnT([0.0, 0.0, 0.0])                    

which can also be called using DxGeo. We can add to the image above by creating extended tangent vectors TVectorE the include their base points

Vx = TVectorE.(ηx,Z)
11-element Array{TVectorE{SnTVector{Float64},SnPoint{Float64}},1}:
 SnT([0.0, 0.4, 0.5])E_Sn([1.0, 0.0, 0.0])                              
 SnT([-0.0563164, 0.355568, 0.493844])E_Sn([0.987688, 0.156434, 0.0])   
 SnT([-0.0988854, 0.304338, 0.475528])E_Sn([0.951057, 0.309017, 0.0])   
 SnT([-0.127117, 0.249482, 0.445503])E_Sn([0.891007, 0.45399, 0.0])     
 SnT([-0.141068, 0.194164, 0.404508])E_Sn([0.809017, 0.587785, 0.0])    
 SnT([-0.141421, 0.141421, 0.353553])E_Sn([0.707107, 0.707107, 0.0])    
 SnT([-0.129443, 0.0940456, 0.293893])E_Sn([0.587785, 0.809017, 0.0])   
 SnT([-0.106921, 0.0544789, 0.226995])E_Sn([0.45399, 0.891007, 0.0])    
 SnT([-0.0760845, 0.0247214, 0.154508])E_Sn([0.309017, 0.951057, 0.0])  
 SnT([-0.0395075, 0.00625738, 0.0782172])E_Sn([0.156434, 0.987688, 0.0])
 SnT([0.0, 0.0, 0.0])E_Sn([6.12323e-17, 1.0, 0.0])                      

and add that as one further set to the Asymptote export.

renderAsymptote("jacobiGeodesicDxGeo.asy",asyExportS2Signals;
    render = asyResolution,
    curves=[geodesicCurve], points = [ [x,y], Z], tVectors = [Vx],
    colors=Dict(
        :curves => [black],
        :points => [TolVibrantOrange,TolVibrantCyan],
        :tvectors => [TolVibrantCyan]
    ),
    dotSizes = [3.5,2.], lineWidth = 0.75, cameraPosition = (1.,1.,.5)
)

A Jacobi field for \$D_xg(t,x,y)[\\eta]\$

If we further move the end point, too, we can derive that Differential in direction

ξy = SnTVector([0.2,0.,-0.5])
ηy = DyGeo.(Ref(M),Ref(x),Ref(y),T,Ref(ξy))
Vy = TVectorE.(ηy,Z)
11-element Array{TVectorE{SnTVector{Float64},SnPoint{Float64}},1}:
 SnT([0.0, -0.0, 0.0])E_Sn([1.0, 0.0, 0.0])                              
 SnT([0.00312869, -0.0197538, -0.0782172])E_Sn([0.987688, 0.156434, 0.0])
 SnT([0.0123607, -0.0380423, -0.154508])E_Sn([0.951057, 0.309017, 0.0])  
 SnT([0.0272394, -0.0534604, -0.226995])E_Sn([0.891007, 0.45399, 0.0])   
 SnT([0.0470228, -0.0647214, -0.293893])E_Sn([0.809017, 0.587785, 0.0])  
 SnT([0.0707107, -0.0707107, -0.353553])E_Sn([0.707107, 0.707107, 0.0])  
 SnT([0.097082, -0.0705342, -0.404508])E_Sn([0.587785, 0.809017, 0.0])   
 SnT([0.124741, -0.0635587, -0.445503])E_Sn([0.45399, 0.891007, 0.0])    
 SnT([0.152169, -0.0494427, -0.475528])E_Sn([0.309017, 0.951057, 0.0])   
 SnT([0.177784, -0.0281582, -0.493844])E_Sn([0.156434, 0.987688, 0.0])   
 SnT([0.2, 0.0, -0.5])E_Sn([6.12323e-17, 1.0, 0.0])                      

and we can look at the total effect, where the TVectorEs even verify that only tangent vectors are added that have a common base point

Vb = Vx .+ Vy
11-element Array{TVectorE{SnTVector{Float64},SnPoint{Float64}},1}:
 SnT([0.0, 0.4, 0.5])E_Sn([1.0, 0.0, 0.0])                             
 SnT([-0.0531877, 0.335814, 0.415627])E_Sn([0.987688, 0.156434, 0.0])  
 SnT([-0.0865248, 0.266296, 0.32102])E_Sn([0.951057, 0.309017, 0.0])   
 SnT([-0.0998779, 0.196021, 0.218508])E_Sn([0.891007, 0.45399, 0.0])   
 SnT([-0.0940456, 0.129443, 0.110616])E_Sn([0.809017, 0.587785, 0.0])  
 SnT([-0.0707107, 0.0707107, 0.0])E_Sn([0.707107, 0.707107, 0.0])      
 SnT([-0.0323607, 0.0235114, -0.110616])E_Sn([0.587785, 0.809017, 0.0])
 SnT([0.0178201, -0.00907981, -0.218508])E_Sn([0.45399, 0.891007, 0.0])
 SnT([0.0760845, -0.0247214, -0.32102])E_Sn([0.309017, 0.951057, 0.0]) 
 SnT([0.138276, -0.0219008, -0.415627])E_Sn([0.156434, 0.987688, 0.0]) 
 SnT([0.2, 0.0, -0.5])E_Sn([6.12323e-17, 1.0, 0.0])                    
renderAsymptote("jacobiGeodesicResult.asy",asyExportS2Signals;
   render = asyResolution,
   curves=[geodesicCurve], points = [ [x,y], Z], tVectors = [Vx,Vy,Vb],
   colors=Dict(
       :curves => [black],
       :points => [TolVibrantOrange,TolVibrantCyan],
       :tvectors => [TolVibrantCyan,TolVibrantCyan,TolVibrantTeal]
  ),
  dotSizes = [3.5,2.], lineWidth = 0.75, cameraPosition = (1.,1.,0.)
)

A Jacobi field for the effect of two differentials (blue) in sum (teal)

Literature