Stefan Sommer

Tom Fletcher

Xavier Pennec

published or not. The documents may come   L'archive ouverte pluridisciplinaire

Introduction to differential and Riemannian geometry

Stefan Sommer , * , Tom Fletcher * * and Xavier Pennec †

Introduction

When data exhibit non-linearity, the mathematical description of the data space must often depart from the convenient linear structure of Euclidean vector spaces. Nonlinearity prevents global vector space structure but we can nevertheless ask which mathematical properties from the Euclidean case can be kept while still preserving the accurate modelling of the data. It turns out that in many cases, local resemblance to a Euclidean vector space is one such property. In other words, up to some second order approximation, the data space can be linearized in limited regions while forcing a linear model on the entire space would introduce too much distortion.

The concept of local similarity to Euclidean spaces brings us exactly to the setting of manifolds. Topological, differential, and Riemannian manifolds are characterized by the existence of local maps, charts, between the manifold and a Euclidean space. These charts are structure preserving: They are homeomorphisms in the case of topological manifolds, diffeomorphisms in the case of differential manifolds, and, in the case of Riemannian manifolds, they carry local inner products that encode the non-Euclidean geometry.

The following sections describe these foundational concepts and how they lead to notions commonly associated with geometry: curves, length, distances, geodesics, curvature, parallel transport, volume form. In addition to the differential and Riemannian structure, we describe one extra layer of structure: Lie groups that are manifolds equipped with smooth group structure. Lie groups and their quotients are examples of homogeneous spaces. The group structure provides relations between distant points on the group and thereby additional ways of constructing Riemannian metrics and deriving geodesic equations.

Topological, differential, and Riemannian manifolds are often covered by separate graduate courses in mathematics. In this much briefer overview, we describe the general concepts, often sacrificing mathematical rigor to instead provide intuitive reasons for the mathematical definitions. For a more in-depth introduction to geometry, the interested reader may for example refer to the sequence of books by John M. Lee on topological, differentiable, and Riemannian manifolds [START_REF] John | Introduction to topological manifolds[END_REF][START_REF] John | Introduction to smooth manifolds[END_REF][START_REF] John | Riemannian manifolds[END_REF] or to the book on Riemannian geometry by do Carmo [dC92]. More advanced references include [START_REF] Kolář | Natural Operations in Differential Geometry[END_REF], [START_REF] Jost | Riemannian Geometry and Geometric Analysis[END_REF] and [START_REF] Mikhailovich | Geometry VI: Riemannian geometry[END_REF].

Manifolds

A manifold is a collection of points that locally, but not globally, resembles Euclidean space. When the Euclidean space is of finite dimension, we can without loss of generality relate it to R d for some d > 0. The abstract mathematical definition of a manifold specifies the topological, differential, and geometric structure by using charts, maps between parts of the manifold and R d , and collections of charts denoted atlases. We will discuss this construction shortly, however, we first focus on the case when the manifold is a subset of a larger Euclidean space. This viewpoint is often less abstract and closer to our natural intuition of a surface embedded in our surrounding 3D Euclidean space.

Let us exemplify this by the surface of the earth embedded in R 3 . We are constrained by gravity to live on the surface of the earth. This surface seems locally flat with two dimensions only, and we use two dimensional maps to navigate the surface. When traveling far, we sometimes need to change from one map to another. We then find charts that overlap in small parts, and we assume that the charts provide roughly the same view of the surface in those overlapping parts. For a long time, the earth was even considered to be flat because its curvature was not noticeable at the scale at which it was observed. When considering the earth surface as a two dimensional restriction of the 3D ambient space, the surface is an embedded submanifold of R 3 . On the other hand, when using maps and piecing the global surface together using the compatibility of the overlapping parts, we take the abstract view using charts and atlases.

Embedded Submanifolds

Arguably the simplest example of a 2-dimensional manifold is the sphere S 2 . Relating to the previous example, when embedded in R 3 , we can view it as an idealized model for the surface of the earth. The sphere with radius 1 can be described as the set of unit vectors in R 3 , i.e. the set

S 2 = {(x 1 , x 2 , x 3 ) ∈ R 3 | (x 1 ) 2 + (x 2 ) 2 + (x 3 ) 2 = 1} .
(1.1)

Notice from the definition of the set that all points of S 2 satisfy the equation (x 1 ) 2 + (x 2 ) 2 + (x 3 ) 2 -1 = 0. We can generalize this way of constructing a manifold to the following definition:

Definition 1 (Embedded manifold). Let F : R k → R m be a differentiable map such that the Jacobian matrix dF(x) = ( ∂ ∂x j F i (x)) i j has constant rank kd for all x ∈ F -1 (0). Then the zero-level set M = F -1 (0) is an embedded manifold of dimension d. R k → R m . Here, F : R 3 → R is given by the sphere equation x → (x 1 ) 2 + (x 2 ) 2 + (x 3 ) 2 -1, and the manifold M = S 2 is of dimension 3 -1 = 2.

The map F is said to give an implicit representation of the manifold. In the above example, we used the definition with F(x) = (x 1 ) 2 + (x 2 ) 2 + (x 3 ) 2 -1.

The fact that M = F -1 (0) is a manifold is often taken as the consequence of the submersion level set theorem instead of a definition. The theorem states that with the above assumptions, M has a manifold structure as constructed with charts and atlases. In addition, the topological and differentiable structure of M is in a certain way compatible with that of R k letting us denote M as embedded in R k . For now, we will be somewhat relaxed about the details and use the construction as a working definition of what we think of as a manifold.

The map F can be seen as a set of m constraints that points in M must satisfy. The Jacobian matrix dF(x) at a point in x ∈ M linearizes the constraints around x and its rank kd indicates how many of them are linearly independent. In addition to the unit length constraints of vectors in R 3 defining S 2 , additional examples of commonly occurring manifolds that we will see in this book arise directly from embedded manifolds or as quotients of embedded manifolds.

Example 1. d-dimensional spheres S d embedded in R d+1 . Here we express the unit length equation generalizing (1.1) by

S d = {x ∈ R n+1 | x 2 -1 = 0} .
(1.2)

The squared norm x 2 is the standard squared Euclidean norm on R d+1 .

Example 2. Orthogonal matrices O(k) on R k have the property that the inner products U i , U j of columns U i , U j of the matrix U ∈ M (k,k) vanish for i j and equal 1 for i = j. This gives k 2 constraints, and O(k) is thus an embedded manifold in M (k,k) by the equation

O(k) = U ∈ M (k,k) UU -Id k = 0 (1.3)
with Id k being the identity matrix on R k We will see in Section 1.7.3 that the rank of the map

F(U) = UU -Id k is k(k+1)
2 on O(k), and it follows that O(k) has dimension

k(k-1)
2 .

Charts and Local Euclideaness

We now describe how charts, local parametrizations of the manifold, can be constructed from the implicit representation above. We will use this to give the more abstract definition of a differentiable manifold. When navigating the surface of the earth, we seldom use curved representations of the surface but instead rely on charts that give a flat, 2D representation of regions limited in extent. It turns out that this analogy can be extended to embedded manifolds with a rigorous mathematical formulation.

Definition 2. A chart on the d-dimensional manifold M is a diffeomorphic mapping φ : U → Ũ from an open set U ⊂ M to an open set Ũ ⊆ R d .
The definition exactly captures the informal idea of representing a local part of the surface, the open set U, with a mapping to a Euclidean space, in the surface case R 2 .

When using charts, we often say that we work in coordinates. Instead of accessing points on M directly, we take a chart φ : U → Ũ and use points in φ(U) ⊂ R d instead. This gives us the convenience of having a coordinate system present. However, we need to be aware that the choice of the coordinate system affects the analysis, both theoretically and computationally. When we say that we work in coordinates x = (x 1 , . . . , x d ), we implicitly imply that there is a chart φ such that φ -1 (x) is a point on M.

φ ψ U Ũ Ṽ V M R d R d
It is a consequence of the implicit function theorem that embedded manifolds have charts. Proving it takes some work but we can sketch the idea in the case of the implicit representation map F : R k → R m having Jacobian with full rank m. Recall the setting of the implicit function theorem (see e.g. [START_REF] John | Introduction to smooth manifolds[END_REF]): Let F : R d+m → R m be continuously differentiable and write (x, y) ∈ R d+m such that x denotes the first d coordinates and y the last m coordinates. Let d y F denote the last m columns of the Jacobian matrix dF, i.e. derivatives of F taken with respect to variations in y. If d y F has full rank m at a point (x, y) where F(x, y) = 0, then there exists an open neighborhood Ũ ⊆ R d of x and a differentiable map g : Ũ → R m such that F(x, g(x)) = 0 for all x ∈ Ũ.

The only obstruction to using the implicit function theorem directly to find charts is that we may need to rotate the coordinates on R d+m in order to find coordinates (x, y) and a submatrix d y F with full rank. With this in mind, the map g ensures that F(x, g(x)) = 0 for all x ∈ Ũ, i.e. the points (x, g(x)), x ∈ Ũ are included in M. Setting U = g( Ũ), we get a chart φ : U → Ũ directly by the mapping (x, g(x)) → x.

Abstract Manifolds and Atlases

We now use the concept of charts to define atlases as collections of charts and from this the abstract notion of a manifold.

Definition 3 (Atlas). An atlas of a differentiable set M is a family of charts (φ i ) i=1,...,N with φ i : U i → Ũi such that • φ i cover M: For each x ∈ M, there exists i ∈ {1, . . . , N} such that x ∈ U i , • φ i are compatible: For each pair i, j ∈ {1, . . . , N} where U i ∩ U j is non-empty, the composition

φ i • φ -1 j : φ j (U i ∩ U j ) → R d is a differentiable map.
An atlas thus ensures the existence of at least one chart covering a neighborhood of each point of M. This allows the topological and differential structure of M to be given by definition from the topology and differential structure of the image of the charts, i.e. R d . Intuitively, the structure coming from the Euclidean spaces R d is pulled back using φ i to the manifold. In order for this construction to work, we must ensure that there is no ambiguity in the structure we get if the domain of multiple charts cover a given point. The compatibility condition ensures exactly that.

Definition 4 (Manifold). Let M be a set with an atlas (φ i ) i=1,...,N with φ i :

U i → Ũi , Ũi ⊆ R d . Then M is a manifold of dimension d.
Remark 1. Until now, we have been somewhat loose in describing maps as being "differentiable". The differentiability of maps on a manifold comes from the differential structure which in turn is defined from the atlas and the charts mapping to R d . The differential structure on R d allows derivatives up to any order but the charts may not support this when transferring the structure to M. To be more precise, in the compatibility condition, we require the compositions φ i • φ -1 j to be C r as maps from R d to R d for some integer r. This gives a differentiable structure on M of the same order. In particular, when r ≥ 1, we say that M is a differentiable manifold, and M is smooth if r = ∞. One may also require only r = 0 in which case M is a topological manifold with no differentiable structure.

