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

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

*NB:* a version of SageMath at least equal to 7.5 (8.2 for the SymPy part) is required to run this notebook.

First we set up the notebook to use LaTeX display:

In [1]:

```
%display latex
```

Manifolds are constructed via the global function `Manifold`

:

In [2]:

```
M = Manifold(2, 'M')
print(M)
```

By default, `Manifold`

returns a manifold over $\mathbb{R}$:

In [3]:

```
M.base_field()
```

Out[3]:

Other base fields must be specified with the optional keyword `field`

, like

```
M = Manifold(2, 'M', field='complex')
```

We may check that $M$ is a topological space:

In [4]:

```
M in Sets().Topological()
```

Out[4]:

Actually, $M$ belongs to the following categories:

In [5]:

```
M.categories()
```

Out[5]:

As we can see, by default, `Manifold`

constructs a smooth manifold.
If one would like to stick to the topological level, one should write

```
M = Manifold(2, 'M', structure='topological')
```

Smooth manifolds are implemented by the Python class `DifferentiableManifold`

:

In [6]:

```
type(M)
```

Out[6]:

The type of `M`

appears as `DifferentiableManifold-with-category`

because it is actually a subclass of `DifferentiableManifold`

, which is dynamically generated by SageMath's category mechanism:

In [7]:

```
isinstance(M,
sage.manifolds.differentiable.manifold.DifferentiableManifold)
```

Out[7]:

The class `DifferentiableManifold`

inherits from `TopologicalManifold`

:

In [8]:

```
isinstance(M, sage.manifolds.manifold.TopologicalManifold)
```

Out[8]:

We declare a chart, along with the symbols used to denote the coordinates (here $x=x^0$ and $y=x^1$) by

In [9]:

```
U = M.open_subset('U')
XU.<x,y> = U.chart()
XU
```

Out[9]:

Open subsets are implemented by a (dynamically generated) subclass of
`DifferentiableManifold`

, since they are manifolds in their own:

In [10]:

```
isinstance(U,
sage.manifolds.differentiable.manifold.DifferentiableManifold)
```

Out[10]:

Points on $M$ are created by passing their coordinates in a given chart:

In [11]:

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

The syntax `U(...)`

used to create $p$ as an element of $U$ reflects the
**parent/element pattern** employed in SageMath; indeed $U$ is the parent of $p$:

In [12]:

```
p.parent()
```

Out[12]:

Points are implemented by the class `ManifoldPoint`

.
The principal attribute of this class is a Python dictionary storing the point's coordinates
in various charts:

In [13]:

```
p._coordinates
```

Out[13]:

Of course, we can recover the point's coordinates by letting the chart act on the point:

In [14]:

```
XU(p)
```

Out[14]:

Let us introduce a second chart on $M$:

In [15]:

```
V = M.open_subset('V')
XV.<xp,yp> = V.chart("xp:x' yp:y'")
XV
```

Out[15]:

and declare that $M$ is covered by only two charts, i.e. that $M=U\cup V$:

In [16]:

```
M.declare_union(U,V)
```

In [17]:

```
M.atlas()
```

Out[17]:

We define the transition map `XU`

$\to$ `XV`

on $W = U\cap V$ as follows:

In [18]:

```
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)
XU_to_XV.display()
```

Out[18]:

The value of the argument `restrictions1`

means that
$W = U\setminus \{S\}$, where $S$ is the point of coordinates $(x,y)=(0,0)$,
while the value of `restrictions2`

means that
$W = V\setminus \{N\}$, where $N$ is the point of coordinates $(x',y')=(0,0)$.
Since $M=U\cup V$, we have then
$$
U = M \setminus \{N\},\qquad
V = M \setminus \{S\},\quad\mbox{and}\quad
W = M \setminus \{N, S\} .
$$

The transition map `XV`

$\to$ `XU`

is obtained by computing the inverse of the one defined above:

In [19]:

```
XU_to_XV.inverse().display()
```

Out[19]:

At this stage, the smooth manifold $M$ is fully specified, being covered by
one atlas with all transition maps specified. One may have recognized that
$M$ is nothing but the **2-dimensional sphere**:
$$
M = \mathbb{S}^2 ,
$$
with `XU`

(resp. `XV`

) being
the chart of \defin{stereographic coordinates}\index{stereographic!coordinates}
from the North pole $N$ (resp. the South pole $S$).

Since the transition maps have been defined, we can ask for the coordinates $(x',y')$ of the point $p$:

In [20]:

```
XV(p)
```

Out[20]:

This operation has updated the dictionary `_coordinates`

:

In [21]:

```
p._coordinates
```

Out[21]:

As a example of smooth map, let us consider the canonical embedding of $M=\mathbb{S}^2$ into $\mathbb{R}^3$. We need first to introduce $\mathbb{R}^3$ as a 3-dimensional smooth manifold, endowed with a single chart:

In [22]:

