Möbius Geometry of Spheres

by

Rolf Sulanke

March 1, 1999, last revision November 3, 2000

Summary

This notebook treates 2-spheres within the 3-sphere,or the euclidean 3-space, as objects of Möbius Geometry. First we introduce some basic objects of n-dimensional Möbius geometry, which is the conformal geometry of the n-sphere. We emphasize the case n=3, in which Möbius geometry can be visualized in the Euclidean 3-space by stereographical projection from the north pole. For this purpose we construct a version of the stereographical projection,and its inversion,which relates isotropic vectors of the pseudo-euclidean 5-space and points of the 3-space.The bijectivity between subspheres of the 3-sphere and one-dimensional euclidean subspaces of the 5-dimensional pseudo-euclidean vector space is established. Functions describing this bijective relation are constructed. Finally we find and visualize the geodesics of the sphere space.
The concepts developed in this notebook are collected in the package mspher.m.
This notebook continues the notebook eusphere.nb, which contains the metric geometry of spheres. It uses tools of pseudo-euclidean linear algebra, developed in the notebook pseuklid.nb.
For the concepts of Möbius geometry applied here see the preprint [3] and my homepage [4].

References

[1] Coxeter, H. S. M.  Inversive distance. Ann. Mat. Pura Apli. 4, 71 (1966), 73 - 83. 
[2] Kuska, Jens-Peer. http://phong.informatik.uni-leipzig.de/ ~kuska/mview3d.html.

[3] Sulanke,  R.  Möbius invariants for pairs of spheres ((Subscript[S, 1]^m), (Subscript[S, 2]^l)) in the Möbius space (S^n), Beiträge zur Algebra und Geometrie, 41 (2000), No.1, 233-246. Preprint No. 373, SFB 288, Berlin 1999.  http://xxx.lanl.gov/abs/math.MG/9902131.

[4] Sulanke,  R.  Homepage:  http://www-irm.mathematik.hu-berlin.de/~sulanke/

Keywords

Generalized angles between spheres. Conformal invariant. Inversive distance.  Spheres defined by spacelike vectors. Random spacelike unit vectors. Random spheres. Spacelike vectors corresponding to spheres and planes.  Spheres through four points. Stereographical projection of isotropic vectors.  Space of all spheres. Geodesics in the sphere space; spacelike, timelike, and isotropic geodesics.

Address

sulanke@mathematik.hu-berlin.de

The author will be very thankful for hints, corrections, comments.

Initialization

The needed packages

If you want tu use all the notebooks of the item "Spheres", ensure that the packages
euvec.m
pseuvec.m
mspher.m
mcirc.m
spiralsf.m
and the master package init.m are in your working directory. These packages are contained in the item "Spheres".
 For the present notebook we only need the packages euvec.m, pseuvec.m, mspher.m.

3D graphics, animations

In the notebook command lines using ThreeScript, or MVShow3D appear. They may used for 3D-representations and animations of Mathematica Graphics3D objects. To apply them you must have the necessary software on your computer. Look in the last section.

For Windows

For the commands in this subsection:
DO NOT CONFIRM THE EVALUATION OF THE INITIALIZATION CELLS !!!

If the result of the following command

Directory[]

does not coincide with the working directory, e.g. "C:\Eigene Dateien\mymath\spheres.v4", containing the needed packages, insert the Path of your working directory in the next command line and activate it first, without calling the initialization:

SetDirectory["C:\\Eigene Dateien\\mymath\\spheres"]

If you want to use the package euvec.mpseuvec.m,mspher.m only,  carry out the next command.

Needs["mspher`"]

If you want  that all commands contained in the packages of the item "Spheres" are available, activate

Get["init.m"]

Now: Activite Kernel/Evaluation/Evaluate Initialization from the menu!

For Linux (and other  systems: not tested!)

For the next command:
DO NOT CONFIRM THE EVALUATION OF THE INITIALIZATION CELLS !!!

If you want to use the package euvec.mpseuvec.m,mspher.m only, remove the package init.m from the working directory and carry out the next command.

Needs["mspher`"]

Now: Activite Kernel/Evaluation/Evaluate Initialization from the menu!

1  Introduction

1.1 Pseudo-euclidean vector algebra and other prerequisites

We consider the n-dimensional sphere, shortly n-sphere, as the standard hyperellipsoid in the (n+1)-dimensional projective space. In homogeneous coordinates x ↦ (x[1], ..., x[n+2])  its equation is

(∑_(i = 1)^(n + 1))(((x[i]))^2) - (((x[n + 2]))^2) = 
      0

By Liouville's theorem (and some elementary consideration about functions of one complexe variable in the case n=2) the conformal group of the n-sphere (n>1) coincides with the group of  all projective transformations preserving this hyperellipsoid. Elementary linear algebra can be used to show, that this group can be represented as the pseudo-orthogonal group  (factored by ∓ identity) of the (n+2)-dimensional real vector space underlying the (n+1)-dimensional projective space. Therefore Möbius geometry is treated by tools of pseudo-euclidean linear algebra. The basic concepts of pseudo-euclidean linear algebra are developed in the notebook pseuklid.nb and contained in the package pseuvec.m. They are introduced automatically as part of the master package init.m of this directory.  Before importing and using them one has to fix the dimension dim of the underlying vector space and the index  ind of the (non-degenerate, symmetric) scalar product. These basic constants are now defined by

dim  =  5;  ind  =  1;

because in this notebook we are mainly interested in the 3-dimensional sphere. You can change these constants only before starting the evaluation of the initialization cells  by changing the last command line. But remember that these constants enter in almost all of the treated concepts; therefore they should be defined before using the packages.

The basic invariant of pseudo-eulidean linear algebra is the scalar product:

?pssp

The matrix of the scalar products in any orthonormed basis is

MatrixForm[Table[pssp[stb[i],  stb[j]], {i, 5}, {j, 5}]]

Terminology: Spheres are usual 2-dimenional spheres if we consider the 3-dimensional spherical geometry; more generally one may consider m-spheres in the n-sphere, where m, n with 0≤m≤n denote the dimensions. For  m = n-1 the (n-1)-spheres are called  hyperspheres.

In 1.2 you can find the new symbols and evoke the corresponding usages.

1.2 The usages

After a package has been loaded you can get information about the  listed concepts in the usual way, e.g.

?orthonorm

For getting the usage and the complete definition, enter e.g.

?? spt

Since the names are prepended by the ContextNames the definition obtained in this way are hard to read. Therefore you will find some of them in the original form in unevaluable cells.  It is recommended to look in the packages.

Here you can get information about the commands which are defined in mspher.m. For the functions and mofules defined in the packages euvec.m, pseuvec.m  see the notebooks euspher.nb, pseuklid.nb.

?spt

?stepro

?invstepro

?confinv

?gencos

?plotsphere

?euklidsphere

?euklidsphereplot3D

?spunit

?pstg5frame

?vradius

?vcenter

?spacevec

?sphvec

?planev

2  Points: Adaption of the stereographical projection

2.1 Stereographical projection of isotropic vectors

In the projective context, the points of the 3-sphere are the generating lines of the isotropic cone

pssp[x,x]  = x[1]^2 + x[2]^2 + x[3]^2 + x[4]^2 - x[5]^2= 0,

analogously in higher dimensions. Our aim is to describe the point of the euclidean 3-space obtained by stereographical projection sterproj (defined in pseuvec.m, see also eusphere.nb) of the 3-sphere onto its equatorial hyperplane as a function of the isotropic vector x. There exists a natural bijection between the  lines intersecting the hyperplane

- pssp[stb[5], x] = x[5]  = 1 						(H)

of the 5-space and the points of the euclidean 4-space  defined by

?spt

\(spt[v_List]\  := \ 
      Module[{w, u, d = Length[v]}, 
        w\  = v[\([d]\)]; 
			If[
          w === 0, \ \ Return[
            Print[v\ , " has last component ", w]], 
				u = 
            Take[v/w, d - 1]; Return[u]]]\)

Now remember that the  3-sphere appears as  intersection of the hyperplane (H) with the isotropic cone. If we restrict spt onto the isotropic cone and combine it with the general stereographical projection

?sterproj

we obtain the wanted map

?stepro

stepro[v_List]  :=  sterproj[spt[v]]