Because of the implicit function theorem, embedded submanifolds in the sense of Definition 1 have charts and atlases. Embedded submanifolds are therefore particular examples of abstract manifolds. In fact, this goes both ways: The Whitney embedding theorem states that any d-dimensional manifold can be embedded in R k with k ≤ 2d so that the topology is induced by the one of the embedding space. For Riemannian manifolds defined below, this theorem only provides a local C 1 embedding and not a global smooth embedding.

Example 3. The projective space P d is the set of lines through the origin in R d+1 . Each such line intersects the sphere S d in two points that are antipodal. By identifying such points, expressed by taking the quotient using the equivalence relation x ∼ -x, we get the representation P d S d / ∼. Depending on the properties of the equivalence relation, the quotient space of a manifold may not be a manifold in general (more details will be given in Chapter 9 of this book). In the case of the projective space, we can verify the above abstract manifold definition. Therefore, the projective space cannot be seen as an embedded manifold directly, but it can be seen as the quotient space of an embedded manifold.

Tangent Vectors and Tangent Space

As the name implies, derivatives lies at the core of differential geometry. The differentiable structure allows taking derivatives of curves in much the same way as the usual derivatives in Euclidean space. However, spaces of tangent vectors to curves behave somewhat differently on manifolds due to the lack of the global reference frame that the Euclidean space coordinate system gives. We here discuss derivatives of curves, tangent vectors, and tangent spaces.

[0,T]

γ φ R k R d
• γ in R d . If M is embedded, γ is in addition a curve in R k .
The derivative γ(t) is a tangent vector in the linear tangent space T γ(t) M. It can be written in coordinates using φ as γ = γi ∂ x i . In the embedding space, the tangent space T γ(t) M is the affine d-dimensional subspace γ(t) + ker dF(γ(t)) of R k .

Let γ : [0, T ] → R k be a differentiable curve in R k parametrized on the interval [0, T ]. For each t, the curve derivative is

d dt γ(t) = γ =             d dt γ 1 (t) . . . d dt γ k (t)             . (1.4)
This tangent or velocity vector can be regarded as a vector in R k , denoted the tangent vector to γ at t. If M is an embedded manifold in R k and γ(t) ∈ M for all t ∈ [0, T ], we can regard γ as a curve in M. As illustrated on Figure 1.3, the tangent vectors of γ are also tangential to M itself. The set of tangent vectors to all curves at p = γ(t) span a d-dimensional affine subspace of R k that approximates M to the first order at x. This affine space has an explicit realization as x + ker dF(x) where x = γ(t) is the foot-point and ker dF denotes the kernel (null-space) of the Jacobian matrix of F. The space is called the tangent space T x M of M at the point x. In the embedded manifold case, tangent vectors thus arise from the standard curve derivative, and tangent spaces are affine subspaces of R k . On abstract manifolds, the definition of tangent vectors becomes somewhat more intricate. Let γ be a curve in the abstract manifold M, and consider t ∈ [0, T ]. By the covering assumption on the atlas, there exists a chart φ : U → Ũ with γ(t) ∈ U. By the continuity of γ and openness of U, γ(s) ∈ U for s sufficiently close to t. Now the curve γ = φ • γ in R d is defined for such s. Thus we can take the standard Euclidean derivative γ of γ. This gives a vector in R d . In the same way as we define the differentiable structure on M by definition to be that inherited from the charts, it would be natural to let a tangent vector of M be γ by definition. However, we would like to be able to define tangent vectors independently of the underlying curve. In addition, we need to ensure that the construction does not depend on the chart φ.

One approach is to define tangent vectors from their actions on real valued functions on M.

Let f : M → R be a differentiable function. f • γ is then a function from R to R whose derivative is d dt f • γ(t) . (1.5) This operation is clearly linear in f in the sense that d dt ((α f + βg) • γ) = α d dt ( f • γ) + β d dt (g • γ)
when g is another differentiable function and α, β ∈ R. In addition, this derivative satisfies the usual product rule for the derivative of the pointwise product f • g of f and g. Operators on differentiable functions satisfying these properties are called derivations, and we can define tangent vectors and tangent spaces as the set of derivations. I.e. v ∈ T x M is a tangent vector if it defines a derivation v( f ) on functions f ∈ C 1 (M, R). It can now be checked that the curve derivative using a chart above defines derivations. By the chain rule, is can be seen that these derivations are independent of the chosen chart.

The construction of T x M as derivations is rather abstract. In practice, it is often most convenient to just remember that there is an abstract definition and otherwise think of tangent vectors as derivatives of curves. In fact, tangent vectors and tangent spaces can also be defined without derivations using only the derivatives of curves. However, in this case one must define a tangent vector as an equivalence class of curves because multiple curves can result in the same derivative. This construction, although in some sense more intuitive, therefore has its own complexities.

The set {T x M | x ∈ M} has a structure of a differentiable manifold in itself. It is called the tangent bundle T M. It follows that tangent vectors v ∈ T x M for some x ∈ M are also elements of T M. T M is a special case of a fiber bundle (a local product of spaces whose global topology may be more complex). We will later see other examples of fiber bundles, e.g. the cotangent bundle T * M and the frame bundle FM.

A local coordinate system x = (x 1 , . . . x d ) coming from a chart induces a basis

∂ x = (∂ x 1 , . . . ∂ x d ) of the tangent space T x M. Therefore, any v ∈ T x M can be expressed as a linear combination of ∂ x 1 , . . . ∂ x d . Writing v i for the i th entry of such linear combi- nations, we have v = d i=1 v i ∂ x i .
Remark 2 (Einstein summation convention). We will often use the Einstein summation convention that dictates an implicit sum over indices appearing twice in lower and upper position in expressions, in particular in coordinate expressions and tensor calculations. For example, in the coordinate basis mentioned above, we have v = v i ∂ x i where the sum d i=1 is implicit because the index i appears in upper position on v i and lower position on ∂ x i . Just as a Euclidean vector space V has a dual vector space V * consisting of linear functionals ξ : V → R, the tangent spaces T x M and tangent bundle T M have dual spaces, the cotangent spaces T *

x M and cotangent bundle T * M. For each x, elements of the cotangent space T *

x M are linear maps from T x M to R. The coordinate basis (∂ x 1 , . . . ∂ x d ) induces a similar coordinate basis (dx 1 , . . . dx d ) for the cotangent space. This basis is defined from evalution on ∂ x i by dx j (∂ x i ) = δ i j where the delta-function δ i j is 1 if if i = j and 0 otherwise. The coordinates v i for tangent vectors in the coordinate basis had upper indices above. Similarly, coordinates for cotangent vectors conventionally have lower indices such that ξ = ξ i dx i for ξ ∈ T * x M again using the Einstein summation convention. Elements of T * M are called covectors. The evaluation ξ(v) of a covector ξ on a vector v is sometimes written (ξ|v) or ξ, v . Note that the latter notation with brackets is similar to the notation for inner products used below.

Differentials and Pushforward

The interpretation of tangent vectors as derivations above allows taking derivatives of functions. If X is a vector field on M, we can use this pointwise to define a new function on M by taking derivatives at each point, i.e. X( f )(x) = X(x)( f ) using that X(x) is a tangent vector in T x M and hence a derivation that acts on functions. If instead f is a map between two manifolds f : M → N, we get the differential d f : T M → T N as a map between the tangent bundle of M and N. In coordinates, this is d f (∂ x i ) j = ∂ x i f j with f j being the j th component of f . The differential d f is often denoted the pushforward of f because it uses f to map, i.e. push, tangent vectors in T M to tangent vectors in T N. For this reason, the pushforward notation f * = d f is often used. When f is invertible, there exists a corresponding pullback operation

f * = d f -1 .
As a special case, consider a map f between M and the manifold R. Then f * = d f is a map from T M to T R. Because R is Euclidean, we can identify the tangent bundle with R itself, and we can consider

d f a map T M → R. Being a derivative, d f | T x M is linear for each x ∈ M and d f (x) is therefore a covector in T *
x M. Though the differential d f is also a pushforward, the notation d f is most often used because of its interpretation as a covector field.

Riemannian Manifolds

So far, we defined manifolds as having topological and differential structure, either inherited from R k when considering embedded manifolds, or via charts and atlases with the abstract definition of manifolds. We now start including geometric and metric structures.

The topology determines the local structure of a manifold by specifying the open sets and thereby continuity of curves and functions. The differentiable structure allowed us to define tangent vectors and differentiate functions on the manifold. However, we have not yet defined a notion of how "straight" manifold valued curves are. To obtain such a notion, we need to add a geometric structure, called a connection, which allows to compare neighboring tangent spaces and characterizes the parallelism of vectors at different points. Indeed, differentiating a curve on a manifold gives tangent vectors belonging at each point to a different tangent vector space. In order to compute the second order derivative, the acceleration of the curves, we need a way to map the tangent space at a point to the tangent space at any neighboring point. This is the role of a connection ∇ X Y, which specifies how the vector field Y(x) is derived in the direction of the vector field X(x) (Figure 1.4). In the embedding case, tangent spaces are affine spaces of the embedding vector space, and the simplest way to specify this mapping is through an affine transformation, hence the name affine connection introduced by Cartan [START_REF] Cartan | Sur les variétés à connexion affine et la théorie de la relativité généralisée (première partie)[END_REF]. A connection operator also describes how a vector is To define the acceleration as the difference of neighboring tangent vectors, we need to specify a mapping to connect a tangent space at one point to the tangent spaces at infinitesimally close points. In the embedding case, tangent spaces are affine spaces of the embedding vector space, and the simplest way to specify this mapping is through an affine transformation.

transported from a tangent space to a neighboring one along a given curve. Integrating this transport along the curve specifies the parallel transport along this curve. However, there is usually no global parallelism as in Euclidean space. As a matter of fact, transporting the same vector along two different curves arriving at the same point in general leads to different vectors at the endpoint. This is easily seen on the sphere where traveling from north pole to the equator, then along the equator for 90 degrees and back to north pole turns any tangent vector by 90 degrees. This defect of global parallelism is the sign of curvature.

By looking for curves that remain locally parallel to themselves, i.e. such that ∇ γ γ = 0, one defines the equivalent of "straight lines" in the manifold: geodesics. One should notice that there exists many different choices of connections on a given manifold which lead to different geodesics. However, geodesics by themselves do not quantify how far away from each other two points are. For that purpose, we need an additional structure: a distance. By restricting to distances which are compatible with the differential structure, we enter into the realm of Riemannian geometry.

Riemannian Metric

A Riemannian metric is defined by a smoothly varying collection of scalar products ., . x on each tangent space T x M at points x of the manifold. For each x, each such scalar product is a positive definite bilinear map ., . x :

