jncf18_vector.ipynb

Kernel: SageMath (stable)

Smooth manifolds, vector fields and tensor fields

This woksheet accompanies the lecture Symbolic tensor calculus on manifolds at JNCF 2018.

Click here to download the worksheet file (ipynb format). To run it, you must start SageMath with the Jupyter notebook, via the command sage -n jupyter

%display latex

Let us restore the 2-sphere manifold constructed in the preceeding worksheet:

M = Manifold(2, 'M')
U = M.open_subset('U')
XU.<x,y> = U.chart()
V = M.open_subset('V')
XV.<xp,yp> = V.chart("xp:x' yp:y'")
M.declare_union(U,V)
XU_to_XV = XU.transition_map(XV, 
                             (x/(x^2+y^2), y/(x^2+y^2)), 
                             intersection_name='W',
                             restrictions1= x^2+y^2!=0, 
                             restrictions2= xp^2+yp^2!=0)
XV_to_XU = XU_to_XV.inverse()
M.atlas()

We also reconstruct the point $p$:

p = U((1,2), chart=XU, name='p')
print(p)

Point p on the 2-dimensional differentiable manifold M

and the embedding $\mathbb{S}^2 \to \mathbb{R}^3$:

R3 = Manifold(3, 'R^3', r'\mathbb{R}^3')
XR3.<X,Y,Z> = R3.chart()
Phi = M.diff_map(R3, {(XU, XR3): 
                       [2*x/(1+x^2+y^2), 2*y/(1+x^2+y^2),
                        (x^2+y^2-1)/(1+x^2+y^2)],
                      (XV, XR3): 
                       [2*xp/(1+xp^2+yp^2), 2*yp/(1+xp^2+yp^2),
                        (1-xp^2-yp^2)/(1+xp^2+yp^2)]},
                  name='Phi', latex_name=r'\Phi')
Phi.display()

graph = XU.plot(chart=XR3, mapping=Phi, number_values=25, 
                label_axes=False) + \
        XV.plot(chart=XR3, mapping=Phi, number_values=25, 
                color='green', label_axes=False) + \
        p.plot(chart=XR3, mapping=Phi, label_offset=0.05)
show(graph, viewer='threejs', online=True)

Finally we reconstruct the scalar field $f$:

f = M.scalar_field({XU: 1/(1+x^2+y^2), XV: (xp^2+yp^2)/(1+xp^2+yp^2)},
                   name='f')
f.display()

and assign the Python variable CM to the algebra of scalar fields:

CM = M.scalar_field_algebra()
CM

Tangent vectors

The tangent space at the point $p$ introduced above is generated by

Tp = M.tangent_space(p)
Tp

It is a vector space over $\mathbb{R}$, which is represented by Sage's Symbolic Ring:

print(Tp.category())

Category of finite dimensional vector spaces over Symbolic Ring

The dimension of $T_p M$ is the same as that of $M$:

dim(Tp)

Tangent spaces are implemented as a class inherited from TangentSpace via the category framework:

type(Tp)

The class TangentSpace actually inherits from the generic class FiniteRankFreeModule, which, in SageMath, is devoted to free modules of finite rank without any distinguished basis:

isinstance(Tp, FiniteRankFreeModule)

Two bases of $T_p M$ are already available: those generated by the derivations at $p$ along the coordinates of charts XU and XV respectively:

Tp.bases()

None of these bases is distinguished, but one if the default one, which simply means that it is the basis to be considered if the basis argument is skipped in some methods:

Tp.default_basis()

A tangent vector is created as an element of the tangent space by the standard SageMath procedure new\_element = parent(...), where ... stands for some material sufficient to construct the element:

vp = Tp((-3, 2), name='v')
print(vp)

Tangent vector v at Point p on the 2-dimensional differentiable manifold M

Since the basis is not specified, the pair $(-3,2)$ refers to components with respect to the default basis:

vp.display()

We have of course

vp.parent()

vp in Tp

As other manifold objects, tangent vectors have some plotting capabilities:

graph += vp.plot(chart=XR3, mapping=Phi, scale=0.5, color='gold')
show(graph, viewer='threejs', online=True)

Module of vector fields

The $C^\infty(M)$-module of vector fields on $M$, $\mathfrak{X}(M)$, is obtained as

YM = M.vector_field_module()
YM

YM.category()

YM.base_ring() is CM

$\mathfrak{X}(M)$ is not a free module (at least its SageMath implementation does not belong to the class FiniteRankFreeModule):

isinstance(YM, FiniteRankFreeModule)

This is because $M=\mathbb{S}^2$ is not a parallelizable manifold:

M.is_manifestly_parallelizable()

Via the category mechanism, the module $\mathfrak{X}(M)$ is implemented by a dynamically-generated subclass of the class VectorFieldModule, which is devoted to modules of vector fields on non-parallelizable manifolds:

type(YM)

On the contrary, the set $\mathfrak{X}(U)$ of vector fields on $U$ is a free module of finite rank over the algebra $C^\infty(U)$:

YU = U.vector_field_module()
isinstance(YU, FiniteRankFreeModule)

YU.base_ring()

This is because the open subset $U$ is a parallelizable manifold:

U.is_manifestly_parallelizable()

being a coordinate chart domain:

U.is_manifestly_coordinate_domain()

We can check that in the atlas of $U$, at least one chart has $U$ for domain:

U.atlas()

The rank of $\mathfrak{X}(U)$ is the manifold's dimension:

rank(YU)

Via the category mechanism, the free module $\mathfrak{X}(U)$ is implemented by a dynamically-generated subclass of the class VectorFieldFreeModule, which is devoted to modules of vector fields on parallelizable manifolds:

type(YU)

The class VectorFieldFreeModule is itself a subclass of the generic class FiniteRankFreeModule:

class_graph(
sage.manifolds.differentiable.vectorfield_module.VectorFieldFreeModule
).plot()

Since $U$ is a chart domain, the free module $\mathfrak{X}(U)$ is automatically endowed with a basis: the coordinate frame associated to the chart:

YU.bases()

Let us denote by eU this frame. We can set eU = YU.bases()[0] or alternatively

eU = YU.default_basis()
eU

Another equivalent instruction would have been eU = U.default_frame().

Similarly, $\mathfrak{X}(V)$ is a free module, endowed with the coordinate frame associated to stereographic coordinates from the South pole, which we denote by eV:

YV = V.vector_field_module()
YV.bases()

eV = YV.default_basis()
eV

If we consider the intersection $W=U\cap V$, we notice its module of vector fields is endowed with two bases, reflecting the fact that $W$ is covered by two charts: $(W,(x,y))$ and $(W,(x',y'))$:

W = U.intersection(V)
YW = W.vector_field_module()
YW.bases()

Let us denote by eUW and eUV these two bases, which are actually the restrictions of the vector frames eU and eV to $W$:

eUW = eU.restrict(W)
eVW = eV.restrict(W)
YW.bases() == [eUW, eVW]

The free module $\mathfrak{X}(W)$ is also automatically endowed with automorphisms connecting the two bases, i.e. change-of-frame operators:

W.changes_of_frame()

The first of them is

P = W.change_of_frame(eUW, eVW)
P

It belongs to the general linear group of the free module $\mathfrak{X}(W)$:

P.parent()

and its matrix is deduced from the Jacobian matrix of the transition map XV $\to$ XU:

P[:]

An example of vector field

We introduce a vector field $v$ on $M$ by

v = M.vector_field(name='v')
v[eU, 0] = f.restrict(U)
v[eU, 1] = -2
v.display(eU)

Notice that at this stage, we have defined $v$ only on $U$, by setting its components in the vector frame eU, either explicitely as scalar fields, like the component $v^0$ set to the restriction of $f$ to $U$ or to some symbolic expression, like for the component $v^1$: the $-2$ will be coerced to the constant scalar field of value $-2$. We can ask for the scalar-field value of a components via the double-bracket operator, since eU is the default frame on $M$, we don't have to specify it:

v[[0]]

v[[0]].display()

Note that the single bracket operator returns a chart function of the component:

v[0]

The restriction of $v$ to $W$ is of course

v.restrict(W).display(eUW)

Since we have a second vector frame on $W$, namely eVW, and the change-of-frame automorphisms are known, we can ask for the components of $v$ with respect to that frame:

v.restrict(W).display(eVW)

Notice that the components are expressed in terms of the coordinates $(x,y)$ since they form the default chart on $W$. To have them expressed in terms of the coordinates $(x',y')$, we have to add the restriction of the chart $(V,(x',y'))$ to $W$ as the second argument of the method display():

v.restrict(W).display(eVW, XV.restrict(W))

We extend the expression of $v$ to the full vector frame XV by continuation of this expression:

v.add_comp_by_continuation(eV, W, chart=XV)

We have then

v.display(eV)

At this stage, the vector field $v$ is defined in all $M$. According to the hairy ball theorem\index{hairy ball theorem}, it has to vanish somewhere. Let us show that this occurs at the North pole, by first introducing the latter, as the point of stereographic coordinates $(x',y')=(0,0)$:

N = M((0,0), chart=XV, name='N')
print(N)

Point N on the 2-dimensional differentiable manifold M

As a check, we verify that the image of $N$ by the canonical embedding $\Phi: \mathbb{S}^2 \to \mathbb{R}^3$ is the point of Cartesian coordinates $(0,0,1)$:

XR3(Phi(N))

The vanishing of $\left. v\right| _N$:

v.at(N).display()

On the other hand, $v$ does not vanish at the point $p$ introduced above:

v.at(p).display()

We may plot the vector field $v$ in terms of the stereographic coordinates from the North pole:

v.plot(chart=XU, chart_domain=XU, max_range=2, 
       number_values=5, scale=0.4, aspect_ratio=1)

or in term of those from the South pole:

v.plot(chart=XV, chart_domain=XV, max_range=2, 
       number_values=9, scale=0.05, aspect_ratio=1)

Thanks to the embedding $\Phi$, we may also have a 3D plot of $v$ atop of the 3D plot already obtained:

graph_v = v.plot(chart=XR3, mapping=Phi, chart_domain=XU, 
                 number_values=7, scale=0.2) + \
          v.plot(chart=XR3, mapping=Phi, chart_domain=XV, 
                 number_values=7, scale=0.2)
show(graph + graph_v, viewer='threejs', online=True)

Note that the sampling, performed on the two charts XU and XV is not uniform on the sphere. A better sampling would be acheived by introducing spherical coordinates.

Some details about the implementation of vector fields

Since $M$ is not parallelizable, the fundamental representation of a vector field on $M$ is via its restrictions to parallelizable open subsets; they are stored in the dictionary _restrictions, whose keys are the open subsets:

v._restrictions

Let us consider one of these restrictions, for instance the restriction to $U$:

vU = v._restrictions[U]
vU is v.restrict(U)

Since $U$ is a parallelizable open subset, the fundamental representation of $v$ on it is via its components with respect to (possibly various) vector frames on $U$; they are stored in the dictionary _components, whose keys are the vector frames:

vU._components

v._restrictions[W]._components

Let us perform some algebraic operation on $v$:

w = v + f*v
w

The code for the addition is defined is accessible via

v.__add__??

As for the addition of scalar field (cf. this worksheet), we see that __add__() is implemented at the level of the generic class Element from which VectorField inherits. When both operands of the addition have the same parent, as here, __add__() invokes the method _add_() (note the single underscore on each side of add). This operator is implemented at the level of TensorField, as it can be checked from the source code:

v._add_??

vU._add_??

Action of a vector field on a scalar field

vf = v(f)
vf

vf.display()

Tensor fields

Let us start with a 1-form, namely the differential of $f$:

df = f.differential()
df

print(df)

1-form df on the 2-dimensional differentiable manifold M

df.tensor_type()

v.tensor_type()

We construct a tensor of type $(1,1)$ by taking the tensor product $v\otimes \mathrm{d}f$:

t = v * df
t

t.display()

t.display(eV)

t.display_comp()

t.parent()

t._restrictions