Its inverse is generated by the vector

invstepro[v_List]  :=  Append[invsterproj[v], 1]

?invstepro

?invsterproj

2.2 Examples

spt[{2, 2, 2, 2, 2, 2, 2, 2, 4}]

vec[x]

spt[vec[x]]

spt[stb[1]+stb[2]]

spt[stb[5]]

spt[stb[5] + stb[2] Cos[t] + stb[3] Sin[t]]

stepro[stb[5] + stb[2] Cos[t] + stb[3] Sin[t]]

invstepro[%]

Simplify[%]

If one applies stepro on non-isotropic vectors one obtains incorrect results without any warning! If there are any doubts one should control the character ch of the arguments! Here is such an example:

vf[s_, t_]  := 
      stb[1] + (stb[2] Cos[t] + stb[3] Sin[t]) Cos[s]  + 
        stb[5] Sin[s]

vf[s, t]

ch[vf[s, t]]

since ch[vf[s,t]]depends on s,t!

pssp[vf[s, t], vf[s, t]]

fct[s_]  =  Simplify[%]

Plot[fct[s],  {s, 0,  Pi}]

spt[vf[s, t]]

stepro[vf[s, t]]

Thus, we obtained a result for all values of t, but it corresponds to  a point of a the sphere only for solutions of Cos[s] == 0.

vf[Pi/2, t]

sterproj[spt[vf[Pi/2,  t]]]

invstepro[%]

spt[nul[]]

invstepro[nul[3]]

Remove[vf,  fct]

2.3 A formula for the scalar products

As an example for symbolic calculation we derive a formula, which expresses the pseudo-euclidean scalar products of the images by investepro of two points in terms of the euclidean scalar products of their radius vectors. The  points are

vx  = vec[x, 3];  vy  =  vec[y, 3];

{vx, vy} // MatrixForm

The preimages are the rows of the matrix

{invstepro[vx],  invstepro[vy]}

The preimages  are isotropic:

Simplify[Together[pssp[%[[1]], %[[1]]]]]

Their pseudo-euclidean scalar product is

pssp[invstepro[vx],  invstepro[vy]]

Simplify[Together[%]]

Simplify[% /. {x[1]^2 + x[2]^2 + x[3]^2 ->  nvx^2, 
          y[1]^2 + y[2]^2 + y[3]^2 -> nvy^2}]

Therefore, the following formula is true:

Simplify[
      pssp[invstepro[vx],  
          invstepro[
            vy]] ==  -2*(vx - vy) . (vx - 
                  vy)/((1 + vx . vx)*(1 + vy . vy))]

Thus we obtained the wanted formula, which symbolically can be written as

pssp[invstepro[x], invstepro[y]]
== | 
      x - y | ^2/((1+ | 
                x | ^2)*(1+ | 
                y | ^2))

Remark.  The isotropic vector

invstepro[{x1, x2, x3}]

is a special representant of the corresponding point

stepro[invstepro[{x1, x2, x3}]]

defined by the condition

pssp[invstepro[{x1, x2, x3}],  stb[5]] ==  -1

Thus, the calculated scalar product does not have any invariant meaning in Möbius geometry. In Möbius geometry all pairs of (distinct) points are equivalent.

3 Spheres

In this section we shall establish the basic bijection between spheres (including planes as spheres of infinite radius) of the euclidean 3-space and one-dimensional euclidean subspaces of the 5-dimensional pseudo-euclidean vector space, see also subsection 3.5.

3.1 Spheres through four points and corresponding spacelike vectors: vradius, vcenter, sphereplot3D

The definitions of the spacelike vector functions sphvec, spacevec

The following Module yields an unnormed spacelike vector which defines the sphere through four points in general position (i.e. do not belong to a circle or a line):

?sphvec

\(\(\(sphvec[p1_List, p2_List, p3_List, 
        p4_List]\)\(:=\)\(Module[{v = nul[5]}, 
        If[Length[p1] ≠ 3 || Length[p2] ≠ 3 || 
            Length[p3] ≠ 3 || Length[p4] ≠ 3, 
          Print["An argument does not have the correct
dimension!"] Return[]]; \ 	v = 
          Simplify[
            pscross[{invstepro[p1], invstepro[p2], invstepro[p3], 
                invstepro[p4]}]]; 
If[v === nul[5], 
          Print["The points are not in general
position!"]; Return[], Return[v], Return[v]]]\)\(	\)\)\)

Applying this procedure to four points lying on the sphere with center{x,y,z} and radius r, we obtain after some simplifications (see below) the formula

?spacevec

vn  = spacevec[x, y, z, r]

pssp[vn, vn]

Simplify[Together[%]]

Examples, applications

We use euklidsphereplot3D to plot the corresponding sphere or plane. For a detailed discussion see subsection 3.3. 
Here we get the standard unit sphere:

sphvec[stb[1, 3],  stb[2, 3],  
      stb[3, 3],  -stb[3, 3]]

euklidsphereplot3D[%]

In Möbius geometry, planes appear as stereographical projections of spheres through the pole; sphvec includes this case:

sphvec[nul[3],  stb[2, 3],  -stb[2, 3],  stb[3, 3]]

euklidsphereplot3D[%]

sphvec[stb[1, 3],  -stb[1, 3],  
      stb[2, 3],  -2 stb[2, 3]]

euklidsphereplot3D[%]

If the four points ly on a circle, one gets a corresponding message:

sphvec[stb[1, 3],  -stb[1, 3],  
      stb[2, 3],  -stb[2, 3]]

{invstepro[stb[1, 3]],  invstepro[-stb[1, 3]], 
      invstepro[ stb[2, 3]],  invstepro[-stb[2, 3]]}

pscross[%]

The spacelike vector  as a funtion of center and radius of the sphere: spacvec

We demonstrate the developing of the procedure spacevec mentioned above. This is the spacelike vector corresponding to the sphere with center {x,y,z}and radius r:

\(sphvec[{x, y, z} + r*stb[1, 3], {x, y, z} - 
        r*stb[1, 3], 
		{x, y, z} + r*stb[2, 3], \ {x, y, z} + 
        r*stb[3, 3]]\)

%*(1 + (r^2) - 2 r x + (x^2) + (y^2) + (z^2)) (1 + (r^2) + 
          2 r x + (x^2) + (y^2) + (z^2)) (1 + (r^2) + (x^2) + 2 r y + 
          (y^2) + (z^2)) (1 + (r^2) + (x^2) + (y^2) + 2 r z + 
          (z^2))

%/(16  r^3)

Since the norming is unessential, we  defined in the package mspher.m

spacevec[x_, y_, z_, 
        r_]  :=  {2 x, 2 y, 2 z, -1 - (r^2) + (x^2) + (y^2) + (z^2),  
          1 - (r^2) + (x^2) + (y^2) + (z^2)}/(2*r)

\(Simplify[
      pssp[{2\ x, 2\ y, 2\ z, \(-1\) - r\^2 + x\^2 + y\^2 + z\^2, 
          1 - r\^2 + x\^2 + y\^2 + z\^2}, \ 
		{2\ x, 2\ y, 
          2\ z, \(-1\) - r\^2 + x\^2 + y\^2 + z\^2, 
          1 - r\^2 + x\^2 + y\^2 + z\^2}]]\)

spacevec[x, y, z, r]

Simplify[pssp[spacevec[x, y, z, r],  spacevec[x, y, z, r]]]

euklidsphereplot3D[spacevec[0, 0, 10, 2]]

Here one gets error messages:

spacevec[0, 0, 1, 0]

The following yields a spacelike unit vector of the sphere with radius r and center {x,0,0}:

spherevec[x_, r_] :=  spacevec[x, 0, 0, r]

spherevec[x,  r]

euklidsphere[spherevec[x, r]][u, v]

% /. {x  ->  3, r ->  2}

ParametricPlot3D[%,  {u,  -1,  1}, {v,  -1,  
        1}]

spacevec[0, 0, 0, 1]

spacevec[0, 0, 0, -r]

pssp[%,  %]

Simplify[%]

Inversion of spacevec, and an euclidean definition of the sphere corresponding to a spacelike vector

From the formula

vn  = spacevec[x, y, z, r]

we conclude

Simplify[Together[vn[[5]] -  vn[[4]]]]

LEMMA. The sphere corresponding to the spacelike vector vn is a plane if and only if one has vn[[4]] = vn[[5]].
Proof: The relation between spheres (including planes) and one-dimensional euclidean subspaces of the 5-dimensional pseudo-euclidean space is bijective. Since the north pole of the 3-sphere corresponds to the isotropic vector stb[4]+stb[5], our condition is equivalent to

pssp[stb[4]  +  stb[5],  vec[v]] ==  0

COROLLARY. The euclidean sphere which corresponds to the spacelike unit  vector

vec[v]

with v[4] ≠ v[5]  is the sphere with center {x, y, z}

x  =  v[1]/(v[5] - v[4]); y  =  v[2]/(v[5] - v[4]); 
    z  =  v[3]/(v[5] - v[4]);

and radius

r  =  1/(v[5] - v[4]);

Remember that vn is normed: pssp[vn,vn] ==1. Thus we conclude for radius and center of the sphere corresponding to the spacelike unit vector v:

vradius[wx_List]  :=  
      Module[{s, q,  vx}, 
q  =  
            pssp[wx, wx]; 
If[q ≤  0, 
            Print[wx,  " is not spacelike!"]; 
            Return[]]; 
vx  =  
            wx/Sqrt[q]; 
s = 
            vx[[5]] -  vx[[4]]; 
If[
            s == 0,  Return[Infinity], Return[1/s], Return[1/s]]] /; 
        Length[wx]  ==  5

vcenter[vx_List]  :=  
      Module[{s}, 
If[pssp[vx, vx] ≤  0, 
            Print[vx,  " is not spacelike!"]; 
            Return[]]; 
s = 
            vx[[5]] -  vx[[4]]; 
If[
            s == 0,  Return[{Infinity, Infinity, Infinity}]]; 
          Return[{vx[[1]]/s, vx[[2]]/s, vx[[3]]/s}]] /; 
        Length[vx] == 5

WARNINGS!
1. These definitions make sense for spacelike vectors only!

vcenter[{1, 0, 0, 1, 3}]

vradius[nul[]]

2. If the corresponding sphere is a plane  center and radius will become infinite,  e.g.

vcenter[stb[1]]

vcenter[{1, 0, 0, 1, 1}]

This case is treated in the next subsection.

Orientation

A negative radius yields the negatively oriented sphere.

spacevec[x, y, z, r]

spacevec[x, y, z, -r]

For v[[4]] > v[[5]] the radius of the sphere becomes negative. Changing the sign of v one may always reach r > 0; the center is independent of the sign. The sign of the radius can be interpreted as the orientation of the sphere. If one considers the sphere corresponding to v as subsphere of (S^3) (that means: not carry out the stereographic projection), the sphere divides the 3-sphere into two disks: the inner and the outer disk. Let us define the inner disks as the set of all points (= isotropic 1-spaces) w of the 3-sphere  (S^3) with pssp[v,w] < 0. Then changing the sign of v means permuting inner and outer disk. It follows that the spacelike unit vectors correspond bijectively to the oriented subspheres of  (S^3).

The sphere corresponding to a spacelike unit vector

Applying the standard euclidean parameter representation of a sphere ( see the notebook eusphere.nb) we get immediately the sphere corresponding to the spacelike unit vector  v:

?sphere

sphere[vcenter[vn], vradius[vn]][u, v]

Simplify[Together[%]]

Some examples

vradius[{1, 0, 0, 2, 1}]

vcenter[{1, 0, 0, 1, 3}]

vradius[stb[4]*5]

vcenter[stb[4]*5]

?sphereplot3D

This is the standard unit sphere with center in the origin:

sphereplot3D[vcenter[stb[4]],  vradius[stb[4]]]

Table[
      sphereplot3D[vcenter[{0, 0, 0, Sqrt[1 + j^2],  j}],  
        vradius[{0, 0, 0, Sqrt[1 + j^2],  j}],  Pi,  2*Pi, -Pi/2, 
        Pi/2], {j, 0, 4}]

Show[%]

3.2 Planes

Before constructing  modules which yield the sphere or plane corresponding to a spacelike vector v we find the spacelike vector corresponding to a plane of the euclidean 3-space.

The definition planev

?planev

The following Module yields an unnormed spacelike vector which defines the plane, fulfilling the equation w.x==p:

\(planevec[w_List, \ p_]\  := \ 
      Module[{v1, v2, v3, u4\  = stb[4] + stb[5], 
		u1, u2, u3}, 
		If[
          Length[w] != 3 || \ w\  === nul[3], 
			Print[
            w, \ " is not a normal vector!"]; Return[]]; 
		v1 = 
          w*p/w . w; 
		If[w[\([1]\)] === 0, 
          v2\  = v1 + stb[1, 3]; 
			If[w[\([2]\)] === \ 0, \ 
            v3\  = \ v1\  + \ stb[2, 3], 
            v3\  = {0, \ p/w[\([2]\)], 
                0}], 
			v2\  = \ {p/w[\([1]\)]\  - \ 
                w[\([2]\)]/w[\([1]\)], 1, 
              0}; 
			v3\  = \ {p/w[\([1]\)]\  - \ w[\([3]\)]/w[\([1]\)], 
              0, 1}]; 
		u1 = \ invstepro[v1]; 
		u2 = \ 
          invstepro[v2]; 
		u3 = \ invstepro[v3]; 

        pscross[{u1, u2, u3, u4}]]\)

The general evaluation with succeeding simplification

Simplify[planevec[{x, y, z}, h]]

\(%*\((\((h^2\  + \ 2*x^2\  - \ 2*h*y\  + \ y^2)\)*\((h^2\  + \ x^2\  + \ 
              y^2\  + \ z^2)\)*
\ \ \ \ \ \ \ \((h^2\  + \ 2*x^2\  - \ 
              2*h*z\  + \ z^2)\))\)\)

%/((8*x^3)*(-x^2 - y^2 - z^2 + 
              h*(y + z)))

Simplify[%]

finally leads to the following definition, contained in the package mspher.m:

\(planev[a_, b_, c_, \ p_]\  := \ 
      If[Not[a^2 + b^2 + c^2 === 0], 
	Return[{a, \ b, \ c, \ p, \ p}], 
        Print["Normal vector null is not admissible!"]; 
        Return[]]\)

planev[a, b, c, p]

planev[0, 0, 0, p]

Clearly, if we norm planev[a,b,c,p] by

dotnorme[{a, b, c}]

{a, b, c} . {a, b, c} ==  1

or, what means the same,

normalize[{a, b, c, p, p}]

pssp[{a, b, c, p, p}, {a, b, c, p, p}]  == 1

we have a bijection between the manifold of all oriented planes of the euclidean 3-space and the manifold of all spacelike unit vectors v with the property v[[4]]==v[[5]]. Usually one names an equation n.x == p with normed 3-vector n the Hesse normal form of the plane equation. The following module yields a parameter representation of the plane with this normal form, which we translate into the 5-vector vc = {a,b,c,p,p}:

plane[vc_List][u_, v_]  :=  
      Module[{q, p, w, ob}, 
q = 
          pssp[vc, vc]; 
If[
          q ≤ 0 ||  vc[[4]] ≠  vc[[5]], 
          Print[vec,  " incorrect argument!"]; 
          Return[]];  
p  =  
          vc[[5]]/Sqrt[q]; 
w  =  
          Take[vc, 3]; 
ob  =  
          orthobase[w]; 
Return[
          ob[[3]]*p  +  ob[[1]]*u  +  ob[[2]]*v]]

Examples

nx = planev[1, 1, 1, 1]

pssp[nx, nx]

plane[nx][u, v]

ParametricPlot3D[%,  {u,  -3, 3},  {v,  -3, 3}]

nx  = planev[0, 0, 1, 0]

plane[%][u, v]

ParametricPlot3D[%,  {u,  -3, 3},  {v,  -3, 3}]

Clear[nx]

The plane through three points

Through three points in general position there exists exactly one plan in the euclidean 3-space. The following Module gives a spacelike 5-vector corresponding to this plane:

spvplane3pts[p1_List, p2_List, p3_List] :=  
      Module[{q1, q2, q3, vq,  wrong}, 
If[
          Not[Length[p1] == 3 && Length[p2] == 3 && 
              Length[p3] == 3],  
          Goto[wrong]]; 
q1  =  invstepro[p1]; 
        q2 =  invstepro[p2]; 
        q3  =  invstepro[p3]; 
vq  =  
          pscross[{q1, q2, q3, stb[4] + stb[5]}]; 
If[
          vq  ==  nul[5],  Goto[wrong], Return[vq], 
          Return[vq]]; 
Label[
          wrong]; 
Print["Incorrect arguments, or points \
not in general position!."]; 
Return[]]

spvplane3pts[nul[3],  stb[1, 3],  -stb[1, 3]]

spvplane3pts[nul[3],  stb[1, 3],  stb[2, 3]]

plane[%][u, v]

ParametricPlot3D[%,  {u,  -3, 3},  {v,  -3, 3}]

spvplane3pts[stb[3, 3],  stb[1, 3],  stb[2, 3]]

plane[%][u, v]

ParametricPlot3D[%,  {u,  -3, 3},  {v,  -3, 3}]

This is a random point in the unit cube of the euclidean 3-space:

randomv[3]

A random plane as the plane through three random points:

spvplane3pts[randomv[3], randomv[3], randomv[3]]

plane[%][u, v]

ParametricPlot3D[%,  {u,  -3, 3},  {v,  -3, 3}]

3.3 The sphere or plane corresponding to a spacelike vector

In this section we join the results of the two foregoing sections: We construct a procedure, which gives the parameter representation of the sphere or plane corresponding to a spacelike vector, and apply it in a plot command.

The definition of the parameter representation euklidsphere

euklidsphere[va_List][u_, v_]  :=  
      Module[{q,}, 
q  =  
            pssp[va, va]; 
If[q ≤  0,  
            Print[va,  " is not spacelike!"];  
            Return[]]; 
If[
            va[[4]]  ==  va[[5]],  Goto[pla]]; 
          Return[sphere[vcenter[va], vradius[va]][u*Pi, 
              v*Pi/2]]; 
Label[pla]; 
 
          p  =  va[[5]]/Sqrt[q]; 
w  =  
            Take[va, 3]; 
ob  =  
            orthobase[w]; 
Return[
            ob[[3]]*p  +  ob[[1]]*u  +  
              ob[[2]]*v]
] /; 
        Length[va] ==  5

Examples

euklidsphere[stb[5]][u, v]

euklidsphere[stb[4]][u, v]

euklidsphere[stb[1]][u, v]

euklidsphere[stb[2]][u, v]

euklidsphere[stb[3]][u, v]

Evaluate repeated the following cell; "almost surely" you will get a random halfsphere.

vrand  =  randomv[5]; 
    ParametricPlot3D[euklidsphere[vrand][u, v],  {u, -1, 1}, {v,  
        0,  1}]

Clearly, the probability for vrand[4] = vrand[5] is zero, therefore planes do happen with probability zero. Furthermore remember, that the probability for a random vector in the 5-dimensional unit cube to be not spacelike, is small; therefore the probability of error messages is small. Indeed, as already shown in the notebook pseuklid.nb, we have for the empirical probability nearly 6.3 %.

( nr  =  0;  
      For[i  =  1,  i ≤  10000, i++, 
        If[ch[randomv[5]] <  1,  nr  =  nr + 1]];  
      nr/10000)

N[%]

Remove[nr, vrand]

euklidsphereplot3D and plotsphere

Now it is easy to build the plot command for euklidsphere:

euklidsphereplot3D[vn_List, u0_: - 1, u1_:  1,  v0_: - 1, v1_:  1, 
        opts___] := 
      ParametricPlot3D[
        Evaluate[euklidsphere[vn][u, v]], {u, u0, u1},  {v, v0, v1}, 
        opts]

The package mspher.m contains a similar procedure plotsphere:

?plotsphere

Now we test euklidsphereplot3D  and compare it with plotsphere.

3.4 Examples, comparision

Example: equidistant unit spheres

Timing[
      Table[plotsphere[spacevec[x, 0, 0, 1]],  {x,  0,  12,  
          3}]]

Show[%[[2]]]

Timing[
      Table[euklidsphereplot3D[spacevec[x, 0, 0, 1]],  {x,  0,  12,  
          3}]]

Show[%[[2]]]

grser  = %

ThreeScript["spherser.ts", grser]

With plotsphere 3.41 s, and with euklidsphereplot3D 1.15 s have been needed. There is yet another difference: with euklidsphereplot3D the geographical coordinates are defined with respect to the fixed orthonormed standard frame, and with plotsphere the coordinates depend on the considered sphere. Furthermore, since in the definition of plotsphere we used the stereographical projection, and in euklidsphereplot3D not, the coordinates and the representations generated by plotsphere are often distorted.

Example: The spacelike frame vectors

Timing[Table[plotsphere[stb[i]],  {i, 4}]]

Show[%[[2]]]

Timing[Table[euklidsphereplot3D[stb[i]],  {i, 4}]]

Show[%[[2]]]

?spherevec

Timing[
      Table[plotsphere[spherevec[2  t, t/2]],  {t,  0.5,  9.5,  
          3}]]

Show[%[[2]],  PlotRange  ->  All]

grser2 = %

ThreeScript["spherse2.ts", grser2]

Timing[
      Table[euklidsphereplot3D[spherevec[2  t, t/2]],  {t,  0.5,  9.5,  
          3}]]

Show[%[[2]]]

We define a function which constructs a certain monotone sequence from a given sequence of positive numbers:

ms[rs_List]  :=  
      Table[2*Sum[rs[[j]], {j, 1,  i - 1}] +  rs[[i]],  {i, 1, 
          Length[rs]}]

Here is a random  sequence of positive real numbers:

rd  =  randomv[]

Now we obtain a sequence of spheres with random radii rd less than 1 and centers on the x-axis in such a way that to consecutive spheres are tangent to each other; the sequence ms[rd] defines the centers of the spheres.

ms[rd]

Timing[
      Table[plotsphere[spherevec[ms[rd][[i]], rd[[i]]]],  {i, 1, 
          5}]]

Show[%[[2]]]

Timing[
      Table[euklidsphereplot3D[
          spherevec[ms[rd][[i]], rd[[i]]]],  {i, 1, 
          5}]]

Show[%[[2]]]

Timing[plotsphere[stb[1]*Cos[Pi/5]  +  stb[4]*Sin[Pi/5]]]

Timing[
      euklidsphereplot3D[stb[1]*Cos[Pi/5]  +  stb[4]*Sin[Pi/5]]]

\(Timing[
      Table[plotsphere[
          stb[1]*Cos[i*Pi/5]\  + \ stb[4]*Sin[i*Pi/5]], 
	{i, \ 1, \ 
          5}]]\)

Show[%[[2]]]

\(Timing[
      Table[euklidsphereplot3D[
          stb[1]*Cos[i*Pi/5]\  + \ stb[4]*Sin[i*Pi/5]], 
	{i, \ 1, \ 
          5}]]\)

Show[%[[2]]]

Here we obtain some exemplars of the bundle of planes containing the y-axis:

\(Timing[
      Table[plotsphere[
          stb[1]*Cos[i*Pi/5]\  + \ stb[3]*Sin[i*Pi/5]], 
	{i, \ 1, \ 
          5}]]\)

Show[%[[2]]]

\(Timing[
      Table[euklidsphereplot3D[
          stb[1]*Cos[i*Pi/5]\  + \ stb[3]*Sin[i*Pi/5]], 
	{i, \ 1, \ 
          5}]]\)

Show[%[[2]]]

euklidsphereplot3D is essentially faster than plotsphere, and gives nicer pictures!

Evaluation problems: Floating point numbers

WARNING!
By numerical reasons, the results may be distorted, or completely false. E.g. we met the following situation

vd[t_]  :=  stb[1]*Cosh[t]  + stb[5]*Sinh[t]

euklidsphere[vd[t]][u, v]

Timing[
      Table[euklidsphereplot3D[vd[t*(2*Pi)]],  {t, 0, 4, 
          0.5}]]

On the other hand, Mathematica may calculate and Plot the last two spheres of the sequence, if we calculate with higher precision::

euklidsphereplot3D[N[{Cosh[7*Pi], 0, 0, 0, Sinh[7*Pi]}, 20]]

N[{Cosh[8*Pi], 0, 0, 0, Sinh[8*Pi]}, 30]

euklidsphereplot3D[N[{Cosh[8*Pi], 0, 0, 0, Sinh[8*Pi]}, 30]]

The reason is that Mathematica, by the 0.5 as the running index step, applied floating point numbers in the calculations of the sphere sequence! With lower precision, the vectors seem to be isotropic; they are very near to the asymptotic (isotreopic) cone. With  .5 replaced by 1/2 the same happens:

Timing[
      Table[euklidsphereplot3D[vd[t*(2*Pi)]],  {t, 0, 4, 
          1/2}]]

Besides of this, Mathematica does not know the simplest trigonometric rules:

Cosh[6*Pi]^2 - Sinh[6*Pi]^2

N[%]

TrueQ[Cosh[6*Pi]^2 - Sinh[6*Pi]^2 == 1]

TrueQ[Cosh[x]^2 - Sinh[x]^2 == 1]

TrueQ[Cos[x]^2 +  Sin[x]^2 == 1]

Using Simplify the equations become True:

TrueQ[Simplify[Cosh[x]^2 - Sinh[x]^2] == 1]

TrueQ[Simplify[Cos[x]^2 +  Sin[x]^2] == 1]

Conclusion: It is often better introducing more Simplify's in the procedures.
But then it may happen that the evaluation needs much more time.

3.5 The manifold of all spheres, random spheres

The sphere manifold

Any (non oriented) sphere is represented by a pseudo-euclidean subspace of dimension 4; the manifold of all spheres is the open orbit of pseudo-euclidean 4-spaces under the action of the Möbius group.  Going to the orthogonal complements one finds that the spacelike unit-hypersphere, the hyper-hyperboloid

pssp[x,x]=x[1]^2 + x[2]^2 + x[3]^2 + x[4]^2 - x[5]^2= 1		(H)

is a twofold covering of the manifold of all non-oriented spheres. We immediately can write down a parameter representation of this hypersurface calculating

Append[unit[3][{s,  t,  u}],  0]*Cosh[v]  +  
      stb[5]*Sinh[v]

?unit

Context[unit]

?spunit

spunit[s, t, u,  v]

Here we used unit[3][{s,t,u}], the parameter representation of the unit 3-sphere in the 4-dimensional euclidean subspace orthogonal to stb[5], see euvec.m, or eusphere.nb. Obviously, spunit can be generalized easily to arbitrary dimension.

Now we may apply  the euklidsphere functions of  section 3.3 on this vector and obtain a parameter representation or plot of the corresponding sphere. In an old version of this notebook we started with the now deleted function psphere:

?vsphere

\(\(\(	\)\(Timing[\(psphere[spunit[s, t, u, v]]\)[x, y]]\)\)\)

The needed time was 1483.02 s (on a very fast PC).

Mathematica were working a long time calculating an odd looking formula.  Also Simplify did not help very much. FullSimplify I aborted after more than one hour.  With psphere replaced by euklidsphere  we get in less than 5 s a very nice formula:

\(\(\(	\)\(Timing[
      Simplify[\(euklidsphere[spunit[s, t, u, v]]\)[x, 
          y], 
s ∈ Reals && t ∈ Reals && 
          u ∈ Reals && v ∈ Reals && x ∈ Reals && 
          y ∈ Reals]]\)\)\)

A bad definition

Therefore, and by certain evaluation advantages, we would prefer the direct definition

\(\(euklidspherefct[s_, t_, u_, v_]\)[x_, 
        y_] := 
{\(Cos[π\ x]\ Cos[\(π\ y\)\/2] + Cos[s]\ Cos[t]\ Cos[u]\ Cosh[v]\)\/\(\(-Cosh[v]\)\ Sin[u] + Sinh[v]\), \(Cos[t]\ Cos[u]\ Cosh[v]\ Sin[s] + Cos[\(π\ y\)\/2]\ Sin[π\ x]\)\/\(\(-Cosh[v]\)\ Sin[u] + Sinh[v]\), \(Cos[u]\ Cosh[v]\ Sin[t] + Sin[\(π\ y\)\/2]\)\/\(\(-Cosh[v]\)\ Sin[u] + Sinh[v]\)}\)

Obviously, we then get infinite vectors for u==v==0:

euklidspherefct[s, 0, 0, 0][x, y]

On the other hand, we get the correct result

euklidsphere[spunit[s, 0, 0, 0]][x, y]

Something must be wrong!

euklidsphere[spunit[s, t, u, v]][x, y]

% /. u -> 0

% /. v -> 0

The reason must be that under the general calculation of euklidsphere[spunit[s,t,u,v]][x,y] an error appears.
Indeed, looking at the definition of  euklidsphere in mspher.m one sees that for

spunit[s, t, u, v]

the procedure, of course, en general would refuse

%[[4]]  == %[[5]]

Therefore vectors as

spunit[s, t, 0, 0]

are not recognized as those which  generate planes. The attempt to evaluate the command line Return[xtest]

va  = spunit[s, t, 0, 0];  xtest  =  sphere[vcenter[va], vradius[va]][u*Pi, 
            v*Pi/2]

leads to infinite vectors what is reflected in the calculation of the general formula. Thus, it does not make sense to evaluate for values with u = v = 0, which correspond to planes, the function

euklidspherefct[s, t, u, v][x, y]

A better definition

Since the derived formula is not general enough, we set

Clear[euklidspherefct]

euklidspherefct[s_, t_, u_, v_][x_, 
        y_] :=  euklidsphere[spunit[s, t, u, v]][x, y]

and obtain the correct result:

euklidspherefct[s, 0, 0, 0][x, y]

euklidsphere[spunit[s, 0, 0, 0]][x, y]

Timing[euklidspherefct[s, 0, 0, 0][x, y]]

spunit[Pi/4, Pi,  Pi/4, 0]

euklidsphereplot3D[spunit[Pi/4, Pi,  Pi/4, 0]]

euklidspherefct[Pi/4, Pi,  Pi/4, 0][x, y]

Since

pssp[spunit[2*Pi*Random[],  Pi*Random[]/2,  0, 0], 
        stb[4] + stb[5]] ==  0

Commands like the following always create a set of planes through the origin:

Table[
      euklidsphereplot3D[
        spunit[2*Pi*Random[],  Pi*Random[]/2,  0, 0]], {20}]

Show[%]

A random spacelike unit vector

sprandomv  :=  
      Module[{re,  rs,  rt, rv,  i}, 
re[i_]  :=  
          Sign[0.5  -  Random[]]; 
rs[]  :=  
          Random[]*2*Pi; 
rt[i_]  :=  
          Random[]*Pi/2; 
rs[]  :=  
          re[3]*Random[]*10; 
spunit[rs[], re[1]*rt[1],  
          re[2]*rt[2],  rs[]] ]

WARNING: The distribution of the random spacelike unit vector and the corresponding random sphere is not invariant under the action of the Möbius group! The corresponding random sphere is:

vte  =  sprandomv

pssp[vte, vte]

Timing[euklidsphere[vte][x,  y]]

Timing[euklidsphereplot3D[vte]]

vte  =  Table[sprandomv, {j, 10}]

Table[pssp[%[[j]], %[[j]]],  {j, 10}]

Timing[
      Table[euklidsphere[vte[[j]]][x,  y], {j, 10}]]

vte  =  Table[sprandomv, {j, 8}]

Table[
      euklidsphereplot3D[vte[[j]],  0, 1,  -1, 1], {j, 
        8}]

Show[%,  PlotRange  ->  All]

Since the random radii and centers are wide spread  the union of the halfspheres give not a nice picture

Table[{vradius[vte[[j]]],  vcenter[vte[[j]]]}, {j, 
        8}]

Oriented spheres

Already in section 3.1 we mentioned oriented spheres. Let v be a spacelike vector. The sphere S(v) defined by v corresponds to the isotropic vectors x solving pssp[v,x] == 0. This defines a linear pseudo-euclidean 4-space dividing  the isotropic cone into two parts. On the 3-sphere these parts are the "inner ball", defined by pssp[v, y] < 0, and the outer ball, defined by pssp[v,y] > 0. Clearly, if we replace v by -v, the corresponding sphere does not change, but inner and outer ball exchange there roles. Therefore the hyper-hyperboloid (H) represents the manifold of all oriented spheres contained in  the 3-sphere bijectively.

vte  =  sprandomv

pssp[vte, vte]

Timing[euklidsphere[vte][x,  y]]

Timing[euklidsphereplot3D[vte]]

Timing[euklidsphere[-vte][x,  y]]

Timing[euklidsphereplot3D[-vte]]

vte

vradius[vte] == -vradius[-vte]

vcenter[vte]  ==  vcenter[-vte]

vec[n]

euklidsphere[vec[n]][x, y] + euklidsphere[-vec[n]][x, 
        y]

Simplify[Together[%]]

It follows that the spheres are identical as sets; the radius vectors with the same parameter values with respect to the centers are opposite: The sum of the radius vectors is twice the center.

4 The generalized angle

4.1 Definition. Cases

The consideration of this section are valid for hyperspheres in any dimension (with exception of plotting commands, of course). We consider two spheres, represented by the spacelike vectors v, w. Since the map v -> sphere[v] is equivariant, the generalized cosinus gencos of these vectors already yields a complete invariant system  in the space of oriented sphere pairs, as can be proved easily.

?gencos

gencos[v_List,  w_List]  :=  
      Simplify[pssp[v,  w]/Sqrt[pssp[v, v] pssp[w, w]]]

We have to distinguish three cases: The linear hull of {v, w} can be 1) euclidean, 2) pseudo-euclidean, 3) isotropic. Here and in the following we assume that the spheres are distinct, i.e. the vectors v, w are linearely independent.

?spacevec

spacevec[x, y, z, r]

pssp[spacevec[x, y, z, r], spacevec[X, Y, Z, R]]

Together[%]

Since the distance square of the centers of the spheres in the euclidean space is

d^2 = { X - x, Y - y, Z - z} . { X - x, Y - y, Z - z}

we obtain the conformal invariant gencos of two spheres in terms of their euclidean invariants:

?confinv

confinv[r, R, d]

If the spheres intersect this invariant is the cosinus of the intersection angle; If they do not intersect the invariant

ArcCosh[confinv[r, R, d]]

has been considered by H. S. M. Coxeter in [1]; he called it the inversive distance of the spheres.

4.2 Intersecting spheres

The spheres intersect in a circle (resp. (n-2)-sphere) iff v, w span an euclidean subspace. Then gencos[v,w] obviously is the cosinus of the  usual angle between the vectors v, w, and this is valid iff -1 < gencos[v,w] < 1:

phi[v_,  w_]  :=  ArcCos[gencos[v, w]]

is the angle between the normal vectors (= radii vectors) of the spheres along the intersection. This follows easily from the cos-theorem.

Example

vx  =  spherevec[0,  1]

euklidsphereplot3D[vx]

sph1  =  %

euklidsphereplot3D[spherevec[0.7,  0.8]]

sph2  =  %

Show[{sph1, sph2}]

gencos[vx,  spherevec[0.7,  0.8]]

wx =  spherevec[x,  r]

pssp[%,  %]

Simplify[%]

gencos[vx,  wx]

The absolute value of this invariant equals 1 in the cases x = ±(1-r) or x = ±(1+r), and this means that the spheres are tangent to each other, see also section 4.4 .

Intersecting planes

A simple example of symbolic calculations: let us take two general planes:

na =  planev[x, y, z, h]

pssp[na, na]

nb =  planev[X, Y, Z, H]

Simplify[(pssp[na, nb])^2/(pssp[na, na]*
            pssp[nb, nb])]

And this is the cos^2 of the intersection angle, if the planes are not parallel, and 1 else, as expected!If the planes are parallel we have

nc  =   planev[x, y, z, k]

Simplify[(pssp[na, nc])^2/(pssp[na, na]*
            pssp[nc, nc])]

Parallel planes are tangent to each other at infinity.

4.3 Non-intersecting spheres

In the pseudo-euclidean case the othogonal complement of the linear hull[v, w] is euclidean, thus, it does not contain isotropic vectors, and the intersection of the spheres is empty, and vice versa. We have

fct[x1_,  r1_, x2_,  r2_]  :=  
      gencos[spherevec[x1,  r1],  spherevec[x2,  r2]]

fct[x1, r1, x2, r2]

Since d = x1 - x2 is the euclidean distance between the centers of the spheres, we obtain the conformal invariant of two spheres with radii r1, r2 and distance d of its centers as

?confinv

confinv[r1_,  r2_,  d_]  :=  
      (r1^2  +  r2^2  - d^2)/(2 r1 r2)

More generally we derived this expression already in 4.1.

In the non-intersecting case one has Abs[confinv[r1, r2, d]]>1, and vice versa.

confinv[r1,  r2,  d]

4.4 Tangential spheres

The spheres are tangent to each other at the common point, which correspond to the isotropic vector of the linear hull[v,w], iff one has Abs[confinv[r1, r2, d]]=1.
We have two subcases:

a)  d = |r1 + r2|

confinv[r1,  r2,  r1  +  r2]

Simplify[%]

spherevec[1 + r,  r]

Simplify[2 r %]

Table[
      euklidsphereplot3D[spherevec[1 + r,  r]],  {r,  0.2,  1,  
        0.2}]

tangser  =  %

Show[%]

euklidsphereplot3D[stb[4]]

spst = %

Append[tangser, spst]

taser  =  %

Show[taser]

tanser  = %

MVShow3D[tanser]

ThreeScript["tangsph.ts", tanser]

The common point of all these spheres is

stepro[{1, 0, 0, 0, 1}]

b)  d = |r1  - r2|

confinv[r1,  r2,  r1  -  r2]

Simplify[%]

spherevec[1 - r,  r]

Simplify[2 r %]

Table[
      euklidsphereplot3D[spherevec[1 - r,  r]],  {r,  0.2,  1,  
        0.2}]

tangser2  =  %

Show[tangser2]

tanser2  = %

MVShow3D[tanser2]

ThreeScript["tangsph2.ts", tanser2]

The common point of all these spheres is again

stepro[{1, 0, 0, 0,  1}]

4.5 Intersection of a sphere and a plane

Now we consider the intersection of a plane with a sphere. For the conformal invariant of a sphere and a plane we calculate

pssp[spacevec[x, y, z, r],  planev[X, Y, Z, p]]

Together[%]

\(%^2/\((pssp[spacevec[x, y, z, r], spacevec[x, y, z, r]]*
			pssp[\ 
            planev[X, Y, Z, p], \ planev[X, Y, Z, p]])\)\)

Simplify[%]

If plane and sphere really intersect, this is the cosinus square of the intersection angle. In the general case it corresponds to the conformal invariant confinv of two spheres. The plane is tangent to the sphere iff this expression equals one.

Example

For the intersection of a unit sphere with center {x,0,0} and the y,z-plane through the origin we obtain

sph[x_]  :=  spherevec[x,  1]

sph[x]

pl  =  planev[1, 0, 0, 0]

gencos[pl,  sph[x]]

pssp[sph[x],  sph[x]]

Simplify[%]

{euklidsphereplot3D[stb[1]], 
      euklidsphereplot3D[{1, 0, 0, -1/2, 1/2}]}

Show[%]

MVShow3D[%]

5 Geodesics of the sphere space

5.1 Introduction

In subsection 3.5 we identified the space of oriented spheres with the hypersurface  H: pssp[x,x] = 1 of spacelike unit vectors in the pseudo-euclidean vector space: a "pseudo-hypersphere".  It is a pseudo-Riemannian symmetric space of constant curvature. As in Riemannian sphere geometry, the geodesics of  H are the intersections of the pseudo-hypersphere H with 2-planes through the center o. But now we  have three kinds of such planes, and the corresponding geodesics. We take the equatorial sphere stb[4] as starting point of  the geodesics. To find the types of geodesic we may restrict to the linear hull[stb[5],stb[1]], as a subspace of the tangent space of this hyper-hyperboloid at stb[4]: look at the isotropy action, the action of the subgroup of the Möbius group preserving stb[4] on the tangential space of H at stb[4]. Obviously, this is the usual standard action of the pseudo-orthogonal group O(3,1) on the 4-dimensional pseudo-euclidean vector space of index 1. 
We remember the following situation:

se1 = euklidsphereplot3D[stb[3]]

se2  = euklidsphereplot3D[stb[4]]

Show[{se1, se2}, PlotRange   -> {-1, 1}]

The intersection of the spheres defined by stb[4] and stb[3] is the unit circle in the x,y-plane around o. The spheres containing this circle are exactly represented by the unit vectors of the linear hull[stb[4], stb[3]]. All the spheres corresponding to spacelike vectors of the linear hull[stb[1], stb[2], stb[5]] intersect all the spheres of this family orthogonally.

5.2 Timelike geodesics

A timelike geodesics is the intersection of H with a pseudo-Euclidean plane, say, spanned by stb[4], stb[5]. This 2-plane contains two isotropic directions, i.e. points of the 3-sphere. The geodesic consists of all spheres which intersect all the spheres through these two points  in an angle Pi/2. If we take one of this point as the north pole (corresponding to stb[4,4]) and the other as the origin o , and the "spheres" through these points are the planes through the origin (under stereographical projection).Therefore  the geodesic is the set of all spheres with center o. The tangential vector of this geodesic at t=0 is stb[5]. Therefore the geodesic is timelike.

Five spheres on this timelike geodesic (for better visibility we only plot the halfspheres):

Table[euklidsphereplot3D[{0, 0, 0, Cosh[t], Sinh[t]},-1,1,0,1],
{t, -1, 1, 1/2}]

grr  = %

Show[grr]

grr1  = %

MVShow3D[grr1]

ThreeScript["spheres1.ts", grr1]

Replacing stb[4] by stb[1] gives a Möbius equivalent, also again timelike geodesic. In the euclidean picture, this geodesic starts at the "sphere" defined by stb[1], that is the y,z-plane. It is defined as the one-parameter sphere family corresponding to the spacelike vectors {Cosh[t], 0, 0, 0, Sinh[t]}:

Table[
      euklidsphereplot3D[{Cosh[t], 0, 0, 0, Sinh[t]}, -1, 0, -1, 
        1], {t, -5, 5}]

Show[%,  PlotRange  ->  All]

tgr  =  %

One sees that we obtained a one parameter family of pairwise non intersecting spheres,  containing the y,z-plane and contracting to two points, which correspond to the origin and the infinite point in our first example.These points are

{stepro[stb[1]  + stb[5]], stepro[stb[1]  - stb[5]]}

The points are not arbitrary; they are in involution with respect to the starting sphere (here: the y,z-plane) of the geodesic. Consider all circles through the two points: they cut the y,z-plane orthogonally. Since the geodesic is generated moving the starting sphere by a 1-parameter group of Möbius transformations  preserving these two points (the pseudo-orthogonal group of the subspace spanned by stb[1], stb[5]):

{{Cosh[t], 0,  0,  0,  Sinh[t]}, 
{0,  1,  0,  
          0,  0}, 
{0,  0,  1,  0,  
          0}, 
{0,  0,  0,  1,  
          0}, 
{Sinh[t], 0,  0,  0,  Cosh[t]}} . {1, 0, 
        0, 0, 0}

{{Cosh[t], 0,  0,  0,  Sinh[t]}, 
{0,  1,  0,  
          0,  0}, 
{0,  0,  1,  0,  
          0}, 
{0,  0,  0,  1,  
          0}, 
{Sinh[t], 0,  0,  0,  Cosh[t]}} . {1, 0, 
        0, 0, 1}

{{Cosh[t], 0,  0,  0,  Sinh[t]}, 
{0,  1,  0,  
          0,  0}, 
{0,  0,  1,  0,  
          0}, 
{0,  0,  0,  1,  
          0}, 
{Sinh[t], 0,  0,  0,  Cosh[t]}} . {1, 0, 
        0, 0, -1}

all the spheres of the geodesic intersect the circle bundle through the two points orthogonally, and therefore the geodesic is the family of orthogonal trajectories of the bundle. Thus, the recipe to conctruct a timelike geodesic through a given sphere is the following: Take a point outside the sphere, reflect it at the sphere, consider the bundle of all circles through the point pair in involution obtained this way: The geodesic is the family of orthogonal trajectories of this bundle, which are spheres. All these geodesics are Möbius equivalent. The set of these geodesics is a six-dimensional homogeneous space under the action of the Möbius group.

5.3 Spacelike geodesics

The spacelike geodesic corresponding to the intersection with the linear hull[stb[1], stb[4]] is the "circle"
{Sin[t], 0, 0, Cos[t],0}.
Five spheres on this spacelike geodesic (again we plot halfspheres only):

Table[
      euklidsphereplot3D[
{Sin[t], 0, 0, Cos[t], 
          0}, -1, 0, -1, 1], 
{t, -1/3, 1/3, 
        1/6}]

Show[%]

spacege  =  %

MVShow3D[spacege]

ThreeScript["spacege4.ts", spacege]

These spheres intersect in the unit circle centered in o in the y,z-plane :

ParametricPlot3D[{0, Cos[t],  Sin[t]}, {t, 0, 2 π}]

circ1yz = %

Show[{spacege, circ1yz}]

spacegec  = %

MVShow3D[spacegec,  MVLineTubeSize  ->  .01]

The whole geodesic consists of all the spheres through the fixed circle defined by the orthogonal complement of
linearhull{{Cos[t], 0, 0, Sin[t], 0}, { Cos[t + Pi/2], 0, 0, Sin[t + Pi/2],0}} = linearhull{stb[1], stb[4]}. 
As a special case of a spacelike geodesic we get the bundle of planes through a line. Since the space of all circles is six-dimensional, we can state: The space of all spacelike geodesics in the sphere space is a six-dimensional homogeneous space under the action of the Möbius group, equivariant isomorphic to the space of all circles.

Table[
      euklidsphereplot3D[
{Sin[t], 0, 0, Cos[t], 
          0}, -1, 0, -1, 1], 
{t, 0, Pi, 
        Pi/8}]

Show[%]

spacegeo1  =  %

MVShow3D[spacegeo1]

Show[{spacegeo1, circ1yz}]

spacegeo2  =  %

MVShow3D[spacegeo2,  MVLineTubeSize  ->  .05]

5.4 Isotropic geodesics

An isotropic geodesic through stb[4] is the intersection of H with the isotropic plane through the origin, spanned by the vectors stb[4] and stb[1] + stb[5]. This is a generating line contained in the hyperboloid H, with parameter representation
 {t,0,0,1, t}.
The starting point: stb[4], the unit sphere with centre o, and any other sphere of this geodesic have the isotropic span:
linearhull[{1,0,0,0,1}, {0,0,0,1,0}], and therefore they intersect in the point

cpt  = {1, 0, 0, 0, 1};
stepro[cpt]

Because spheres which intersect in exactly one point must have a common tangential plane at this point,  the isotropic geodesic under consideration is the 1-parameter family of spheres through the point {1,0,0} having there the same tangent plane as the unitsphere with centre o, and this is the plane cutting the x-axis orthogonally. This tangent plane corresponds to the spacelike vector

planev[1, 0, 0, 1]

and belongs to the geodesic. If vec[n] is a spacelike unit vector the corresponding sphere of which is tangent to the unit sphere at the point {1,0,0}, it must be a solution of the following system

\(Solve[{pssp[cpt, {x1, x2, x3, x4, x5}] == 
          0, 
		pssp[stb[4], {x1, x2, x3, x4, x5}] == 1}, 
	{x1, x2, x3, 
        x4, x5}]\)

The general solution is

vs =  {t, a, b, 1, t}

pssp[vs, vs]

Since vs should be a spacelike unit vector it follows a=b=0: The isotropic geodesic consists exactly of all spheres tangent to the unit sphere at the point {1,0,0}.

Here are 5 spheres on this isotropic geodesic:

Table[
      euklidsphereplot3D[{t, 0, 0, 1, t}, 0, 1,  -1, 
        1, 
PlotPoints ->  20], {t, -3, 3}]

Show[%]

isogeod1 = %

MVShow3D[isogeod1]

ThreeScript["isoge1.ts", isogeod1]

Another isotropic geodesic through the unit sphere in the linear hull under consideration is given by the spacelike vectors

{t, 0, 0, 1, -t}

Table[
      euklidsphereplot3D[{t, 0, 0, 1, -t}, -1, 0,  -1, 
        1, 
PlotPoints ->  20], {t, -3, 3}]

Show[%]

The geodesic consists of all spheres having at the common point

stepro[{1, 0,  0,  0,  -1}]

the plane parallel to the y,z-plane as its common tangential plane. Obviously, there exists a 2-parameter family of isotropic geodesics through stb[4], corresponding to the points of the sphere, i.e. to the generators of the isotropic cone in the 4-dimensional tangential hyperplane of the hyper-hyperboloid of spacelike unit vectors at stb[4]. The set of all isotropic geodesics is a five-dimensional homogeneous manifold under the action of the Möbius group.

6 ThreeScript, MathGL3d

If you want to generate *.ts files for animations outside of Mathematica, activate the next section. To use the animation program MathGL3d in this notebook, activate section 6.2.

6.1 Threescript

The ThreeScript commands generate *.ts files appropriate  for certain animation programs, e.g. the  freeware program  MathGL3d  created by Jens-Peer Kuska [2] (for Windows95, Linux and other platforms).

Under Windows the commercial 3D-Animation programs MathLive and DynamicVisualizer are available. To produce .ts-files one

Needs["Graphics`ThreeScript`"]

?ThreeScript

The standard extension of these files is ts: "filename.ts" or "Path/filename.ts", e.g.

ThreeScript["tangr.ts",  tangr3]

6.2 MathGL3d

The free animation program MathGL3d by Jens-Peer Kuska [2] can be used from a command line and within a Mathematica session. A new version of this program is in preparation. Jens-Peer Kuska has sent me a new startup file OpenGLViewer.m, and friendly permitted to add it to the item "Spheres". Replace the old file OpenGLViewer.m by this new OpenGLViewer.m. If you did this, changes of Mathematica's ParametricPlot3D, especially the sometimes disturbing generating of MeshGraphics3D objects, will not take place.
Establish the link to the current Mathematica session with the command:

Get["MathGL3d`OpenGLViewer`"]

If you want to close the link, enter

MVClose[]

The following command clears the Viewer window; this and other MathGL3d commands you may enter also with  palettes.

MVClear[];

MVShow3D[\[Placeholder], MVPolygonShading -> MVSmooth, 
      MVGrayBackground -> False, MVTexture -> None, 
      MVTextureMapType -> MVMeshUVMapping, 
      MVScaleTexture -> {1, 1}, 
      MVTranslateTextureTarget -> Automatic, 
      MVRotateTextureTarget -> Automatic, 
      MVScaleTextureTarget -> Automatic]

MVShow3D[gre1
, MVPointSphereSize -> 0.8,  
      MVLineTubeSize  ->  0.03]

MVShow3D[tgr3,  MVPointSphereSize  ->  1.8]

MVShow3D[ex123, MVPointSphereSize -> 1.8]

MVShow3D[isogeod1]

Options[MVShow3D]

MVClose[]

Appendix

Obsolete: plotsphere

A procedure which directly plots the sphere corresponding to the spacelike vector nv is defined in mspher.m as follows:

?plotsphere

\(plotsphere[nv_List, \ opts___]\  := \ 
      Module[{a, \ p, \ r, \ m4\  = \ nul[4], \ z\  = \ nul[4], \ 
	vx, \ 
          tt\  = \ IdentityMatrix[4], \ u, \ v}, 
	
        a\  = \ pssp[nv, \ nv]; 
	
        If[Length[nv]\  != \ 5\  || \ a\  <= \ 0, \ 
          Return[Print[
				"Vector ", \ 
              nv, \ ": incorrect argument!"]]]; 
	
        m4\  = \ Take[nv, \ 4]; \ a\  = \ m4\  . \ m4; \ 
        m4\  = \ m4/Sqrt[a]; 
	p\  = \ nv[\([5]\)]/Sqrt[a]; \ 
        r\  = \ Sqrt[1\  - \ p^2]; 
	z\  = \ p\ m4; \ 
	
        tt\  = \ orthoframe[m4]; 
	
        vx[u_, \ v_]\  := \ 
          z\  + \ \((\((tt[\([1]\)]\ Cos[u]\  + \ tt[\([2]\)]\ Sin[u])\)\ Cos[
                      v]\  + \ 
		tt[\([3]\)]\ Sin[v])\)\ r; 
	
        ParametricPlot3D[
          sterproj[vx[u, \ v]], \ \ {u, \ 0, \ 2\ Pi}, \ {v, \ \(-Pi\)/2, \ 
            Pi/2}, \ opts]]\)

?orthoframe

Moving frame on the sphere manifold:pstg5frame

We construct a pseudo-orthogonal moving 5-frame of the hypersurface (H), using its parameter representation spunit. We follow the procedure which us led to the definition of tg4frame, see the notebook euspher.nb:

pstg1[s_, t_, u_, v_]  :=  
      normalize[D[spunit[x,  t, u, v], x] /. x  ->  s]

pstg1[s, t, u, v]

pstg2[s_, t_, u_, v_]  :=  
      normalize[D[spunit[s,  x, u, v], x] /. x  ->  t]

pstg2[s, t, u, v]

pstg3[s_, t_, u_, v_]  :=  
      normalize[D[spunit[s,  t, x, v], x] /. x  ->  u]

pstg3[s, t, u, v]

pstg4[s_, t_, u_, v_]  :=  
      normalize[D[spunit[s,  t, u, x], x] /. x  ->  v]

pstg4[s, t, u, v]

spunit[s, t, u, v]

Simplifying these vectors we defined in mspher.m the moving orthonormal tangentially adapted  5-frame of the spacelike unitsphere in the pseudo-euclidean 5-space,  analogous to tg4frame:

?pstg5frame

\(pstg5frame[s_, \ t_, \ u_, \ 
        v_]\  := 
	{{\(-\ Sin[s]\), Cos[s], 0, 0, 
          0}, 
		{\(-Cos[s]\)\ Sin[t], \(-\ Sin[s]\)\ Sin[t], Cos[t], 0, 
          0}, 
		{\(-Cos[s]\)\ Cos[t]\ \ Sin[u], \(-Cos[t]\)\ \ Sin[
              s]\ Sin[u], \(-\ Sin[t]\)\ Sin[u], Cos[u], 
          0}, 
		{Cos[s]\ Cos[t]\ Cos[u]\ Sinh[v], 
          Cos[t]\ Cos[u]\ Sin[s]\ Sinh[v], Cos[u]\ Sin[t]\ Sinh[v], 
          Sin[u]\ Sinh[v], Cosh[v]}, 
		{Cos[s]\ Cos[t]\ Cos[u]\ Cosh[v], 
          Cos[t]\ Cos[u]\ Cosh[v]\ Sin[s], Cos[u]\ Cosh[v]\ Sin[t], 
          Cosh[v]\ Sin[u], Sinh[v]}}\)

pst5  = pstg5frame[s,  t,  u,  v];

Table[
      Simplify[pssp[pst5[[i]],  pst5[[j]]]],  {i, 5}, {j, 
        5}]

Simplify[Det[pst5]]

spunit[s, t, u, v]  +  
      spunit[s + Pi, -t, -u, -v]

We conclude that changing the sign of v = spunit[s,t,u,v] changes the orientation of the orthogonal pseudo-euclidean subspace, if we consider the pseudo-euclidean 5-space as oriented.

The isotropic vectors in the tangential space spanned by (pstg5frame[[i]], i = 1,...,4) are:

\(\(isotr[s_, t_, u_, v_]\)[x_, 
        y_]\  := 
\(\(\(pstg5frame[s, \ t, \ u, \ v]\)[\([1]\)] 
          Cos[x]\ Cos[
            y]\  + \ 
\(pstg5frame[s, \ t, \ u, \ v]\)[\([2]\)]\ Sin[x]\ Cos[
            y]\  + 
\(pstg5frame[s, \ t, \ u, \ v]\)[\([3]\)]\ Sin[
            y]\  + 
\(pstg5frame[s, \ t, \ u, \ 
            v]\)[\([4]\)]\)\(
\)\(	\)\)\)

iso  =  Simplify[isotr[s, t, u, v][x, y]]

Simplify[pssp[iso, iso]]