T x M × T x M → R, see Fig- ure 1.5. The inner product gives a norm • x : T x M → R by v 2 = v, v x .
In a given chart, we can express the metric by a symmetric positive definite matrix g(x). The i j th entry of the matrix is denoted g i j (x) and given by the dot product of the coordinate basis for the tangent space: g i j (x) = ∂ x i , ∂ x j x . This matrix is called the local representation of the Riemannian metric in the chart x and the dot products of two vectors v and w in T x M is now in coordinates v, w x = v g(x)w = v i g i j (x)w j . The components g i j of the inverse g(x) -1 of the metric defines a metric on covectors by ξ, η x = ξ i g i j η j . Notice how the upper indices of g i j fit the lower indices of the covector in the Einstein summation convention. This inner product on T *

x M is called a cometric.

Curve Length and Riemannian Distance

If we consider a curve γ(t) on the manifold, we can compute at each abscissa t its velocity vector γ(t) and its norm γ(t) , the instantaneous speed. For the velocity vector, we only need the differential structure, but for the norm, we need the Riemannian metric at the point γ(t). To compute the length of the curve, the norm is integrated along the curve:

L(γ) = γ(t) γ(t) dt = γ(t), γ(t) γ(t) 1 2 dt .
The integrals here are over the domain of the curve, e.g. [0, T ]. We write L b a (γ) = b a γ(t) γ(t) dt to be explicit about the integration domain. This gives the length of the curve segment γ(a) to γ(b).

The distance between two points of a connected Riemannian manifold is the minimum length among the curves γ joining these points:

dist(x, y) = min γ(0)=x,γ(1)=y L(γ). (1.6)
The topology induced by this Riemannian distance is the original topology of the manifold: open geodesic balls constitute a basis of open sets. The Riemannian metric is the intrinsic way of measuring length on a manifold. The extrinsic way is to consider the manifold as embedded in R k and compute the length of a curve in M as for any curve in R k . In section 1.2.4, we identified the tangent spaces of an embedded manifold with affine subspaces of R k . In this case, the Riemannian metric is the restriction of the dot product on R k to the tangent space at each point of the manifold. Embedded manifolds thus inherit also their geometric structure in the form of the Riemannian metric from the embedding space.

Geodesics

In Riemannian manifolds, locally length-minimizing curves are called metric geodesics. The next subsection will show that these curves are also auto-parallel for a specific connection, so that they are simply called geodesics in general. A curve is locally length minimizing if for all t and sufficiently small s, L t+s t (γ) = dist(γ(t), γ(t + s)). This implies that small segments of the curve realize the Riemannian distance. Finding such curves is complicated by the fact that any time-reparameterization of the curve is authorized. Thus one often defines geodesics as critical points of the energy functional E(γ) = 1 2 γ 2 dt. It turns out that critical points for the energy also optimize the length functional. Moreover, they are parameterized proportionally to their arc length removing the ambiguity of the parametrization.

We now define the Christoffel symbols from the metric g by

Γ k i j = 1 2 g km ∂ x i g jm + ∂ x j g mi -∂ x m g i j .
(1.7)

Using the calculus of variations, it can be shown that the geodesics satisfy the second order differential system γk + Γ k i j γi γ j = 0 (1.8)

We will see the Christoffel symbols again below for the connection.

Levi-Civita Connection

The fundamental theorem of Riemannian geometry states that on any Riemannian manifold, there is a unique connection which is compatible with the metric and which has the property of being torsion-free. This connection is called the Levi-Civita connection. For that choice of connection, shortest curves have zero acceleration and are thus geodesics in the sense of being "straight lines". In the following, we only consider the Levi-Civita connection unless explicitly stated.

The connection allows us to take derivatives of a vector field Y in the direction of another vector field X expressed as ∇ X Y. This is also denoted the covariant derivative of Y along X. The connection is linear in X and obeys the product rule in Y so that

∇ X ( f Y) = X( f )Y + f ∇ X Y for a function f : M → R with X( f )
being the derivative of f in the direction of X using the interpretation of tangent vectors as derivations. In a local coordinate system, we can write the connection explicitly using the Christoffel symbols by

∇ ∂ x i ∂ x j = Γ k i j ∂ x k .
With vector fields X and Y having coordinates X(x) = v i (x)∂ x i and Y(x) = w i (x)∂ x i , we can use this to compute the coordinate expression for derivatives of Y along X:

∇ X Y = ∇ v i ∂ x i (w j ∂ x j ) = v i ∇ ∂ x i (w j ∂ x j ) = v i (∂ x i w j )∂ x j + v i w j ∇ ∂ x i ∂ x j = v i (∂ x i w j )∂ x j + v i w j Γ k i j ∂ x k = v i (∂ x i w k ) + v i w j Γ k i j ∂ x k = X(w k ) + v i w j Γ k i j ∂ x k .
Using this, the connection allows us to write the geodesic equation (1.8) as the zero acceleration constraint:

0 = ∇ γ γ = γ(γ k ) + γi γ j Γ k i j ∂ x k = γk + γi γ j Γ k i j ∂ x k .
The connection also defines the notion of parallel transport along curves. A vector v ∈ T γ(t 0 ) M is parallel transported if it is extended to a t-dependent family of vectors with v t ∈ T γ(t) M and ∇ γ(t) v t = 0 for each t. Parallel transport can thereby be seen as a map P γ,t : T γ(t 0 ) M → T γ(t) M linking tangent spaces. The parallel transport inherits linearity from the connection. It follows from the definition that γ is a geodesic precisely if γ(t) = P γ,t (γ(t 0 )).

It is a fundamental consequence of curvature that parallel transport is dependent on the curve along which the vector is transported: With curvature, the parallel transports P γ,T and P φ,T along two curves γ and φ with the same end-points γ(t 0 ) = φ(t 0 ) and γ(T ) = φ(T ) will differ. The difference is denoted holonomy, and the holonomy of a Riemannian manifold vanishes only if M is flat, i.e. has zero curvature.

Completeness

The Riemannian manifold is said to be geodesically complete if the definition domain of all geodesics can be extended to R. This means that the manifold has no boundary nor any singular point that we can reach in a finite time. For instance, R d -{0} with the usual metric is not geodesically complete because some geodesics will hit 0 and thus stop being defined in finite time. On the other hand, R d is geodesically complete. Other examples of complete Riemannian manifolds include compact manifolds implying that S d is geodesically complete. This is a consequence of the Hopf-Rinow-De Rham theorem that also states that geodesically complete manifolds are complete metric spaces with the induced distance, and that there always exists at least one minimizing geodesic between any two points of the manifold, i.e. a curve whose length is the distance between the two points.

From now on, we will assume that the manifold is geodesically complete. This assumption is one of the fundamental properties ensuring well-posedness of algorithms for computing on manifolds.

Exponential and Logarithm Maps

Let x be a point of the manifold that we consider as a local reference point and v a vector of the tangent space T x M at that point. From the theory of second order differential equations, it can be shown that there exists one and only one geodesic γ (x,v) (t) starting from that point x = γ (x,v) (0) with tangent vector v = γ(x,v) (0). This geodesic is at first defined in a sufficiently small interval around zero but since the manifold is assumed geodesically complete, its definition domain can be extended to R. Thus, the points γ (x,v) (t) are defined for each t and each v ∈ T x M. This allows to map vectors in the tangent space to the manifold using geodesics: the vector v ∈ T x M can be mapped to the point of the manifold that is reached after a unit time t = 1 by the geodesic γ (x,v) (t) starting at x with tangent vector v. This mapping Exp x :

T x M -→ M v -→ Exp x (v) = γ (x,v) (1)
is called the exponential map at point x. Straight lines passing 0 in the tangent space are transformed into geodesics passing the point x on the manifold and distances along these lines are conserved (Fig. 1.6). When the manifold is geodesically complete, the exponential map is defined on the entire tangent space T x M but it is generally one-to-one only locally around 0 in the tangent space corresponding to a local neighborhood of x on M. In the sequel, we denote by -→ xy or Log x (y) the inverse of the exponential map where the inverse is defined: this is the smallest vector as measured by the Riemannian metric such that y = Exp x ( -→ xy). In this chart, the geodesics going through x are represented by the lines going through the origin: Log They have zero acceleration ∇ γ γ and their velocity vectors are parallel transported γ(t) = P γ,t (γ(t 0 )). Geodesics locally realize the Riemannian distance so that dist(x, γ(t)) = t v for sufficiently small t. (right) The tangent space T x S 2 and Exp x give an exponential chart mapping vectors v ∈ T x S 2 to points in S 2 by Exp x (v). The cut locus of x is its antipodal point, and the injectivity radius is π. Note that the equator is the set Exp

x (v) v = π 2 .
Thus, the exponential chart at x gives a local representation of the manifold in the tangent space at a given point. This is also called a normal coordinate system or normal chart if it is provided with an orthonormal basis. At the origin of such a chart, the metric reduces to the identity matrix and the Christoffel symbols vanish. Note again that the exponential map is generally only invertible locally around 0 ∈ T x M and Log x y is therefore only locally defined, i.e. for points y near x.

The exponential and logarithm maps are commonly referred to as the Exp and Logs maps.

Cut Locus

It is natural to search for the maximal domain where the exponential map is a diffeomorphism. If we follow a geodesic γ (x,v) (t) = Exp x (t v) from t = 0 to infinity, it is either always minimizing for all t, or it is minimizing up to a time t 0 < ∞. In this last case, the point z = γ (x,v) (t 0 ) is called a cut point and the corresponding tangent vector t 0 v a tangential cut point. The set of all cut points of all geodesics starting from x is the cut locus C(x) ∈ M, and the set of corresponding vectors the tangential cut locus C(x) ∈ T x M. Thus, we have C(x) = Exp x (C(x)), and the maximal definition domain for the exponential chart is the domain D(x) containing 0 and delimited by the tangential cut locus.

It is easy to see that this domain is connected and star-shaped with respect to the origin of T x M. Its image by the exponential map covers the manifold except the cut locus and the segment [0, -→ xy] is transformed into the unique minimizing geodesic from

x to y. Hence, the exponential chart has a connected and star-shaped definition domain that covers all the manifold except the cut locus C(x):

D(x) ⊆ R d ←→ M -C(x) - → xy = Log x (y) ←→ y = Exp x ( - → xy) .
From a computational point of view, it is often interesting to extend this representation to include the tangential cut locus. However, we have to take care of the multiple representations: Points in the cut locus where several minimizing geodesics meet are represented by several points on the tangential cut locus as the geodesics are starting with different tangent vectors (e.g. antipodal points on the sphere and rotation of π around a given axis for 3D rotations). This multiplicity problem cannot be avoided as the set of such points is dense in the cut locus. The size of D(x) is quantified by the injectivity radius i(M, x) = dist(x, C(x)), which is the maximal radius of centered balls in T x M on which the exponential map is one-to-one. The injectivity radius of the manifold i(M) is the infimum of the injectivity over the manifold. It may be zero, in which case the manifold somehow tends towards a singularity (e.g. think to the surface z = 1/ x 2 + y 2 as a sub-manifold of R 3 ).

Example 4. On the sphere S d (center 0 and radius 1) with the canonical Riemannian metric (induced by the ambient Euclidean space R d+1 ), the geodesics are the great circles and the cut locus of a point x is its antipodal point x = -x. The exponential chart is obtained by rolling the sphere onto its tangent space so that the great circles going through p become lines. The maximal definition domain is thus the open ball D = B d (π). On its boundary ∂D = C = S d-1 (π), all the points represent x, see Figure 1.6.

For the real projective space P d (obtained by identification of antipodal points of the sphere S d ), the geodesics are still the great circles, but the cut locus of the point {x, -x} is now the equator of the two points, with antipodal points identified (thus the cut locus is P d-1 ). The definition domain of the exponential chart is the open ball

D = B d ( π
2 ), and the tangential cut locus is the sphere ∂D = S d-1 ( π 2 ) where antipodal points are identified.

Elements of Analysis in Riemannian Manifolds

We here outline further constructions on manifolds relating to taking derivatives of functions, the intrinsic Riemannian measure, and defining curvature. These notions will be used in following chapters of this book for instance for optimization algo-rithms.

Gradient and Musical Isomorphisms

Let f be a smooth function from M to R. Recall that the differential d f (x) evaluated at the point x ∈ M is a covector in T *

x M. Therefore, contrary to the Euclidean situation where derivatives are often regarded as vectors, we cannot directly interpret d f (x) as a vector. However, thanks to the Riemannian metric, there is a canonical way to identify the linear form d f ∈ T *

x M with a unique vector v ∈ T x M. This is done by defining v ∈ T x M to be a vector satisfying d f (w) = v, w x for all vectors w ∈ T x M. This mapping corresponds to the transpose operator that is implicitly used in Euclidean spaces to transform derivatives of functions (row vectors) to column vectors. On manifolds, the Riemannian metric must be specified explicitly since the coordinate system used may not be orthonormal everywhere.

The mapping works for any covector, and is often denoted the sharp map : T * M → T M. It has an inverse in the flat map : T M → T * M. In coordinates, (ξ ) i = g i j ξ j for a covector ξ = ξ j dx j and (v ) i = g i j v j for a vector v = ∂ x j v j . The maps and are denoted musical isomorphisms because they raise or lower the indices of the coordinates.

We can use the sharp map to define the Riemannian gradient as a vector:

grad f = (d f ) .
This definition corresponds to the classical gradient in R k using the standard Euclidean inner product as Riemannian metric. Using the coordinate representation of the sharp map, we get the coordinate form (grad f ) i = g i j ∂ x j f of the gradient.

Hessian and Taylor Expansion

The covariant derivative of the gradient, the Hessian, arises from the connection ∇:

Hess f (X, Y) = ∇ X ∇ Y f = (∇ X (d f ))Y = ∇ X grad f , Y .
Here, the two expressions on the right are given using the action of the connection on the differential form d f (a covector), or the vector field grad f = (d f ) . Its expression in a local coordinate system is

Hess f = ∇d f = (∂ x i x j f -Γ k i j ∂ k f )dx i dx j .
Let now f x be the expression of f in a normal coordinate system at x. Its Taylor expansion around the origin in coordinates is:

f x (v) = f x (0) + d f x v + 1 2 v H f x v + O( v 3 )
where d f x = (∂ x i f ) is the Jacobian matrix of first order derivatives and H f x = (∂ x i x j f ) is the Euclidean Hessian matrix. Because the coordinate system is normal, we have

f x (v) = f (Exp x (v))
. Moreover, the metric at the origin reduces to the identity: d f x = (grad f ) T , and the Christoffel symbols vanish so that the matrix of second derivatives H f x corresponds to the Hessian Hess f . Thus, The Taylor expansion can be written in any coordinate system:

f Exp x (v) = f (x) + grad f (v) + 1 2 Hess f (v, v) + O( v 3 ) (1.9)

Riemannian Measure or Volume Form

In a vector space with basis A = (a 1 , . . . a n ), the local representation of the metric is given by g = A A where A = [a 1 , . . . a n ] is the matrix of coordinates change from A to an orthonormal basis. Similarly, the measure or the infinitesimal volume element is given by the volume of the parallelepiped spanned by the basis vectors: dV = |A| dx = |g| dx with | • | denoting the matrix determinant. In a Riemannian manifold M, the Riemannian metric g(x) induces an infinitesimal volume element on each tangent space, and thus a measure on the manifold that in coordinates has the expression

dM(x) = |g(x)| dx .
The cut locus has null measure, and we can therefore integrate indifferently in M or in any exponential chart. If f is an integrable function of the manifold and f x (v) = f (Exp x (v)) is its image in the exponential chart at x, we have

x∈M f (x)dM(x) = v∈D(x) f x (v) |g(Exp x (v))| dv ,

Curvature

The curvature of a Riemannian manifold measures its deviance from local flatness. We often have a intuitive notion of when a surface embedded in R 3 is flat or curved: For example, a linear subspace of R 3 is flat while the sphere S 2 is curved. This idea of curvature is expressed in the Gauss curvature. However, for high dimensional spaces, the mathematical description becomes somewhat more intricate. We will see below several notions of curvature capturing aspects of the nonlinearity of the manifold with varying details. Important to note is that while vanishing curvature implies local flatness of the manifold, this is not the same as the manifold being globally Euclidean. An example is the torus T 2 which can both be embedded in R 3 inheriting nonzero curvature and be embedded in R 4 in a way in which it inherits a flat geometry. In both cases, the periodicity of the torus remains, which prevents it from being a vector space.

The curvature of a Riemannian manifold is described by the curvature tensor R : T M × T M × T M → T M. It is defined from the covariant derivative by evaluation on vector fields X, Y, Z:

R(X, Y)Z = ∇ X ∇ Y Z -∇ Y ∇ X Z -∇ [X,Y] Z .
(1.10)

The bracket [X, Y] denotes the anticommutativity of the fields X and Y. If f is a differentiable function on M, the new vector field produced by the bracket is given by its application to f : [X, Y] f = X(Y( f )) -Y(X( f )). The curvature tensor R can intuitively be interpreted at x ∈ M as the difference between parallel transporting the vector Z(x) along an infinitesimal parallelogram with sides gives by X(x) and Y(x), see Figure 1.7 (left). As noted earlier, parallel transport is curve dependent, and the difference between transporting infinitesimally along X(x) and then Y(x) as opposed to along Y(x) and then X(x) is a vector in T x M. This difference can be calculated for any vector z ∈ T x M. The curvature tensor when evaluated at X, Y, i.e. R(X, Y), is the linear map T x M → T x M given by this difference. The reader should note that two different sign conventions exist for the curvature tensor: the definition (1.10) is used in a number of reference books in physics and mathematics [START_REF] Misner | Gravitation[END_REF][START_REF] John | Riemannian manifolds[END_REF][START_REF] Klingenberg | Riemannian Geometry[END_REF][START_REF] Mikhailovich | Geometry VI: Riemannian geometry[END_REF][START_REF] Jost | Riemannian Geometry and Geometric Analysis[END_REF]. Other authors use a minus sign to simplify some of the tensor notations [Spi79, O'N83, dC92, GHL93, Ber03] and different order conventions for the tensors subscripts and/or a minus sign in the sectional curvature defined below (see e.g. the discussion in [GQ16, p. 399]).

The curvature can be realized in coordinates from the Christoffel symbols:

R(∂ x i , ∂ x j ) ∂ x k = R m ki j ∂ x m = (Γ l jk Γ m il -Γ l ik Γ m jl + ∂ i Γ m jk -∂ j Γ m ik ) ∂ x m .
(1.11)

The sectional curvature κ measures the Gaussian curvature of 2D submanifolds of M, the Gaussian curvature of each point being the product of the principal curvatures of curves passing the point. The 2D manifolds arise as the Gaussian curvature of the geodesic spray of a 2D linear subspace of T x M, see Figure 1.7 (right). Such a 2-plane can be represented by basis vectors u, v ∈ T x M in which case the sectional curvature can be expressed using the curvature tensor by

κ(u, v) = R(u, v)v, u u 2 v 2 -u, v 2 .
(1.12)

The curvature tensor gives the notion of Ricci and scalar curvature that both provides summary information of the full tensor R. The Ricci curvature Ric is the trace over the first and last indices of R with coordinate expression

R i j = R k ki j = R k ik j = g kl R ik jl . (1.13)
Taking another trace, we get the scalar valued quantity the scalar curvature S :

S = g i j R i j . (1.14)
Note that the cometric appears to raise one index before taking the trace.

Lie groups and Homogeneous Manifolds

A Lie group is a manifold equipped with additional group structure such that the group multiplication and group inverse are smooth mappings. Many of the interesting transformations used in image analysis, translations, rotations, affine transforms, etc., form Lie groups. We will in addition see examples of infinite dimensional Lie groups when doing shape analysis with diffeomorphisms as described in chapter 4. We begin by reviewing the definition of an algebraic group.

Definition 5 (Group). A group is a set G with a binary operation, denoted here by concatenation or group product, such that 1. (xy)z = x(yz), for all x, y, z ∈ G, 2. there is an identity element e ∈ G satisfying xe = ex = x, for all x ∈ G,

3. each x ∈ G has an inverse, x -1 ∈ G, satisfying xx -1 = x -1 x = e.
A Lie group is simultaneously a group and a manifold, with compatibility between these two mathematical concepts.

Definition 6 (Lie Group). A Lie group G is a smooth manifold that also forms a group, where the two group operations,

(x, y) → xy : G × G → G product x → x -1 : G → G inverse
are smooth mappings of manifolds.

Example 5. The space of all k × k non-singular matrices forms a Lie group called the general linear group, denoted GL(k). The group operation is matrix multiplication, and GL(k) can be given a smooth manifold structure as an open subset of R k 2 . The equations for matrix multiplication and inverse are smooth operations in the entries of the matrices. Thus, GL(k) satisfies the requirements of a Lie group in Definition 6. A matrix group is any closed subgroup of GL(k). Matrix groups inherit the smooth structure of GL(k) as a subset of R k 2 and are thus also Lie groups.

Example 6. The k × k rotation matrices form a closed matrix subgroup of GL(k) and thus a Lie group. This group is called the special orthogonal group. It is defined as

SO(k) = {R ∈ GL(k) : R T R = Id k and det(R) = 1}.
This space is a closed and bounded subset of R k 2 and thus compact.

Example 7. Classical geometric transformation groups used in image registration such as rigid-body transformations, similarities and affine transformations can also be looked upon as matrix groups via their faithful representation based on homogeneous coordinates.

For each y in a Lie group G, the following two diffeomorphisms of G are denoted left-and right-translations by y:

L y : x → yx (left multiplication) R y : x → xy (right multiplication)
The differential or pushforward (L y ) * of the left translation maps the tangent space T x G to the tangent space T yx G. In particular, (L y ) * maps any vector u ∈ T e G to the vector (L y ) * u ∈ T y G thereby giving rise to the vector field ũ(y) = (L y ) * u. Such a vector field is said to be left-invariant since it is invariant under left multiplication: ũ • L y = (L y ) * ũ = ũ for every y ∈ G. Right-invariant vector fields are defined similarly. A leftor right-invariant vector field is uniquely defined by its value on the tangent space at the identity, T e G.

Recall that vector fields on G can be seen as derivations on the space of smooth functions, C ∞ (G). Thus, two vector fields u and v can be composed to form another operator uv on C ∞ (G) but the operator uv is not necessarily a derivation as it includes second order differential terms. However, the operator uvvu is a vector field on G. Indeed, we can check by writting this expression in a local coordinate system that the second order terms vanish. This leads to a definition of the Lie bracket of vector fields u, v on G, defined as

[u, v] = uv -vu.
(1.15) This is also sometimes called the Lie derivative

L u v = [u, v]
because it is conceptually the derivative of the vector field v in the direction u(x) generated by u at each point x ∈ G.

Definition 7 (Lie algebra). A Lie algebra is a vector space V equipped with a bilinear product

[•, •] : V × V → V, called a Lie bracket, that satisfies 1. [u, v] = -[v, u] (skew symmetry), 2. [[u, v], w] + [[v, w], u] + [[w, u], v] = 0 (Jacobi identity), for all u, v, w ∈ V.
The tangent space of a Lie group G at the identity element, T e G, typically denoted g, forms a Lie algebra. The Lie bracket on g is induced by the Lie bracket on the corresponding left-invariant vector fields. If u, v are two vectors in g, let ũ, ṽ be the corresponding unique left-invariant vector fields on G. Then the Lie bracket on g is given by [u, v] = [ũ, ṽ](e).

The Lie bracket provides a test for whether the Lie group G is commutative. A Lie group G is commutative if and only if the Lie bracket on the corresponding Lie algebra g is zero, i.e., [u, v] = 0 for all u, v ∈ g.

Example 8. The Lie algebra for Euclidean space R k is again R k . The Lie bracket is zero, i.e., [X, Y] = 0 for all X, Y ∈ R k .

Example 9. The Lie algebra for GL(k) is gl(k), the space of all real k × k matrices. The Lie bracket operation for X, Y ∈ gl(k) is given by

[X, Y] = XY -Y X.
Here the product XY denotes actual matrix multiplication, which turns out to be the same as composition of the vector field operators (compare to (1.15)). All Lie algebras corresponding to matrix groups are subalgebras of gl(k).

Example 10. The Lie algebra for the rotation group SO(k) is so(k), the space of skewsymmetric matrices. A matrix A is skew-symmetric if A = -A T .

One-parameter Subgroups

Let ũ(y) = (L y ) * u be a left-invariant vector field. The solution x(t) to the initial value problem ẋ(t) = ũ(x(t)), x(0) = e is called a one-parameter subgroup because it is a morphing of Lie groups: x(s + t) = x(s)x(t). The Lie exponential map exp : g → G is then given by the value of x(t) at t = 1, i.e. exp(u) = x(1). For matrix groups where the Lie group algebra consists of ordinary matrices, exp corresponds to the matrix exponential. The group exponential should not be confused with the Riemannian exponential as they usually differ, unless the group is provided with a bi-invariant metric.

Actions

Let M be a manifold and G a Lie group. The elements of the group can often be used to produce variations of elements of the manifold, e.g. elements of GL(k) linearly transforms elements of the manifold R k . Similarly, affine transformations applies to change images in image registration. These are examples of actions of G on M. Such actions are usually denoted g. x where g ∈ G, x ∈ M. Because the action involves two manifolds, G and M, we will below use x, y to denote elements of M, and g, h to denote elements of G.

Definition 8 (Action). A left action of a Lie group G on a manifold M is a smooth mapping . : G × M → M satisfying 1. e. x = x, ∀x ∈ M, 2. h.(g. x) = (hg). x, ∀x ∈ M, 3. the map x → g. x is a diffeomorphism of M for each g ∈ G.

We will see examples of Lie group actions throughout the book. For example, Chapter 4 on shape analysis relies fundamentally on actions of the group Diff(Ω) of diffeomorphisms of a domain Ω on shape spaces S.

Through the action, a curve g(t) on the group G acts on a point x ∈ M to give a curve g(t). x in M. In particular, one-parameter subgroups define the curves x v (t) = exp(tv). x in M for Lie algebra elements v ∈ g. The derivative

v M (x) := d dt exp(tv). x in T x v (t)
M is denoted the infinitesimal generator associated to v. Some particularly important actions are the actions of G on itself and on the Lie algebra g. These include the actions by left translation g. h := L g (h) = gh, and the action by conjugation g.h := L g (R g -1 )h = ghg -1 . The pushforward of the conjugation gives the adjoint action g.v = (L g • R g -1 ) * v of G on g. The adjoint action is also denoted Ad g v.

From the adjoint action, we get the adjoint operator ad u v = d dt Ad exp(tu) v for v, u ∈ g. This operator is sometimes informally denoted "little ad", and it is related to the Lie bracket by ad

u v = [u, v].
The actions on the Lie algebra have dual actions as well, denoted co-actions: the co-adjoint action G × g * → g * , g. ξ := Ad * g -1 ξ for ξ ∈ g * where the dual of the adjoint is given by Ad * g -1 ξ(v) = ξ(Ad g -1 v) for all v ∈ g. Using the notation (ξ|v) for evaluation ξ(v) of ξ on v, the definition of the dual of the adjoint is (Ad * g -1 ξ|v) = (ξ|Ad g -1 v). The coadjoint operator ad * : g × g * → g * is similarly specified by (ad * v ξ|u) = (ξ|ad v u) for v, u ∈ g and ξ ∈ g * .

Homogeneous Spaces

Let the group G act on M. If, for any x, y ∈ M, there exists g ∈ G such that g. x = y, the action is said to be transitive. In this case, the manifold M is homogeneous. For a fixed x ∈ M, the closed subgroup H = {g ∈ g | g. x = x} is denoted the isotropy subgroup of G, and M is isomorphic to the quotient G/H. Similarly, a closed subgroup H of G leads to a homogeneous space G/H by quotienting out H. Examples of homogeneous spaces are the spheres S n ≡ SO(n + 1)/SO(n) and the orbit shape spaces described in Chapter 4, e.g. the manifold of landmark configurations.

Invariant Metrics and Geodesics

The left-and right translation maps gives a particularly useful way of defining Riemannian metrics on Lie groups. Given an inner product •, • g on the Lie algebra, we can extend it to an inner product on tangent spaces at all elements of the group by setting u, v g := (L g -1 ) * u, (L g -1 ) * v g . This defines a left-invariant Riemannian metric on G because, (L h ) * u, (L h ) * v hg = u, v g for any u, v ∈ T g G. Similarly, we can set u, v g := (R g -1 ) * u, (R g -1 ) * v g to get a right-invariant metric. In the particular case where the metric is invariant to both left-and right-translation, it called bi-invariant.

Geodesics for bi-invariant metrics are precisely one-parameter subgroups, and the Lie group exponential map exp therefore equals the Riemannian exponential map Exp e . For metrics that are left-or right-invariant, but not bi-invariant, the ordinary geodesic equation (1.8) can be simplified using for example Euler-Poincaré reduction. The resulting Euler-Poincaré equations are discussed further in chapter 4 in the case of right-invariant metrics on the group Diff(Ω).

Elements of Computing on Riemannian Manifolds

The Riemannian Exp and Log maps constitute very powerful atomic functions to express most geometric operations for performing statistical computing on manifolds. The implementation of Log x and Exp x is therefore the algorithmic basis of programming on Riemannian manifolds, as we will see in the following.

In a Euclidean space, exponential charts are nothing but orthonormal coordinates systems translated to each point: In this case, -→ xy = Log x (y) = yx and Exp x (v) =

x + v. This example is more than a simple coincidence. In fact, most of the usual operations using additions and subtractions may be reinterpreted in a Riemannian framework using the notion of bipoint, an antecedent of vector introduced during the 19th century. Indeed, one defines vectors as equivalent classes of bipoints, oriented couples of points, in a Euclidean space. This is possible using the canonical way to compare what happens at two different points by translating. In a Riemannian manifold, we can compare vectors using the parallel transport along curves but the curve dependence of the parallel transport prevents global comparison of vectors as in Euclidean space. This implies that each vector has to remember at which point of the manifold it is attached, as is the case for tangent vectors, which relates back to the Euclidean notion of a bipoint. Conversely, the logarithm map may be used to map almost any bipoint (x, y) into a vector -→ xy = Log x (y) of T x M. This reinterpretation of addition and subtraction using logarithm and exponential maps is very powerful when generalizing algorithms working on vector spaces to algorithms on Riemannian manifolds. This is illustrated in Table 1.1 and in the following sections. The Exp and Log maps are different for each manifold and for each metric. They must therefore be determined and implemented on a case by case basis. In some cases, closed form expressions are known, examples being the spheres S d , rotations and rigid body transformations with left invariant metric [START_REF] Pennec | A Framework for Uncertainty and Validation of 3D Registration Methods based on Points and Frames[END_REF], and covariance matrices (positive definite symmetric matrices, so called tensors in medical image analysis) [START_REF] Pennec | A Riemannian Framework for Tensor Computing[END_REF] and Chapter 3. In cases where closed form solutions are not known, geodesics with given initial velocity can be obtained by numerically solving the geodesic ODE (1.8), or by solving the variational problem of finding a minimum energy curve between two points. Thus computing Exp x (v) may be posed as a numerical integration problem, see e.g. [START_REF] Helmke | Optimization and Dynamical Systems[END_REF][START_REF] Hairer | Geometric numerical integration : structure preserving algorithm for ordinary differential equations[END_REF], and computing -→ xy = Log x (y) an optimal control problem. This opens the way to statistical computing in complex spaces than the spaces we have considered up to now, such as spaces of curves, surfaces, and diffeomorphic transformations, as we will see in the following chapters. Geometric computation frameworks such as Theano Geometry1 [KAS17] and Geomstats2 provide numerical implementations of geometric operations on some commonly used manifolds. Theano Geometry uses automatic differentiation to express and compute the derivatives that are essential for differential geometric computations. This results in convenient code for computing Christoffel symbols, curvature tensors, and fiber bundle operations using the parallel transport.

Euclidean space Riemannian manifold Subtraction

- → xy = y -x - → xy = Log x (y) Addition y = x + v y = Exp x (v) Distance dist(x, y) = y -x dist(x, y) = - → xy x Mean value (implicit) i (x i -x) = 0 i -→ xx i = 0 Gradient descent x t+ε = x t -ε∇ f (x t ) x t+ε = Exp x t (-ε grad f (x t )) Geodesic interpolation x(t) = x 0 + t ---→ x 0 x 1 x(t) = Exp x 0 (t ---→ x 0 x 1 )

Examples

Below, we survey ways to express the Exp and Log maps on selected manifolds at the same time exemplifying how particular structure of the spaces can be used for computations.

The Sphere

Let x be a point in S d . From the embedding of S d in R d+1 , the tangent space T x S d can be identified with the d-dimensional vector space of all vectors in R d+1 orthogonal to x. The inner product between two tangent vectors is then equivalent to the usual Euclidean inner product. The exponential map is given by a 2D rotation of x by an angle given by the norm of the tangent, i.e.,

Exp x (v) = cos θ x + sin θ θ v, θ = v . (1.16)
The log map between two points x, y on the sphere can be computed by finding the initial velocity of the rotation between the two points. Let π x (y) = x x, y denote the identity matrix Id k as neutral element e:

O(k) = U ∈ M (k,k) | U U = Id k .
Because the quadratic constraint is smooth and differentiable, O(k) constitute a Lie group, submanifold of the linear space of square matrices. However, it is not connected: Taking the determinant of the constraint gives det(U) 2 = det( Id k ), so that det(U) = ±1. We see that there are two disconnected component of determinant +1 and -1 that cannot be joined by any continuous curve on the space of matrices. Such a curve would have to go through matrices with determinants between -1 and 1 since the determinant is a continuous function. The component of negative determinant includes symmetries that reverse the orientation of the space. It is not a subgroup because the composition of two such transformations of negative determinant has a positive determinant.

The component of positive determinant preserves the orientation of the space and is a subgroup, the group of rotations, or special orthogonal transformations:

SO(k) = R ∈ M (k,k) | R R = Id k , det(R) = 1 .
Let R(t) = R + t Ṙ + O(t 2 ) be a curve drawn on SO(k), considered as an embedded manifold in the vector space of matrices M (k,k) . The constraint R R = Id k is differentiated into ṘR + ṘR = 0 or R Ṙ + R Ṙ = 0 which means that ṘR and R Ṙ are skew-symmetric matrices. Thus, the tangent space T e SO(k) at identity is the vector space of skew matrices and the tangent space at rotation R ∈ SO(k) is its left or right translation:

T R SO(k) = {X ∈ M(k, k) / R X = -(R X) } = {X ∈ M(k, k) / XR = -(XR ) }
Since k × k skew-symmetric matrices have k(k -1)/2 free components, we also obtain that the dimension of the special orthogonal group is k(k -1)/2. In order to put a metric on this Lie group, we may take a metric on the tangent space at the identity and left translate it to any other point resulting in a left-invariant metric. We may similarly right translate it to obtain a right-invariant metric. Since SO(k) is a submanifold of the Euclidean space of matrices M (k,k) , we may also consider the restriction of the embedding Frobenius dot product Tr(X Y) to the tangent spaces at all points. It is common to rescale the Frobenius metric by 1/2 to compensate the fact that we are counting twice each off diagonal coefficient of the skew-symmetric matrices. This induces the metric

X, Y R = 1 2 Tr(XY )
on the tangent space T R SO(k). This metric is invariant by left and right translation. The existence of this biinvariant metric is a particular case due to the compactness of the SO(k) group. Biinvariant metrics on Lie groups have very special properties that will be described in Chapter 5. In particular, as mentioned earlier, geodesics passing the identity are oneparameter subgroups whose equations are given by the matrix exponential: Exp e (X) = exp(X) = +∞ k=0 X k k! . This series is absolutely convergent so that the matrix exponential always exists. Its inverse, the logarithm, may however fail to exist and is also generally not unique when it exists.

For rotations, the exponential of skew symmetric matrices covers the whole rotation group so that the log always exists, but it is not unique: For k = 2, rotating of an angle θ is the same as rotating of an angle θ + 2lπ where l is an integer. In order to understand the structure of rotations in higher dimensions, we may look at the spectral decomposition of a rotation matrix R: The characteristic polynomial P(λ) = det(R -λ Id k ) is a real polynomial of degree k. Thus the k complex eigenvalues are real or conjugate by pairs, and the polynomial can be factored into at most k/2 quadratic terms, potentially with multiplicity, and real linear terms. The conservation of the norm by the rotation Rx = x shows that the modulus of all the eigenvalues is 1. Thus, eigenvalues are e ±iθ j or 1. Since a rotation is a normal matrix, it can be diagonalized, and we conclude that every rotation matrix, when expressed in a suitable coordinate system, partitions into k/2 independent 2D rotations, called Givens rotations [START_REF] Horn | Matrix analysis[END_REF]: R(θ j ) = cos(θ j )sin(θ j ) sin(θ j ) cos(θ j ) = exp θ j 0 -1 1 0 .

Conversely, each skew symmetric matrix Ω = -Ω decomposes the space R k in a direct sum of mutually orthogonal subspaces, which are all invariant under Ω [START_REF] Horn | Matrix analysis[END_REF]. The decomposition has l (possibly equal to zero) 2-dimensional vector subspaces E j on which Ω acts non trivially, and one single subspace F of dimension k -2l, the orthogonal complement of the span of other subspaces, which is the kernel of Ω. For any E j , there exists an orthonormal basis of E j such that Ω restricted to E j is in this basis of the following matrix form: θ j 0 -1 1 0 where θ j ( 0) is the j th angle of rotation of the n-dimensional rotation exp Ω.

We can now come back to the uniqueness of the Log: When the angles of the above k/2 2D rotations decomposing the rotation R are within ] -π, π[ , the logarithm of R is well-defined. Otherwise, one cannot define a unique logarithm. This is only the case for 2D rotations of 180 degrees, whose two 'smallest' real logarithms are the following: 0 -π π 0 and 0 π -π 0 .

Geodesics starting at any point in the group are left, or right, translation of geodesics starting at identity. For instance, γ (R,Y) (t) = R exp(t R Y) is the unique geodesic starting at R with tangent vector Y. The following reasoning underlies this: To find the geodesic starting at R with tangent vector Y, we first left translate Y by R to the tangent space at identity, compute the geodesic starting at e with tangent vector R Y, and left translate back the result by R. Since the metric is bi-invariant, the same mechanism can be implemented with right translation. The formula for the exponential map Exp R : T R SO(k) → SO(k) at point R is thus:

Exp R (X) = R Exp e (R X) = R exp(R X) = exp(XR ) R.
(1.20)

Likewise, to compute the log map of rotation U at rotation R, we first left translate both rotations by R , take the log map of R U at e, and left translate back to result:

Log R (U) = R Log e (R U) = R log(R U) = log(UR ) R.

(1.21)

Additional References

This very compact introduction to differential geometry, Riemannian manifolds, and Lie groups provides only a brief overview of the underlying deep theory. We here provide some references to further reading. There are many excellent texts on the subjects. The lists below are therefore naturally non-exhaustive.

Figure 1 . 1

 11 Figure 1.1 An embedded manifold arises as the zero level subset M = F -1 (0) of the map F :

Figure 1 . 2

 12 Figure 1.2 Charts φ : U → Ũ and ψ : V → Ṽ, members of the atlas covering the manifold M, from the open sets U, V ⊂ M to open sets Ũ, Ṽ of R d , respectively. The compatibility condition ensures that φ and ψ agree on the overlap U ∩ V between U and V in the sense that the composition ψ • φ -1 is a differentiable map.

Figure 1 . 3

 13 Figure 1.3 The curve γ maps the interval [0, T ] to the manifold. Using a chart φ, we can work in coordinates with the curve φ• γ in R d . If M is embedded, γ is in addition a curve in R k .The derivative γ(t) is a tangent vector in the linear tangent space T γ(t) M. It can be written in coordinates using φ as γ = γi ∂ x i . In the embedding space, the tangent space T γ(t) M is the affine d-dimensional subspace γ(t) + ker dF(γ(t)) of R k .

Figure 1 . 4

 14 Figure1.4 Tangent vectors along the red and blue curves drawn on the manifold belong to different tangent spaces. To define the acceleration as the difference of neighboring tangent vectors, we need to specify a mapping to connect a tangent space at one point to the tangent spaces at infinitesimally close points. In the embedding case, tangent spaces are affine spaces of the embedding vector space, and the simplest way to specify this mapping is through an affine transformation.

Figure 1

 1 Figure 1.5 (left) Vectors along a curve, here velocity vectors γ along the curve γ, live in different tangent spaces and therefore cannot be compared directly. A connection ∇ defines a notion of transport of vectors along curves. This allows transport of a vector γ(t -∆t) ∈ T γ(t-∆t) M to T γ(t) M, and the acceleration ∇ γ γ arises by taking derivatives in T γ(t) M. (right) For each point x ∈ M, the metric g defines a positive, bilinear map g x : T x M × T x M → R. Contrary to the Euclidean case, g depends on the base point, and vectors in the tangent space T y M can only be compared by g evaluated at y, i.e. the map g y : T y M × T y M → R.

Figure 1

 1 Figure1.6 (left) Geodesics starting at x with initial velocity v ∈ T x M are images of the exponential map γ(t) = Exp x (tv). They have zero acceleration ∇ γ γ and their velocity vectors are parallel transported γ(t) = P γ,t (γ(t 0 )). Geodesics locally realize the Riemannian distance so that dist(x, γ(t)) = t v for sufficiently small t. (right) The tangent space T x S 2 and Exp x give an exponential chart mapping vectors v ∈ T x S 2 to points in S 2 by Exp x (v). The cut locus of x is its antipodal point, and the injectivity radius is π. Note that the equator is the set Exp x (v) v = π 2 .

Figure 1

 1 Figure 1.7 (Left) The curvature tensor describes the difference in parallel transport of a vector Z around an infinitesimal parallelogram spanned by the vector fields X and Y (dashed vectors).(Right) The sectional curvature measures the product of principal curvatures in a 2D submanifold given as the geodesic spray of a subspace V of T x M. The principal curvatures arise from comparing these geodesics to circles as for the Euclidean notion of curvature of a curve.

Table 1 .

 1 1 Re-interpretation of standard operations in a Riemannian manifold.

http://bitbucket.org/stefansommer/theanogeometry/

http://geomstats.ai

projection of the vector y onto x. Then, Log x (y) = θ (yπ x (x)) y -π x (y) , θ = arccos( x, y ).

(1.17)

2D Kendall Shape Space

Kendall shape space [START_REF] Kendall | Shape Manifolds, Procrustean Metrics, and Complex Projective Spaces[END_REF] represents a shape as an equivalence class of all translations, rotations, and scalings of a set of k points, landmarks, in the plane. A configuration of k points in the 2D plane is considered a complex k-vector, z ∈ C k . Removing translation by requiring the centroid to be zero projects this point to the linear complex subspace V = {z ∈ C k : z i = 0}, which is isomorphic to the space C k-1 . Next, points in this subspace are deemed equivalent if they are a rotation and scaling of each other, which can be represented as multiplication by a complex number, ρe iθ , where ρ is the scaling factor and θ is the rotation angle. The set of such equivalence classes forms the complex projective space, CP k-2 . We think of a centered shape p ∈ V as representing the complex line L p = {z p : z ∈ C\{0} }, i.e., L p consists of all point configurations with the same shape as p. A tangent vector at L p ∈ V is a complex vector, v ∈ V, such that p, v = 0. The exponential map is given by rotating within V the complex line L p by the initial velocity v:

Likewise, the log map between two shapes p, q ∈ V is given by finding the initial velocity of the rotation between the two complex lines L p and L q . Let π p (q) = p p, q / p 2 denote the projection of the vector q onto p. Then the log map is given by Log p (q) = θ q -π p (q) q -π p (q) , θ = arccos | p, q | p q .

(1.19) We will see an example of a different landmark space equipped with a geometric structure coming from the action of the diffeomorphism group in Chapter 4.

Rotations

The set of orthogonal transformations O(k) on R k discussed in Section 1.2.1 is the subset of linear maps of R k , square matrices U ∈ M (k,k) , that preserve the dot product: U x, Uy = x, y . In particular, they conserve the norm of a vector: U x 2 = x 2 . This means that x (U U -Id)x = 0 for all vectors x ∈ R k , which is possible if and only if the matrix U satisfies the quadratic constraint U U = Id. Thus, the inverse transformation is U -1 = U . The composition of two such maps obviously also preserves the dot product so this forms the group of orthogonal transformations with the


MEDIEVAL LEAD SMELTING ON THE Mont Lozère

France A Southern

Ploquin

Ph Allée

M C Bailly-Maître

S Baron

J L De Beaulieu

J Carignan

S Laurent

M Lavoie

C Mahé -Le Carlier

J Peytavin

M Pulido

Keywords: lead metallurgy, medieval, palynology, slag, isotope, petrology, pollution

During a survey on the western part of the Mont Lozère (Cévennes National Park), about 50 lead smelting sites were mapped. Pluridisciplinary research activities were engaged : archaeology with small excavations to check extension of sites, presence of archaeological levels, … ; survey around the granitic massif on mineralisations and mining traces ; characterization, typology of metallurgical wastes (slags, lining, …) ; lead isotopes on mineralisation, slag, soil, plant, … ; palynology and geochemistry on core collected in peatland. Five radiocarbon dates on charcoal linked to slags gave a medieval age (between 985 and 1277 AD). All investigated sites show the same slag typology. The main facies are a black one, often vitreous, and a grey one, white coated, with bigger Pb prills. Chemical data have to be interpreted as complex genesis for slag, i.e. for smelting process : participation of arena, addition of iron (oxydes ?), lime, ... but the signification of the different slag facies is not yet elucidated. Lead isotope data are in agreement with exploitation of mineralisations located in the S-SW border of Mont Lozère Granite but not with exploitation of mineralisations located in northern one. A new unexpected conclusion resulted from results on a peat core : another period of lead smelting, with the same chemical and Pb isotopic characters occured during late La Tène.

INTRODUCTION

This progress report concern a multidisciplinary work about the lead metallurgy whose traces are located on the western part of the granitic Mont Lozère (Cévennes National Park).

The Mont Lozère Granite is a late-hercynian one at the S-E of the Massif Central. As his southern brother, Mont Aigoual Granite, it intruded the metamorphic basement (Cévennes Schistes). After a long period of weathering and erosion, a new cover was deposed (liasic and jurassic formations, mainly lime and dolomite, which form the Causses). Mineralisations occur all around these granites, both in the basement and in the cover but only in the first decametres of the cover (Fig. 1). Most of mineralisation (vein or stratiform) bear Pb as main metal (galena) but also may contain As, Cu, Zn, Sb, and Ag as minor constituant. Associated minerals are quartz and/or barytine [START_REF]Geological maps[END_REF].

During preliminary survey on the western part of the Mont Lozère about 50 lead smelting sites (Fig. 2) were mapped inside a seven square kilometres area (J. Peytavin, S. Lhuillier). Few sites were added during subsequent surveys, but a new area was defined three kilometres eastward. Altitudes of all slag sites are between 1360 and 1430 m. Most sites show a very poor vegetation cover compared to neihbouring areas with no slags. Five 14 C dates conducted on charcoal fragments linked to slags gave a medieval age : between 985 and 1277 AD, overlaping the five dates with 95% confidence level, but four of them are focused on the first half XIIth. C. Because all sites seem show the same characters and typology, we are enclined to believe that these ages are relevant for the whole set of slag sites. But don't forget that a Galloroman mine (dated by oil-lamp, Le Bleymard-Neyrac area) and Renaissance mining (radiocarbon age near Cocurès) exist (Fig. 4). Complementary ways of research were engaged : archaeology with small excavations to check extension of sites, presence of archaeological levels, diagnostic for excavation, …; survey around the granitic massif on mineralisations and mining traces ; characterization, typology of metallurgical wastes (slags, lining, …) ; lead isotopes on Florac Mende Ramponenche mineralisation, slag, soil, plant, … ; palynology and geochemistry on core collected in peatland (Narses Mortes) ; heavy metals pollution in and around smelting workshops.

Historical context

Medieval Gevaudan had about the same extension as the present Département de la Lozère (48). After migrations times during early Middle Ages, Charlemagne (Carolus Magnus) created Aquitaine kingdom for his son Louis (Ludovicus) the VIIth ; this kingdom included Gevaudan. At about 1100, Gevaudan was binded to the County of Toulouse. During Medieval period, from the VIIth C. until Philippe the IVth (socalled Le Bel, end XIIIth C. / beginning XIVth), silver was the main metal for coinage ; consequently a renewal of interest for mining developed notably in the Cévennes as all along the southern border of the Massif Central [START_REF] Bailly-Maitre | L'argent, du minerai au pouvoir dans la France médiévale[END_REF]. The Coinage Hall in Mende might be created about 1170 (attribution of regalia to the Bishop of Mende by Louis the VIIth.). At about the same time Aldebert, the bishop of Mende, tried to lay hands on mines and mining right in his diocese, against local powers (baronnies, priory, …). Last coinage from Mende dated about from 1316. Relationship of mining in Gevaudan and vicissitudes during the Albigeois Crusade are not yet known but we have to emphasize coincidence with the last part of our present datations span (Meaux-Paris treaty : 1229). During the XIIIth C., the Capetian power has tried progressively to recover mining right. In 1306, an acte de paréage between the Bishop and the king Philippe the IVth fixed a joint-authority on Gevaudan.

Very few archives concern directly the Mont Lozère area and its mining/smelting activity, and they are later than the expected time for it. More precision on the datation of this medieval mining/smelting activity is required but archaeological investigations (excavation, …) have not yet given enough information. Also, very few information were obtained about workshop operation.

FROM ORE TO WORKSHOP

One of the aims of this study is to track the source of Pb ores which supplied the smelting sites using the Pb isotopic composition. The lead isotopic composition of galena collected in most surrounding ores is very similar to that of different slag samples (lead prills) sampled on workshop sites (Fig. 3). On the other hand, the high precision of the results allowed us to select the mineralised areas which were probably the ore sources. The Pb isotopic composition of slags is even more homogeneous ( 208/206 Pb = 2.0920.002 ; 206/207 Pb = 1.1790.001 ; 208/204 Pb = 38.6630.025 ; 207/204 Pb = 15.6650.006 ; 206/204 Pb = 18.4760.023), and will allow source tracing in the environment. If we assume that no fractionnement occur during the metallurgical chaîne opératoire, these values, may be used to characterize the lead which was produced on the Mont Lozère.

The slag samples show an homogeneous group unlike the whole set of ores samples. The ores that have supply the metallurgical activity are issue of North-Eastern Montmirat and the Western and Northern Bondons areas (Fig. 4). The others ones might be to exclude because of the different lead isotopic composition. If further investigations will confirm that, in spite of topological proximity, northern and north-eastern mineralisations did not supply ore to eastern medieval smelting workshops on the Mont Lozère, the place of the Barony of Tournel in the mining game will be difficult to situate.

INFORMATION FROM SLAG

On smelting sites of Mont Lozère, slags are studied by macroscopical and microscopical observations and by chemical analysis :

• pieces of lining and tuyère. They are composed by a clayey matrix with a granitic arena (quartz, potassic felspar, ...)

• two main types of slags, black vitreous slags (with inclusions of quartz) and grey slags wich present a white patina (Fig. 5 and6).

• locally, other types of slag may be found. e.g. on site 12, there are also many black and crystallised slags, some of them are flat and only one centimetre thick. Iron bloomery slags have been found on sites 8 and 61 but no evident traces of bloomery workshop were not yet discovered on this part of Mont Lozère. Because of the very high fragmentation of the samples, it is difficult to give details concerning the shape of the furnaces. However some specimens are enougth big and can be interpreted as moulding a basin (diameter about 30 cm). Following Biringuccio such basin could be placed at before or on the side of the furnace. Tuyère might indicate that there was using forced ventilation, but a structure which look like a furnace have two stony walls as a funnel in front of one of dominant wind, and sites are in cleared and windy area. It is also the case of the sites located in the area of Swaledale in Yorkshire [START_REF] Murphy | Smelting residues from boles and simple smeltmills[END_REF].

For petrological and chemical study of the various types of slag and lining, all the major and traces elements were taken into account in order to highlight various contributions within slags by means of diagrams of correlation. E.g. :

• Si, Al and K are correlated, which may imply that they have the same origin. However lining are made up almost exclusively by these three elements. One can thus consider that the contribution of the molten wall can be represented by the parameter "P = SiO 2 +A1 2 O 3 +K 2 O".

• CaO, P 2 O 5 and MgO present a very good correlation. This suggests a contribution of a material containing these elements. It could be dolomite (Ca, Mg)CO 3 , because this mineral constitutes part of the gangue of mineralisation. BaO does not show a clear correlation with other elements, although there is a small tendency to a positive correlation with CaO+MgO+P 2 O 5 . We can think that BaO comes from another material, and it could act of the barytine BaSO 4 which is also a mineral composing the gangue of the ore. One represents the contribution of the gangue of the ore by the parameter "G = CaO+MgO+ P 2 O 5 +BaO".

•The contribution of the ore itself is given by the content of PbO which comes directly from PbS Pb prill channel crystal. The content of PbO in slags is very variable and disturbs strongly the reading of the diagrams. Also, in order to be able to highlight the various types of additions, we chose to use the chemical analyses recomputed at 100% without lead (Fig. 7) : Fig. 7 : Chemical typology of slags (analysis are recalculated without Pb). In this way, almost all samples are distributed on a line (virtual correlation since the totality of the elements taken into account add up 100% of the analysis). Group 1 represents the composition of burnt or molten lining (without contamination by slag). Group 4 represents the load of the furnace i.e. the elements brought by the ore and the associated gangue. Groups 2 and 3 are in intermediate position between these the first two groups. They could be amalgams between the load of the furnace and the melted lining with a lining component more significant (group 2) or a composition charges dominant (group 3).

The presence of remainders of gangue in the charge seems important. According to Birringuccio, it would be preferable to use a cleanest possible ore. According to the experimental work of F Tereygeol (in progress), one could realize that obtaining a pure ore could be done easily during crushing on skins of animals. However the variable report/ratio of (CaO+MgO)/BaO seems to indicate one "bad" washing of the ore that a voluntary addition of products. In some ore galena may be too scatered to be washed perfectly or a part of gangue was left voluntarily ? Another group 5 is not on this line.They are black and crystallized slags of site 12. Some of these slags are very iron rich. The high content in iron is all the more significant as one moves away from the right-hand side. The diagrams carried out did not show any correlation with the other chemical elements. Iron thus comes here from another contribution that of the walls or from that of, the gangue of the ore. The alignment of slags is done towards another contribution : may be slags characteristic of the reduction of the iron ores in bloomery. As said beforef, bloomery slags have been found on sites 8 and 61. One can thus put forth the assumption that bloomery slags were added to ore. In its treaty, Birringuccio mentioned this feature (XVIth C.).

INFORMATION FROM PEATLAND (the word "peatland" is used instead of "peat-bog", seeing that they are minerotrophic) Palynology :

Several palynological studies have been conducted in the French Massif central over the last decades. The postglacial vegetation history is now well known, and the major human disturbances on the vegetation dynamics are documented at a regional scale [START_REF] Beaulieu | Évolution de la végétation sur la bordure montagneuse cévenole au postglaciaire, d'après les pollens[END_REF][START_REF] Beaulieu | Histoire de la flore et de la végétation du Massif Central (France) depuis la fin de la dernière glaciation[END_REF][START_REF] Reille | Late and postglacial vegetation, climate, and human action in the French Massif Central[END_REF]. However, less attention was given to the role of specific human activities on plant communities during the last thousand years. Several lead smelting sites found on the western part of the Mont Lozère are located near peatlands, allowing the opportunity to examine at a fine temporal scale the role of past metallurgical activities on the forest dynamics using paleobotanical indicators (pollen, plant-macrofossils). These activities used large quantity of wood for combustion, especially beech, a tree species currently absent from the altitudinal range where lead smelting sites were found. When did metallurgical activities begin on Mont Lozère ? What was the impact on the vegetation ? To answer these questions, we have analysed the pollen content of sediments collected from numerous peatlands. For some of these sites, geochemical analyses were also conducted.

Narses Mortes peatland is located at an altitude of 1400 m on the south-western edge of Mont Lozère (44°26'N; 3°36'E). Two sediment cores were collected : one at the centre of the site where the thickness of organic matter was maximum (140 cm), and another one near the edge of the peatland (lateral core; 138 cm). Pollen analyses were conducted on the lateral core at 2 cm intervals for the upper part of the core (0-74 cm), and at 4 cm intervals for the basal part (74-138 cm). Geochemical analyses were measured on the same levels (infra). Seven samples were radiocarbon dated by the AMS method. Dates are presented in calendar years.

The lateral core provides a record of more than 2000 years. Four major periods in the vegetation history were defined based on pollen records. During the first period (Zone 1; >275 BC), a forested environment is deduced from the abundance of birch (Betula), deciduous oak (Quercus) and beech (Fagus). The pollen representation (pollen percentages and concentrations) of all tree species declined rapidly after circa 275 BC (Zone 2). The opening of the vegetation cover is associated to a major increase in Pb anomalies reaching a maximum at depth 112 cm. Subsequently, lead anomalies began to drop progressively at the same time of an increase in the pollen representation of Fagus, Quercus and fir (Abies) indicating an expansion of the forest cover. Zone 3 (after circa 710 AD) is characterized by a second opening of the forest cover characterized by a significant decrease of all tree species and a maximum in the abundance of Poaceae and anthropogenous herbaceous species. Finally, the upper part of the diagram (Zone 4) show high values of Pinus corresponding to its introduction during the last decades on abandonned cultivated fields. High lead anomalies reflect modern atmospheric metal inputs.

Pollen analyses show that past metallurgical activities played a strong and significant role on the vegetation dynamics of Mont Lozère. Two major periods of deforestation occurred during the last 2000 years. The first one was characterized by a decrease of all tree species associated to a peak in lead anomalies. This deforestation is closely associated to metallurgical activities, and occurred during the Antic period, as indicated by the radiocarbon dates and the presence of pollen of planted tree species like chestnut (Castanea) and walnut (Juglans). The second episode happened more recently during the Medieval period leading to the definitive disappearance of beech and oak in the area where lead smelting sites are located. However, this deforestation is firstly attributed to agropastoral activities because of the high representation of Poaceae and ruderal taxa, but it is likely probable that metallurgical activities also played a role. The level line, in between zones 3 and 4 , is also the summit of the chemical profils (infra). After these profils, the antic metallurgy began at the 119 cm level until -95 cm (or 91) level at least. The medieval one might begin at 41 cm but it is overshadowed by industrial pollution from level 35 cm.

Chemical and Pb isotopic investigation :

The "lateral core" was splited into two longitudinal sections, one for palynological investigation (supra) whereas the second half core was retrieved and divided into 58 individual samples, at 2 cm intervals, from level 23 cm to 139 cm (too much recent plant in upper part). The "Narses Mortes" peatland is strongly minerotrophic and contains 8 to 60 % ash. Minerotrophic peatland records both atmospheric deposition, soils leaching and the grounwater influence. The measured metal concentrations are normalised to Al contents of peat samples and the metal/Al ratios are compared to that of the Mont-Lozère granite : relative excess in metal concentrations are found in peat samples. The Narses Mortes peatland presents an enrichment for most of the heavy metals for the upper samples, which might be attributed to atmospheric inputs during the last centuries. Only Pb and As are enriched in some deeper samples suggesting a contribution from the surrounding activities. Hence two phases of metallurgical activities can be distinguished : one beginning at 119 cm depth and an another one at 41 cm depth (Fig. 9-11). The 206 Pb/ 207 Pb vs depth diagram shows two anthropogenic periods illustrated by less radiogenic isotopic composition than the background composition. The oldest period seems to have long activity and may correspond to a single event because of the same isotopic composition during 30 cm. The second smelting activity began around 41 cm depth but the isotopic composition is not the same along the profil. Other inputs are also recorded (Fig. 12 and13). These two smelting activities, at different periods, have used the same ores. According to radiocarbon dates and to palynological data, the medieval smelting activity record began at about 41 cm depth. However, a very important activity, recorded in deeper samples (Latenian level, may be continued in Galloroman time) remains being confirmed by archaeological arguments. Minerotrophic peatlands can be used for proxies of atmospheric inputs. The atmospheric metal inputs are not disturbed by the leaching of surrounding soils or groundwaters because of the high isotopic gradient measured throughout the entire profil.

CONCLUSIONS / PROSPECT

Pollen analyses show that past metallurgical activities played a strong and significant role on the vegetation dynamics of Mont Lozère. Two periods with lead metallurgy are recorded by peatland ; however only the medieval one have yet been determined by the radiocarbon dates ; at this time the Antic one is presently "invisible" on the Mont Lozère. Both were using same ores (or ores with the same Pb-isotopic signal). Where are antic workshops? Are their characters similar to medieval ones ? Now, we have to multiply the radiocarbon dates on the whole area.

Two main types of slag are always associated (dark ones and white coated ones) : did they produced after two different operations , e.g. as in Yorkshire [START_REF] Murphy | Smelting residues from boles and simple smeltmills[END_REF], or during only one operation ? Participation of arena, addition of iron (oxydes ?), lime, ...are proved in some cases, but which ones may be generalised ? What role do such additions ? At this time, no well preserved workshops have been excavated.

Fig. 1 :

 1 Fig. 1 : Geological sketch. Black dots are main Pb (Ag) mineralisations. Black rectangle is the slag area.

Fig. 2 :

 2 Fig. 2 : slag sites in the Barrandon Pond area.

Fig. 3 :

 3 Fig. 3: Diagram 208 Pb/ 204 Pb versus 206 Pb/ 204 Pb. We used the NIST 981 Pb reference material and normalised on NIST 997 Tl (7). The RSD (ppm) is 149 ppm on the 208 Pb/ 206 Pb ratio (0.015%).

Fig. 4 :

 4 Fig. 4 : Relationship between mineralisation and smelting areas after Pb-isotopic signals.

Fig. 5 (

 5 Fig. 5 (left) : black glassy slag, with quartz and feldspar inclusions. Fig. 6 (right) : white coated slag, with Pb prill in the bottom of a channel.

Fig. 8 :

 8 Fig.8: Simplified pollen diagram of Narses Mortes peatland (lateral core). Results are expressed as pollen percentages (curves) and pollen concentrations (grains.cm -3 ; bars). Dots represent less than 1%. Pollen zones correspond to the regional vegetation history. (1) : 14 C determination for this sample (level 31 cm) canot be interpreted as an age because it contains high recent plants content. The level line, in between zones 3 and 4 , is also the summit of the chemical profils (infra). After these profils, the antic metallurgy began at the 119 cm level until -95 cm (or 91) level at least. The medieval one might begin at 41 cm but it is overshadowed by industrial pollution from level 35 cm.

Fig. 9 ,

 9 Fig. 9, 10, 11: Elemental chemical analysis. Fig.9 shows the Pb and As concentrations in the peat core and then their correlation. Fig. 10 shows the Zn, Sb and Cd concentrations. Fig. 11 shows the Enrichment Factor (EF) for Pb. EF correspond to : ( [Pb] / [Al] )sample/( [Pb] / [Al] )granite. The normalisation allows to estimate the excess Pb, then the anthropogenic resulting from metallurgical activities.

Fig. 12 :

 12 Fig. 12 : Diagram 206 Pb/ 207 Pb versus depth. We can see the phases of anthropogenic activities illustrated by the less radiogenic signature than the background value. At 79 and 81 cm depth, a modern contamination due to circulation of waters occurs.

Fig. 13 :

 13 Fig. 13 : Diagram 207 Pb/ 204 Pb versus 206 Pb/ 204 Pb. The same ore have been used but at different time. Even a wide background, it is possible to see the history of this core link with the anthropogenic activity on the Mont-Lozère massif.

We have also to confirm that northern and north-eastern ore deposits did not supplied the Mont Lozère workshops … and to understand this feature. We have also to determine data for estimating Pb (Ag) production. And, last but not least, such ancient pollution in and from smelting workshop may be used as a model for the evolution of present industrial pollution.