```
R3 = Manifold(3, 'R^3', r'\mathbb{R}^3')
XR3.<X,Y,Z> = R3.chart()
XR3
```

Out[22]:

The embedding $\Phi: \mathbb{S}^2 \to \mathbb{R}^3$ is then defined in terms of its coordinate expression in the two charts covering $M=\mathbb{S}^2$:

In [23]:

```
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()
```

Out[23]:

We may use $\Phi$ for graphical purposes, for instance to display the grids of the stereographic charts `XU`

(in red) and `XV`

(in green), with the point $p$ atop:

In [24]:

```
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)
```

Let us define a scalar field, i.e. a smooth map $M\to \mathbb{R}$:

In [25]:

```
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()
```

Out[25]:

Scalar fields are implemented by the class `DiffScalarField`

and the function chart representations are
stored in the attribute `_express`

of this class, which is a
Python dictionary whose keys are the various charts defined on $M$:

In [26]:

```
f._express
```

Out[26]:

One may wonder about the compatibility of the two coordinate expressions
provided in the definition of $f$. Actually, to enforce the compatibility, it is possible to declare the scalar field in a single chart, `XU`

say,
and then to obtain its expression in chart `XV`

by analytic continuation
from the expression in $W=U\cap V$, where both expressions are known, thanks
to the transition map `XV`

$\to$ `XU`

:

In [27]:

```
f0 = M.scalar_field({XU: 1/(1+x^2+y^2)})
f0.add_expr_by_continuation(XV, U.intersection(V))
f == f0
```

Out[27]:

The representation of the scalar field in a given chart, i.e. the public access to the directory `_express`

, is obtained via the method `coord_function()`

:

In [28]:

```
fU = f.coord_function(XU)
fU.display()
```

Out[28]:

In [29]:

```
fV = f.coord_function(XV)
fV.display()
```

Out[29]:

Both `fU`

and `fV`

are instances of the class `ChartFunction`

:

In [30]:

```
isinstance(fU, sage.manifolds.chart_func.ChartFunction)
```

Out[30]:

Mathematically, **chart functions** are real-valued functions on the codomain of the considered chart. They map coordinates to elements of the base field (here $\mathbb{R}$):

In [31]:

```
fU(1,2)
```

Out[31]:

In [32]:

```
fU(*XU(p))
```

Out[32]:

while scalar fields maps manifold points to $\mathbb{R}$:

In [33]:

```
f(p)
```

Out[33]:

Internally, each chart function stores coordinate expressions with respect to various computational engines:

- SageMath symbolic engine, based on the Pynac backend, with Maxima used for some simplifications or computation of integrals;
- SymPy (Python library for symbolic mathematics);
- in the future, more symbolic engines (e.g. Giac) or numerical ones will be implemented

The coordinate expressions are stored in the dictionary `_express`

,
whose keys are strings identifying the computational engines. By default
only SageMath symbolic expressions, i.e. expressions pertaining
to the so-called Symbolic Ring (`SR`

), are stored:

In [34]:

```
fU._express
```

Out[34]:

The public access to the dictionary `_express`

is performed via the
method `expr()`

:

In [35]:

```
fU.expr()
```

Out[35]:

In [36]:

```
type(fU.expr())
```

Out[36]:

Actually, `fU.expr()`

is a shortcut for `fU.expr('SR')`

since
`SR`

is the default symbolic engine. Note that the class
`Expression`

is that devoted to SageMath symbolic expressions.
The method `expr()`

can also be invoked to get the expression in
another symbolic engine, for instance SymPy:

In [37]:

```
fU.expr('sympy')
```

Out[37]:

In [38]:

```
type(fU.expr('sympy'))
```

Out[38]:

This operation has updated the dictionary `_express`

:

In [39]:

```
fU._express
```

Out[39]:

The default calculus engine for chart functions of chart `XU`

can
changed thanks to the method `set_calculus_method()`

:

In [40]:

```
XU.set_calculus_method('sympy')
fU.expr()
```

Out[40]:

We can have better (i.e. LaTeX) rendering of SymPy objects with

```
from sympy import init_printing
init_printing()
```

We don't use it here to make clear the distinction between `SR`

objects and SymPy ones.

Reverting to SageMath's symbolic engine:

In [41]:

```
XU.set_calculus_method('SR')
fU.expr()
```

Out[41]:

Symbolic expressions can be accessed directly from the scalar field,
`f.expr(XU)`

being a shortcut for `f.coord_function(XU).expr()`

:

In [42]:

```
f.expr(XU)
```

Out[42]:

In [43]:

```
f.expr(XV)
```

Out[43]:

The commutative algebra of scalar fields on $M$ is

In [44]:

```
CM = M.scalar_field_algebra()
CM
```

Out[44]:

In [45]:

```
CM.category()
```

Out[45]:

As for the manifold classes, the actual Python class implementing
$C^\infty(M)$ is inherited from `DiffScalarFieldAlgebra`

via the
category mechanism, hence it bares the name `DiffScalarFieldAlgebra-with-category`

:

In [46]:

```
type(CM)
```

Out[46]:

The class `DiffScalarFieldAlgebra-with-category`

is dynamically generated as a subclass of `DiffScalarFieldAlgebra`

with extra functionalities, like
for instance the method `is_commutative()`

:

In [47]:

```
CM.is_commutative()
```

Out[47]:

Let us have a look at the code of this method; it suffices to use the double question mark:

In [48]:

```
CM.is_commutative??
```

We see from the `File`

field that the code belongs to the category part of SageMath, not to the manifold part, where the class `DiffScalarFieldAlgebra`

is defined.

Regarding the scalar field $f$ introduced above, we have of course

In [49]:

```
f in CM
```

Out[49]:

Actually, in Sage language, `CM`

=$C^\infty(M)$ is the parent of `f`

:

In [50]:

```
f.parent() is CM
```

Out[50]:

The zero element of the algebra $C^\infty(M)$ is

In [51]:

```
CM.zero().display()
```

Out[51]:

while its unit element is

In [52]:

```
CM.one().display()
```

Out[52]:

Let us consider some operation in the algebra $C^\infty(M)$:

In [53]:

```
h = f + 2*CM.one()
h.display()
```

Out[53]:

In [54]:

```
h(p)
```

Out[54]:

We are going to investigate how the above addition is performed. For the Python interpreter
`h = f + 2*CM.one()`

is equivalent to `h = f.__add__(2*CM.one())`

,
i.e. the `+`

operator amounts to calling the method `__add__()`

on `f`

.
To have a look at the source code of this method, we use the double question mark:

In [55]:

```
f.__add__??
```

We see that the method `__add__()`

is implemented at the level
of the class `Element`

from which `DiffScalarField`

inherits, via
`CommutativeAlgebraElement`

.
In the present case, `left`

= `f`

and `right`

= `2*CM.one()`

have the same parent, namely the algebra `CM`

, so that the actual
result is computed in the line

```
return (<Element>left)._add_(right)
```

This invokes the method `_add_()`

(note the single underscore on each side of `add`

). This operator is
implemented at the level of `ScalarField`

, as it can be checked from the source code:

In [56]:

```
f._add_??
```

This reflects a general strategy in SageMath: the arithmetic Python operators `__add__`

, `__sub__`

, etc. are implemented at the
top-level class `Element`

, while the specific element subclasses,
like `ScalarField`

here, implement single-underscore methods `_add_`

, `_sub_`

, etc., which perform the actual computation
when both operands have the same parent.

Looking at the code, we notice that the first step is to search
for the charts in which both operands of the addition operator have a coordinate
expression. This is performed by the method `common_charts()`

;
in the current case, we get the two stereographic charts defined on $M$:

In [57]:

```
f.common_charts(2*CM.one())
```

Out[57]:

In general, `common_charts()`

returns the charts for which both operands
have already a known coordinate expression or for which a coordinate
expression can be computed by a known transition map, as we can see on the source code:

In [58]:

```
f.common_charts??
```

Once the list of charts in which both operands have a coordinate expression,
the addition is performed at the chart function level.
The code for the addition of chart functions defined on the same chart
is (recall that `fU`

is the chart function representing $f$ in chart `XU`

):

In [59]:

```
fU._add_??
```

We notice that the addition is performed on the symbolic expression
with respect to the symbolic engine currently at work (SageMath/Pynac, SymPy,...), as returned by
the method `expr()`

.
Let us recall that the user can change the symbolic engine at any time
by means of the method `set_calculus_method()`

, applied either to
a chart or to an open subset (possibly `M`

itself).
Besides, we notice in the code above that the result of the symbolic addition
is automatically simplified, by means of the method `_simplify`

.
It invokes a chain of simplification functions, which depends on the
symbolic engine (See here
for details).

Let us now discuss the second case in the `__add__`

method of
`Element`

, namely the case for which the parents of both operands are
different. This case is treated via SageMath **coercion model**, which allows one to deal with additions like

In [60]:

```
h1 = f + 2
h1.display()
```

Out[60]:

A priori, `f + 2`

is not a well defined operation, since the integer $2$ does not
belong to the algebra $C^\infty(M)$. However SageMath manages to treat it
because $2$ can be coerced (i.e. automatically and unambigously converted) via `CM(2)`

into a element of $C^\infty(M)$, namely the constant scalar field whose value is $2$:

In [61]:

```
CM(2).display()
```

Out[61]:

This happens because there exists a coercion map from the parent of $2$, namely the ring of integers $\mathbb{Z}$
(denoted `ZZ`

in SageMath), to $C^\infty(M)$:

In [62]:

```
2.parent()
```

Out[62]:

In [63]:

```
CM.has_coerce_map_from(ZZ)
```

Out[63]: