Finite Element Systems for vector bundles : elasticity and curvature

We develop a theory of Finite Element Systems, for the purpose of discretizing sections of vector bundles, in particular those arizing in the theory of elasticity. In the presence of curvature we prove a discrete Bianchi identity. In the flat case we prove a de Rham theorem on cohomology groups. We check that some known mixed finite elements for the stress-displacement formulation of elasticity fit our framework. We also define, in dimension two, the first conforming finite element spaces of metrics with good linearized curvature, corresponding to strain tensors with Saint-Venant compatibility conditions. Cochains with coefficients in rigid motions are given a key role in relating continuous and discrete elasticity complexes.


Introduction
In this paper we first generalize the previously introduced framework of finite element systems (FES) [24] [33] so that it can treat, in particular, elasticity problems, and then provide concrete examples of finite element spaces, some old and some new, that fit the framework.
The general framework provides an approach to finite element discretizations of sections of vectorbundles, and complexes thereof, in particular differential forms with values in a given vector bundle. We make some comments about curvature, but most of the paper concerns the case of flat bundles. For applications in elasticity, the fiber can be identified as the space rigid motions.
In space dimension 2, one can distinguish between two differential complexes related to elasticity, which are formal adjoints of each other and give priority to stresses and strains, respectively. For the stress complex (61), we can check that the spaces defined in [53] and [3] fit the framework. For the strain complex (63), we introduce, also within the framework, some new finite element spaces. They model symmetric 2-tensors (metrics) with a good Saint-Venant operator (linearized curvature).
In the obtained finite element complexes, rigid motion like degrees of freedom play a key role, at every index. The FES framework stresses this design principle, and relates it to the interpretation of elasticity in terms of rigid motion valued fields.
Defining discrete spaces of metrics with good curvature in dimension 2 should be useful, in view of the importance of curved surfaces in several branches of mathematics, both pure and applied, whatever the distinction is. Such applications will be explored elsewhere. Another motivation for this work was to prepare the way for similar constructions in higher dimensions, especially 3 (with classical elasticity in mind) and 4 (with general relativity in mind).
Previous work on FES. Until now, the FES framework has been formulated in order to discretize de Rham complexes. It has been used to define finite element complexes of differential forms on polyhedral meshes [24], accommodate upwinded finite element complexes containing exponentials [26] [31], give new presentations of known elements [36] and to define elements with minimal dimension [28] under various constraints (such as containing given polynomials).
The regularity of the differential forms, in the above mentioned works, was L 2 with exterior derivative in L 2 , and the defined finite elements were natural generalizations of, in particular, the Raviart-Thomas-Nédélec (RTN) spaces [67] [60]. The continuity is thus partial, and can be expressed as singlevaluedness of pullbacks to interfaces, corresponding, for vector fields, to continuity in either tangential or normal directions.
In [33] we extended the FES framework so as to be able to impose stronger interelement continuity. For instance, for a conforming discretization of the Stokes equation, one would like to have spaces of fully continuous vector fields, satisfying a commuting diagram with respect to the divergence operator. For de Rham sequences of higher regularity (H 1 and, if desired, exterior derivative in H 1 ), the required continuity can be expressed as singlevaluedness of all components of the differential form and, if desired, of its exterior derivative too, on interfaces. This led us, in [33], to define FE complexes starting with the Clough-Tocher element, which is of class C 1 , instead of, say, Lagrange elements, which are of class C 0 . This provided the first conforming polynomial composite Stokes element in dimension 3 (and higher), with piecewise constant divergence and the degrees of freedom of [15]. The latter seem to be the natural ones for lowest order approximations.
For mixed finite element methods (MFE), pairs of finite element spaces that are compatible in the sense of Brezzi [20] should be identified. A particularly convenient tool for this purpose, has been the so-called commuting diagram property, see for instance [69] page 552 and 570 and compare with [16] §8. 4 and §8.5. It can sometimes be derived from a commutation property of the interpolators associated with the degrees of freedom. In particular, in [61], finite element grad − curl − div complexes were presented with degrees of freedom providing commuting diagrams.
Arbitrary order finite element complexes of differential forms were defined in [50]. Whitney forms [78] [77] and the RTN spaces appear as special cases (lowest order -arbitrary dimension, and arbitrary order -low dimension, respectively). This connection between numerical methods and differential topology was first pointed out in [18]. Computational electromagnetics has been one of the main motivations [19][62] [51]. Its interpretation in terms of differential forms is quite clearcut compared with the case for, say, computational fluid dynamics.
Systematically developing the theory of finite elements in terms of differential complexes equipped with commuting projections was advocated in [1]. Relating de Rham complexes to differential complexes appearing in elasticity, and viewing both as special cases of complexes of Hilbert spaces, has lead to the finite element exterior calculus (FEEC) [6][8] [2].
Stability of numerical methods is, in many cases, equivalent to the existence of projections onto the finite element spaces, satisfying commuting diagrams, and having appropriate boundedness or compactness properties [6] [38]. Uniformly bounded commuting projections can often be obtained from the interpolator associated with degrees of freedom, by a smoothing procedure [70][23] [6][37] [35] (in chronological order of submission).
The FES framework downplays the role of degrees of freedom and stresses that, for a finite element space on a cell, there are implicit finite element spaces on the subcells. The claim is that making these spaces explicit has numerous benefits. For instance, it suggests defining FE spaces recursively, starting with low-dimensional cells, using some extension procedure, such as harmonic extension [24].
The latter technique is also basic to the Virtual Element Method (VEM) [14]. VEM grad − curl − div complexes are constructed in [13] ; for an interpretation as a FES of differential forms, valid in arbitrary dimension, see §2.1 in [28]. Notice that [13] is based on non-homogeneous harmonic extensions, which produces quite large spaces compared with the minimal ones [28]. Harmonic extensions with respect to a modified metric is used in [26] [31]. In [30] harmonic extensions with respect to a flat metric but a non-zero connection form is used. This illustrates that the FES framework can accomodate FE spaces constructed as solutions of certain PDEs, where it is not necessary to have explicit solution formulas in order to be able to compute with them.
FES principles: restrictions, differentials and cochains. In FES, spaces are posited on cells of all dimensions. Interelement continuity is expressed through certain restriction operators, from spaces on cells to spaces on subcells. The spaces on the subcells can also be arranged in complexes, for certain induced differential operators. The restriction operators and the induced differential operators must satisfy commutation relations to define a system. For a given FES, a condition of compatibility (Definition 3.5), expressed as exactness properties of the restrictions and the induced differentials, ensures the existence of good degrees of freedom and in particular that the harmonic interpolator is well defined and commutes with the differential operators.
In [24], concerning de Rham complexes of low regularity (L 2 with exterior derivative in L 2 ), the relevant restriction operators were pullback by inclusion maps, and the induced differential operators were, again, the exterior derivative. In [33], to ensure H 1 regularity, the restriction operators could remember all components of the differential forms on the subcell, and possibly of the exterior derivative as well (if it is required to be H 1 ). The induced differential operators now acted on all this information. Thus appeared some new vectorbundles on subcells, linked by differential operators that were not exactly the exterior derivative on the subcell: they retain additional information about the ambient cell. We therefore, for the framework, considered general complexes of spaces, not just complexes of differential forms.
General degrees of freedom are not essential in FES, but they are certainly accomodated and sometimes very convenient. On the other hand certain degrees of freedom are paramount for the development of the theory. For de Rham complexes these degrees of freedom are the integration of k-forms on the kdimensional cells of the mesh. This gives rise to the de Rham map, which maps from differential forms to (real valued) cellular cochains ; it commutes with the differentials.
For elasticity complexes, we contend that cellular cochains with coefficients in rigid motions are the right analogue. More precisely we introduce, for each cell of each dimension, a space which is naturally isomorphic to the space of rigid motions. Cochains with coefficients in these spaces form a complex. A generalized de Rham map from elasticity fields to such cochains with coefficients, is then defined and shown to commute.
Finite element elasticity complexes. In dimension 2, for elasticity problems, the stress complex is implicitly used in [53] and made explicit in [3]. Since then, many more discrete stress complexes have been defined, both conforming [9] and non-conforming [10] [7]. See [52] and the references therein for more examples. Notice that the stress complexes in [53] and [3] are composite and start with Clough-Tocher elements, whereas those in [9] are polynomial and start with Argyris elements.
For the systematic design of discrete elasticity complexes, a link between de Rham complexes and elasticity complexes, known as the BGG construction [43], has been developed [4] [5]. In [5], known finite element de Rham complexes were tensorized with vectors in order to get vector valued de Rham complexes. Under a surjectivity condition (see page 58), the diagram chase then yielded new spaces for the elasticity complexes.
The finite element complexes defined here behave naturally with respect to the BGG diagram chase. That is, we can define finite element spaces for some interlinked vector valued de Rham complexes, such that the diagram chase at the discrete level works exactly as at the continuous level: isomorphisms at the continuous level correspond to isomorphisms at the discrete level. One thus needs a large supply of discrete de Rham sequences, corresponding to different regularities, that match at different indices. While this can be dispensed of in the presentation of our elasticity elements, it was an important guiding principle towards their design and we have included remarks to this effect.
Finite element de Rham complexes of higher regularity have been constructed [44] [63]. See also [49] [48] for related Stokes elements. A motivation behind [33] [32] was to have enough such sequences to address elasticity through diagram chasing.
For elasticity complexes and related BGG diagram chases, the case of C ∞ regularity is well established in the literature, but the choice of Sobolev spaces is often not explicit. We introduce several Sobolev spaces for our complexes, many of which are not simple tensor products. There are several possible choices for each smooth complex. To obtain the stress complex one can do the chase in (67) or (72). In the latter, the regularity is expressed with a differential operator that acts on columns, whereas the differentials of the complex act by rows. For the strain complex we study two different regularities, corresponding to two different regularities in the diagram chase. Here also the regularity is expressed in terms of differential operators acting on columns as well as rows. A rationale behind our choice of Sobolev spaces is given in Remark 4.9.
Another tool we have developed for the purposes of constructing elasticity elements are Poincaré -Koszul operators for elasticity complexes [34]. The Cesaro -Volterra path integral is but one example.
Numerical methods for curvature problems. For the general framework, the main novelty here, compared with [33], is that we introduce some generalizations of the de Rham maps. We are interested in discretizing sections of vector bundles. These are equipped with a connection. For applications in (linear) elasticity this connection is flat. We have implicitly linearized around the Euclidean metric, for which the Levi -Civita connection is flat, as well as other associated connections. But, for the definition of discrete vector bundles, we have also been mindful of situations where non-zero curvature is centre stage, and inspired by numerical methods developed for such problems: Regge Calculus (RC) [68], a discrete approach to general relavivity, can be interpreted in a finite element context [22][25] [27] and extended to higher orders [58]. One then obtains, in dimension 2, strain complexes of low regularity: they end with discrete spaces containing measures, typically Dirac deltas at vertices. Here, on the contrary, the finite element fields are at least square integrable throughout the complexes. Notice also that, even though the Regge strain complexes are equipped with commuting projections, the degrees of freedom do not contain rigid motions. They are nevertheless natural: Regge metrics are determined by edge lengths (squared), at lowest order.
A discrete Gauss-Bonnet theorem is valid for RC, as can be proved combinatorially, or by a smoothing technique as introduced in [27]. For the elements presented here, Gaussian curvature (linearized or not) is well defined by classical formulas, so that the Gauss-Bonnet theorem is immediate. While the regularity of Regge elements seems adapted to general relativity theory, higher regularity, as achieved here, could be important to other PDEs in Riemannian geometry, such as those treated in [12][21] [40].
Lattice Gauge Theory (LGT) [79], as extended to a finite element context [29], was also at the back of our minds during this work. In LGT one defines discrete connections and curvature, as well as a discrete Yang-Mills functional, but it is less clear what the discrete covariant exterior derivative and Bianchi identity should be. By contrast, for the discrete theory we develop here, both of these are explicit: the former is to some extent the basic building block, and the latter is obtained in Theorem 2.3.
For the flat case, for which we introduce the FES framework, we prove a variant of the de Rham theorem (e.g. [66] §V.3.): the de Rham map induces isomorphisms on cohomology groups, from the space of gobal sections to the cochains with coefficients, see Theorem 3.2.
Exponential fitting. A role for covariant exterior derivatives in a FES context was also indicated in [26], defining upwinded complexes of differential forms, generalizing exponentially fitted methods. Such constructions have applications to PDEs describing for instance convection diffusion problems or band gaps in photonic cristals. See Remark 3.9 for further details.
Sheaves, differential geometry and formal theory of PDEs. The FES framework can be interpreted as a discrete sheaf theory (see Remark 3.1). For our purposes, sheaf theory can be summarized as a framework for gluing fields that are defined piecewise. This is done in Sobolev spaces with the help of, in particular, partitions of unity, and in finite element spaces, by imposing various forms of continuity through interfaces in a mesh. Sheaves help express this analogy. Discrete sheaf theory for various applications is also developed in [46] [41]. Sheaf theory is usually expressed in the language of categories, which also originates in algebraic topology. It has been suggested that this language, sometimes derided as "abstract nonsense" (following Steenrod, with more or less affection), could be useful not only for other branches of mathematics [54], but for sciences in general [75].
What we refer to here as elasticity complexes are sometimes called Calabi complexes in differential geometry, whereas rigid motions correspond to Killing fields. Sheaf theoretic approaches to such, with applications to linearized gravity, have been considered, see [55].
Our discrete Bianchi identity, being combinatorially quite elementary, has precedents, of course. In particular similar identities can be found in synthetic differential geometry [57]. This axiomatic approach to geometry grew out of topos theory [59] [54], introduced by Grothendieck to provide the ideal "double bed" for the espousals of continuous and discrete theories. It remains to see if it will appeal to numerical analysis.
There is a formal theory of partial differential equations, couched in homological algebra and developed by, in particular, Spencer [74]. See e.g. [65] [72]. Perhaps our paper can hint towards a corresponding theory for finite element spaces. For instance: could there be master finite elements for the Spencer sequences, from which all specific examples could be deduced by diagram chases?
Though not necessary in order to express the end product, namely new finite element spaces, the above mentioned philosophies (more specifically: homological approaches to differential operators and categorical approaches to geometry), were inspiring to us for developing them, even from a distance.
Outline. The paper is organized as follows. In §1 we develop a theory for discrete flat vector bundles. Results concerning vector bundles with curvature are relegated to §2. In §3 we detail the framework of finite element systems, for a given discrete flat vector bundle. In §4 we provide background on elasticity, including relevant differential operators, differential complexes, the BGG diagram chase and Poincaré operators. In §5 we detail FES for the stress complex, detailing the main example of the Johnson-Mercier element. In §6 we detail FES for the strain complex, with two different regularities, providing the new examples of finite elements for strain tensors (metrics) with compatible Saint Venant operator (linearized curvature).

.1 Cellular complexes and cochains
Let T be a cellular complex. If T, T are cells in T we write T T to signify that T is a subcell of T . Each cell of dimension at least one is supposed oriented. Given two cells T and T in T , their relative orientation is denoted o(T, T ). It is 0 unless T is a codimension one subcell of T , in which case it is ±1. The subset of T consisting of k-dimensional cells is denoted T k . The space of k-cochains is denoted C k (T ) and consists of the maps from T k to R. In other words a k-cochain assigns a real number to each cell of dimension k. Notice that C k (T ) has a canonical basis indexed by T k .
The cellular cochain complex is denoted C • (T ). Its differential, also called the coboundary map, is denoted δ : C k (T ) → C k+1 (T ). Its matrix in the canonical basis is given by relative orientations.
All complexes considered in this paper are cochain complexes, in the sense that the differential increases the index. • For each T ∈ T we suppose that we have a vectorspace L(T ). We call this a discrete vectorbundle. We call L(T ) the fiber of L at T .

Discrete vectorbundles with connection
• Moreover, when T is a codimension 1 face of T , we suppose that we have an isomorphism t T T : L(T ) → L(T ), called the transport map from T to T . We call this a discrete connection. LGT was initially defined for cubical complexes [79]. An analogue for simplicical complexes was developed in [29]. There, a discrete vector bundle corresponds to a choice of vector space attached to vertices only, whereas we here associate a vectorspace to each cell, of every dimension, in T . Moreover, in LGT, a discrete connection is defined only on edges, as a choice of isomorphim between the vectorspaces attached to its two vertices ; here the discrete connection has many more variables.
Definition 1.2 (Flatness of discrete connections). Whenever T is a codim-2 face of T , if we let T 0 and T 1 be the two codimension 1 faces of T which have T as a common codimension 1 face, we require that the following diagram commutes: Or, if one prefers: For reasons that will appear later we say that a discrete connection having this property is flat.
For instance we can choose a fixed vector space V and let t T T = id V . This defines a discrete vectorbundle with a flat discrete connection. Discrete vectorbundles of this form, for a choice of vectorspace V , will be called trivial discrete vectorbundles. We may speak of the trivial discrete vectorbundle modelled on V , to make the choice of V explicit.
In this setting one defines a cochain complex with coefficients in L, denoted C • (T , L), as follows: Definition 1.3 (Cochains with coefficients).
• The space C k (T , L) is nothing but Π T ∈T k L(T ), whose elements will be families (u(T )) T ∈T k such that for each T ∈ T k , u(T ) ∈ L(T ). Such a u will be called a k-cochain with coefficients in L.
The operator δ • t will be called the discrete covariant exterior derivative. We notice that we do indeed have a complex: Proof. This is a direct consequence of (1), given what we know about relative orientations.
This definition can be used in particular when a vector space V has been chosen and we let L(T ) = V for all T ∈ T and t T T = id V for all T, T . This cochain complex will be denoted C • (T , V ). With this notation we have in particular that C • (T , R) = C • (T ), the standard cellular cochain complex introduced previously. Remark 1.2 (on invertibility of transport operators). In Definition 1.1 we could allow t T T : L(T ) → L(T ) to be just a morphism (not necessarily an isomorphism), and still get at complex C • (T , L) from (3). However, we have in mind situations where the transport operators t T T mimick the parallel transport associated with a connection on a vectorbundle, and these are isomorphisms. Compare with [45] §I.4.7.

Transport along some paths within a cell
We notice that when the discrete connection is flat, transport along paths within a given cell just depends on the endpoints. We will use the following more precise statement: Lemma 1.2. Suppose T, T ∈ T , and that T is a codimension-k face of T , for some k ≥ 2. Then all sequences T T 0 · · · T k−1 T , where each term in this sequence is a codimension-1 subcell of the next, give the same map t T T k−1 · · · t T0T from L(T ) to L(T ). This map will be denoted t T T .
We could therefore make the following alternative description of a discrete vectorbundle with connection: Definition 1.4 (equivalent definition of flat discrete vectorbundles). For each T ∈ T we suppose that we have, as before, a vectorspace L(T ). Moreover, when T T , we suppose we have an isomorphim t T T : L(T ) → L(T ). We require that t T T = id L(T ) and also that whenever T T T we have t T T = t T T t T T .

Discrete gauge transformations
Suppose we have, for each T ∈ T , two choices of vectorspaces denoted L(T ) and L (T ). For T a codim-1 face of T we suppose that we have transport maps t T T : L(T ) → L(T ) as well as t T T : L (T ) → L (T ). Under these circumstances we define an isomorphism from (L, t) to (L , t ) to be a family of isomorphims θ T : L(T ) → L (T ), one for each T ∈ T , such that: Lemma 1.3. Under the above circumstances θ induces an isomorphism of com- Proof. Bijectivity is obvious. We prove that θ is a cochain morphism (i.e. commutes with the differentials). We let δ • t be the cochain- as required.
A family of isomorphisms θ as above may also be referred to as a discrete gauge transformation. This terminology is used in particular when we have one discrete vectorbundle L, but with two different choices of discrete connections t and t ; in this case θ T will be an automorphism of L(T ), for each T .
Proof. The requirement is that the following commutation relation holds for subcells S, S of T such that S S: This can also be written: This identity holds according to Lemma 1.2.
The preceding corrollary is local, in that it concerns a single cell T ∈ T . The cohomology of the global space C • (T , L) could be different from that of C • (T , L(T )) (for any choice of a fixed T ∈ T ). Remark 1.3 (classification of flat discrete vector bundles). Flat vector bundles over a manifold M , modulo gauge transformations, correspond to representations of the fundamental group of M , modulo conjugacy. For a precise statement, see for instance Theorem 13.2 in [76] (in the context of principal bundles). This seems to carry over to the discrete setting. The fundamental group of a cellular complex can be replaced by the so-called edgepath group (described in Chapter 3 in [73] for the case of simplicial complexes).

Discrete vector bundles with curvature 2.1 Cubes in the barycentric refinement
Consider now an n-dimensional cube. For definiteness we consider the unit cube [0, 1] n in R n and denote it as S. We let (e i ) i∈[ [1,n]] be the canonical basis of R n .
The vertices of S can be indexed by the subsets of [ [1, n]]. For any subset U of [[1, n]], we let p U be the vertex of S defined by: Any face S of S is uniquely determined by two vertices p U and p V , with U ⊆ V , such that the vertices of S are exactly those of the form p W for U ⊆ W ⊆ V . Then we also have: Subsets of [ [1, n]] are partially ordered by inclusion, and this uniquely determines a partial ordering of the vertices. Then p U is the smallest vertex of S and p V the largest. Let T be a simplex. For each face T of T we let b T be the isobarycenter of T or, more generally, a point in the interior of T (so a barycenter with respect to some strictly positive weights), referred to as the inpoint of T . Recall that the barycentric refinement of T is the simplicial complex whose k-dimensional simplices are those of the form [b T0 , b T1 , . . . , b T k ] such that the T i are two by two distinct subsimplices of T satisfying T 0 T 1 . . . T k . We call such subsimplices of T barycentric simplices.
The barycentric refinement may be coarsened as follows: for any two subsimplices T and T of T such that T T , we consider the cell S(T , T ) which is the union of all the barycentric simplices When we start with a simplex T , the cells of the form S(T , T ) form a cellular complex, where each cell is, combinatorially, a cube. The same holds true if the cell T we start with is a cube. We may consider that b T is the smallest vertex of S(T , T ) and that b T is the largest. The vector b T − b T points towards the center of T . When T is an n-dimensional simplex, this procedure will divided it into (n + 1) cubes of dimension n, each one of the form S(T , T ), where T is a vertex of T . We call this the cubical refinement of T .

Discrete curvature and Bianchi identity
If we relax condition (1) we model the parallel transport associated with connections with curvature, as opposed to flat connections. We still define the discrete covariant exterior derivative by (3). We no longer have δ k+1 t δ k t = 0, as in Lemma 1.1. When we compute δ k+1 t δ k t we get an operator C k (T , L) → C k+2 (T , L), which one would like to interpret in terms of a curvature.
Definition 2.1 (Discrete curvature). We suppose that T and T are to cells of T , such that T is a codimension 2 face of T , and we let T 0 and T 1 be the two codimension 1 faces of T which have T as a common codimension 1 face. The curvature of t is then c t defined by: which is associated with the square S(T , T ) (element of the cubical refinement of T ), whose set of vertices is associated with {T , T 0 , T 1 , T }. The sign in this definition is given by the orientation of the square, which is chosen such that the orientations of the two transverse cells T and S(T , T ) induce the orientation of T .
The definitions were chosen so as to have the trivial: We have a Bianchi identity in this setting, which we now detail. Recall that the usual Bianchi identity says that the covariant exterior derivative of the curvature 2-form, which is a priori a certain endomorphism valued 3-form, is 0. In this identity, the relevant covariant exterior derivative is the one associated with the induced connection on the bundle of endomorphisms. The discrete identity will assert that certain linear operators attached to the 3-dimensional cubes, in the cubical refinement, is 0.

Definition 2.2 (Discrete bundle of endomorphisms and its connection).
• For T ∈ T we consider the previously introduced cubical refinement of T , whose k-dimensional cells are (combinatorial) cubes of the form S(T , T ), for subcells T and T of T , where T has codimension k in T . We let S denote the cubical refinement of T . We then define End(L) to be the discrete vectorbundle on S, whose fiber at the cube S(T , T ) is the space of linear maps from L(T ) to L(T ).
• The discrete vectorbundle End(L) on S inherits a discrete connection from the discrete connection of L on T , as follows.
Consider a k-dimensional cube S = S(T , T ), where we say that b T is the smallest vertex and b T is the largest. When S is a codimension 1 face of S there are two possibilities : either b T is the smallest vertex of S and then we let b T0 be the largest, or b T is the largest vertex of S and then we let b T0 be the smallest. We define the transport operator on End(L) through: The spaces C • (S, End(L)) are defined as before and the discrete covariant exterior derivative linking these spaces, is defined as in (3) from the given induced discrete connection on End(L).
Remark 2.1. A discrete connection for the discrete bundle L over T is thus an element of C 1 (S, End(L)), where the element of End(L) attached to each edge of S is bijective. Compare with the fact that, in the continuous setting, the difference between two connections is an endomorphism valued 1-form. Theorem 2.3 (discrete Bianchi identity). The curvature of (L, t), which is defined as an element of C 2 (S, End(L)) by (14), has a covariant exterior derivative (element of C 3 (S, End(L))) which is zero.
Proof. That the discrete covariant exterior derivative of the curvature is zero expresses that for each 3-dimensional cube S(T , T ) a certain linear map from L(T ) to L(T ) is zero. This linear map is a sum of maps of the form: where the cells T T 0 T 1 T represent vertices of the cube, in increasing order. The sum consists of two such contributions from the curvature of each of the six faces of the cube. We thus get twelve maps of the form (17). They cancel two by two ; in fact each map t T T1 t T1T0 t T0T appears twice in the sum, with different signs.
Remark 2.2 (consistency). One would like the discrete covariant exterior derivative to be in some sense consistent with a continuous one. Recall that the coboundary operator acting on (realvalued) simplicial cochains is isomorphic to the exterior derivative acting on Whitney form, via the de Rham map. One would like a similar interpretation of the discrete covariant exterior derivative.
For cellular complexes, an analogue of Whitney forms was provided in [24], by solving recursively, the PDE system d du = 0 and d u = 0, or a discrete analogue. One motivation for identifying induced operators on subcells is to extend this construction to other differential complexes, where one wants to find preimages of cochains with coefficients. This connects with a broader theme of defining finite elements as solutions of local PDEs (possibly discretized at a subgrid scale).
This also raises the question of, to which extent, from a discrete vector bundle, one can reconstruct a continuous vector bundle. In the flat case this seems unproblematic. In the presence of curvature, a condition of small curvature might be necessary, to mimick that fibers vary continuously in the continuous setting. For instance, one could require that the maps in Lemma 1.2 should be close to each other, in some sense. As interesting and perhaps simpler special cases, one could consider the reconstruction of line bundles, and bundles over two-dimensional manifolds.
See [56] for similar considerations. Notice that the discrete curvature transforms naturally under discrete gauge tranformations. Indeed, consider a discrete connection t for L and a discrete connection t for L , as well as gauge transformations θ : L → L such that (5) holds. Then we have: Also, given a gauge transformation θ from (L, t) to (L , t ), there is an induced gauge transformation from End(L) to End(L ), equipped with their induced connections. For a cube S(T , T ) it maps u : The above transformation of curvature can be seen as a special case.

Remark 2.4 (Chern classes).
A motivation for developing a discrete Bianchi identity is to develop a corresponding notion of Chern classes. A first Chern class is defined in this spirit in [64]. The topic is discussed in Chapter 7 of [71].

Finite element systems
Definition 3.1. We fix a flat discrete vector bundle L on T in the above sense. A finite element system on T consists of the following data, which includes both spaces and operators: • We suppose that for each T ∈ T , and each k ∈ Z we are given a vector space A k (T ). For k < 0 we suppose A k (T ) = 0.
• For every T ∈ T and k ∈ Z, we have an operator d k T : A k (T ) → A k+1 (T ) called differential. Often we will denote it just as d. We require d k+1 • Given T, T in T with T T we suppose we have restriction maps: subject to the commutation relations (when T T T ): • For any k-dimensional cell T in T we suppose we have an evaluation map e : A k (T ) → L(T ). We suppose that the following formula holds, for see also Remark 3.6 below.
Remark 3.1 (presheaf interpretation). The first three points in the definition can be rephrased in terms of categories and sheaves. First, the cellular complex T may be considered as a category, where the objects are the cells and the morphisms are the inclusion maps. Complexes of vector spaces also constitute a category. The first three points then say that the family A • (T ), indexed by T ∈ T , equipped with restriction operators associated with inclusions maps, is a contravariant functor from T to the category of complexes of vector spaces. It is therefore a presheaf of complexes (according to the definition in [54] Chapter 17).
Remark 3.2 (inverse system interpretation). The cellular complex T can also be interpreted as a partially ordered set, the order being the inclusion relation. Then the family A • (T ), indexed by T ∈ T , equipped with restriction operators, constitutes an inverse system (also called a projective system) of complexes. See for instance [54] §2.1. (20) is a generalization of Stokes' theorem. Indeed Stokes' theorem may be regarded as the special case where the spaces A k (T ) consist of smooth enough realvalued k-forms on T , the discrete vectorbundle is given simply by L(T ) = R and t T T = id R , and the evaluation map

Remark 3.3 (Stokes' theorem). Identity
Definition 3.2 (gluing spaces). If T is a cellular subcomplex of T , we define: We notice that, in the special case where T is a cell and S(T ) denotes the cellular complex consisting of all the subcells of T in T , then the restriction maps provide an isomorphism: In what follows we will usually not distinguish between A • (T ) and A • (S(T )).
Remark 3.4 (continuity). The condition that, for T T we have u T = r T T u T , for elements u of A k (T ), can be interpreted as a continuity condition. Indeed when two cells have a common face, the condition enforces the two restrictions to the common face to be equal.
Remark 3.5 (glued spaces are inverse limits). If T is a cellular subcomplex of T , the spaces A k (T ) with T ∈ T constitute an inverse system, see Remark 3.2. Then A k (T ) is an inverse limit, and this defines it up to unique isomorphism. See for instance [54] Equation (2.1.2) p. 36.
The notation (21) will be used in particular in the following two cases: • T = T . The spaces of the form A k (T ) are those typically used in a Galerkin method to solve a PDE on the set covered by T . Such spaces will be referred to as the global spaces, as opposed to the local spaces A k (T ) for T ∈ T .
• T = ∂T for a given T ∈ T . By that we mean the set of cells in T included in the boundary of a given cell T ∈ T . That is, ∂T denotes the cellular complex consisting of the strict subcells of T .
Remark 3.6 (reformulation of Stokes' identity as a commutation relation).
Another way of formulating (20) is that for any cellular subcomplex T of T , the evaluations e T : A k (T ) → L(T ) (for T ∈ T k ) provide a morphism of complexes: We will later provide conditions under which the evaluation morphism (23) induces isomorphisms on cohomology groups. This would be an analogue of de Rham's theorem which asserts that the de Rham map, from real valued differential forms to real valued cellular cochains, gives isomorphisms between the respective cohomologies. We therefore refer to the map in (23) as the de Rham map.
We denote by A k 0 (T ) the kernel of the induced map r : A k (T ) → A k (∂T ). We consider that the boundary of a point is empty, so that if T is a point Definition 3.3 (Flabby finite element systems). We say that A admits extensions on T ∈ T , if the restriction map induces a surjection: We say that A admit extensions on T or is flabby, if it admits extensions on each T ∈ T .
This notion corresponds to that of flabby sheaves (faisceaux flasques in French [47]), due to the following result Proof. If one can extend from the boundary of a cell to the cell, then one can extend from subcomplexes to complexes, step by step, incrementing dimension by one each time.
Remark 3.7 (surjectivity of restrictions). In particular if A admits extensions, then, when T is a subcell of T , the restriction A • (T ) → A • (T ) is onto. However this is, in general, a strictly weaker condition than the extension property. To see this, consider for instance the finite element spaces A 0 (T ) consisting of P 1 functions on a quadrilateral S, on its edges E and on its vertices V . Then the restriction from A 0 (S) to each edge A 0 (E) is onto, as are the other restrictions from faces to subfaces, but the restriction from A 0 (S) to A 0 (∂S) is not onto, since the latter has dimension 4 but the former had dimension only 3. In practice therefore, finite element spaces on a square therefore include, in addition to the affine functions, a bilinear function.  • We say that A • is exact on a cell T ∈ T when the following, equivalent, conditions hold: -The following sequence is exact: and moreover the map j T : K(T ) → L(T ) defined in Remark 3.8 is an isomorphism. -The de Rham map A • (T ) → C • (T, L) induces isomorphisms on cohomology.
• We say that A • is locally exact on T when A • is exact on each T ∈ T .
Proof. The equivalence holds by Lemma 1.4.
Definition 3.5. We say that A is compatible when it is flabby and is locally exact.
Remark 3.9 (The flat connections of upwinding and band gap computations). Some PDEs, such as convection diffusion equations, can be expressed with covariant derivatives [26], see also [80]. Consider real valued differential forms. Choose A ∈ Ω 1 (U ), called the connection 1-form. Define the covariant exterior derivative d A by: Then we have: Then dA ∈ Ω 2 (U ) is identified as the curvature of A. We suppose dA = 0, so that the operators d A constitute a complex. For any contractible subdomain T of U we may choose φ ∈ Ω 0 (T ) such that dφ = A. Then we have, for u ∈ Ω k (T ): Furthermore φ is uniquely determined, if we impose, in addition, the value of φ(x T ) for some point x T ∈ T . For every cell T ∈ T , choose an interior point x T ∈ T . Let φ T : T → R be the unique function such that φ T (x T ) = 0 and dφ T = A. For u ∈ Ω k−1 (T ) with k = dim T , we notice: When T is k-dimensional, we let e T : Ω k (T ) → R be the map: And, when T is a codimension 1 face of T we put: Then the preceding identity becomes: as required in (20). As a slight variant, consider complex-valued differential forms (but A still real valued), and let: with comparable consequences. One can even restrict attention to A constant. On a torus, A is not in general globally of the form dφ, but of course locally this still holds and can be used in particular on individual cells. Such observations were made in [42] [17] in the context of band gap computations for photonic cristals. The corresponding numerical methods are variants of exponential fitting.

de Rham type theorems.
The following theorem extends Proposition 5.16 in [35]: Proof. Exactness gives that the map A • (T ) → C • (S(T ), L) induces isomorphims on cohomology groups. From there the proof proceeds as in [35].
We also have the following extension of Proposition 5.17 in [35]: Suppose that A has extensions. Then A is compatible if and only if the following condition holds: For each T ∈ T the sequence A • 0 (T ) has nontrivial cohomology only at index k = dim T , and there the induced map: is an isomorphism (it is well defined by (20)).
Proof. (i) For cellular complexes consisting only of vertices, the equivalence trivially holds because, when T is a point, A • 0 (T ) = A • (T ). (ii) We suppose now that m > 0 and that the equivalence has been proved for cellular complexes consisting of cells of dimension at most n < m. Consider a cellular complex consisting of cells of dimension at most m.
Let T ∈ T be a cell of dimension m. We suppose that the finite element system A is compatible on the boundary of T . Since the boundary is (m − 1)dimensional we may apply the de Rham theorem 3.2 there. In other words A • (∂T ) → C • (∂T, L) induces isomorphisms on cohomology.
We write the following diagram: Comments: • The complex C • 0 (S(T ), L) is rather trivial. The terms consist of cochains that are zero on the boundary of T . In other words the only non-zero space in the complex is at index k = dim T , where it is L(T ).
• On the rows, the second map is inclusion and the third arrow restriction.
Both rows are short exact sequences of complexes.
• The vertical maps are the de Rham map.
We write the two long exact sequences corresponding to the two rows, and connect them by the map induced by the de Rham map.
The equivalence is now proved in two steps: • Suppose that (25) is exact. Then the first and fourth vertical maps are isomorphisms. By the induction hypothesis the second and fifth are isomorphisms. By the five lemma, the third one is an isomorphism. This can be stated by the condition formulated in the theorem.
• Suppose that the stated condition holds. One applies again the five lemma to the long exact sequence, and obtains now that A • (T ) has the same cohomology as C • (S(T ), L).

Extensions, dimension counts and harmonic interpolation.
The following proposition almost exactly reproduces Proposition 2.2 in [36].
Proposition 3.4. Suppose that A is an element system and that T ∈ T . We are interested only in a fixed index k ∈ N. Suppose that, for each cell U ∈ ∂T , each element v of A k 0 (U ) can be extended to an element u of A k (T ) in such a way that, r U T u = v and for each cell U ∈ ∂T with the same dimension as U , but different from U , we have r U T u = 0. Then A k admits extensions on T .
Proof. In the situation described in the proposition we denote by ext U v = u a chosen extension of v (from U to T ).
Pick v ∈ A k (∂T ). Define u −1 = 0 ∈ A k (T ). Pick l ≥ −1 and suppose that we have a u l ∈ A k (T ) such that v and u l have the same restrictions on all l-dimensional cells in ∂T . Put w l = v − r ∂T T u l ∈ A k (∂T ). For each (l+1)-dimensional cell U in ∂T , remark that r U ∂T w l ∈ A k 0 (U ), so we may extend it to the element ext U r U ∂T w l ∈ A k (T ). Then put: Then v and u l+1 have the same restrictions on all (l + 1)-dimensional cells in ∂T .
We may repeat until l + 1 = dim T and then u l+1 is the required extension of v.

Proposition 3.5 (flabbyness vs dimension counts).
Let A be a FES on a cellular complex T . Then: • Equality holds in (40) if and only if A k admits extensions on each T ∈ T .
Proof. The proof in [36] works verbatim.
Definition 3.6. Given a FES A on a cellular complex T , a system of degrees of freedom is a choice of subspace F k (T ) ⊆ A k (T ) , for each k ∈ N and T ∈ T . In that situation we can define maps Φ k (T ) : where the brackets denote the canonical bilinear pairing F k (T ) × A k (T ) → R. We say that the system F is unisolvent on A if Φ k (T ) is an isomorphism for each T ∈ T .
We will use the following result: Proposition 3.6 (unisolvence of degrees of freedom and flabbyness). Suppose that A is a FES on a cellular complex T . Suppose that F is a system of degrees of freedom for A. Suppose that for each T ∈ T , the canonical map A k 0 (T ) → F k (T ) is injective. Suppose that T ∈ T is such that: Then F is unisolvent on A on the cellular complex S(T ), A is flabby on S(T ) and equality holds in (42).
Example 3.1. For each cell T ∈ T , equip A k (T ) with a continuous scalar product ·|· , typically a variant of the L 2 product. We define a system of degrees of freedom F k (T ) ⊆ A k (T ) as follows. For each k, we consider the following space of linear forms on A k (T ): where the last space in the direct sum should be included only for k = dim T . We call these the harmonic degrees of freedom. For compatible finite element systems these degrees of freedom are unisolvent.
If one considers that the linear forms in F k (T ) are defined on more general fields, these DoFs yield a commuting interpolator onto A • (T ), which we call the harmonic interpolator. For more on this topic see §2.4 of [36], in particular Proposition 2.8 of that paper. One can obtain such a FES inside any compatible FES by imposing the degrees of freedom (43) to be zero, while (44) are kept free. This generalizes the construction of Whitney forms on cellular complexes given in [24].

Discrete vector bundles: a dual picture
Notice that the degrees of freedom l • e T for l ∈ L(T ) appearing in (43) play a special role. In practice they often appear in a slightly different way, namely as integration against certain fields, forming a space M (T ) which is more tangible than L(T ) (the parallel transport operators acting on M (T ) can be more natural for instance).
We now make some remarks on this alternative point of view. We suppose that we have for each T a vector space M (T ) and a bilinear form ·, · T on A(T ) × M (T ). Moreover, when T has codimension 1 in T we suppose we have a bijective (linear) restriction map s T T : M (T ) → M (T ), subject to the condition that s T T 0 s T 0 T = s T T 1 s T 1 T . The generalized Stokes theorem takes the form, for φ ∈ M (T ) and u ∈ A k−1 (T ): In practice, this formula often arises as follows. The bilinear form ·, · T on A(T ) × M (T ) is the L 2 scalar product (with respect to, say, the standard Euclidean metric). The space M (T ) is the kernel of the formal adjoint of d. Identity (45) is obtained by integration by parts, m times when d is a differential operator of order m. Only boundary terms remain, by definition of the kernel of the formal adjoint. For boundary cells T , the space M (T ) could be obtained as the kernel of a formal adjoint, or could appear as natural restrictions to T of elements of M (T ).
For instance, one integrates the divergence of a vector field against the constants, which constitute the kernel of the gradient. Then the restriction operator on the constants is the usual trace and the restriction operator on the vector field is the trace of the normal component.
We will be interested in more complicated examples, for instance, integration of vector fields against the kernel of the deformation operator (i.e. rigid motions) which is the formal adjoint of the divergence operator acting on symmetric matrices), and integration of functions against the kernel of the Airy operator (i.e. affine functions), which is the formal adjoint of the Saint-Venant compatibility operator.
This data provides the discrete vector bundle defined by L(T ) = M (T ) and t T T = s T T . Moreover it gives evaluation maps e T : A k (T ) → L(T ) defined by u → u, · T .
Conversely, given the data L, t and e, one can define M (T ) = L(T ) , u, φ T = φ(e T (u)) for φ ∈ M (T ) and u ∈ A k (T ) and s T T = t T T , so the two points of view are equivalent.

Spaces
We work in dimension 2.
• We denote by V = R 2 the space of column vectors and by V t the space of row vectors, • We denote by M = R 2×2 the space of matrices.
• We denote by S = R 2×2 sym the space of symmetric 2 × 2 matrices and by K = R 2×2 skew the space of skewsymmetric 2 × 2-matrices. If T is a subset of R 2 (typically a simplex of dimension 0, 1 or 2, or an open subset of R 2 ), and E is a vector space (such as one of R, V, V t , M, S, K, or a product of such spaces) we denote by Γ(T, E) the set of maps from T to E.
These maps should be smooth enough that differential operators and traces (in the sense of restrictions to subsets of the boundary) make sense, so that one can define a suitable finite element system as subspaces. To get flabbyness, the regularity of the spaces should not be too big. Exactness under the differentials is also a regularity dependent issue. We will get back to the problem of getting the regularity of the fields right, to accomodate these constraints.

Differential operators
We will use the following differential operators, on scalar, vector and matrix fields, defined on some connected domain in R 2 .
• The gradient and curl of a scalar field are (row) vectorfields defined by: • The divergence and curl of (row) vectorfields are defined by: • The gradient of a (column) vectorfield is defined by: • The deformation of a (column) vectorfield is the symmetric part of its gradient: def • The divergence of a matrix field is the (column) vectorfield defined by: • The Airy operator acts as follows on a scalar field u: • The Saint Venant operator is the formal adjoint of the Airy operator: It encodes the Saint-Venant compatibility conditions for being a deformation tensor. The sven operator may also by interpreted as the linearization of scalar curvature around the standard metric.
• We define the matrix J as: It encodes a direct rotation around the origin by π/2. Then the Hessian and the Airy operator are linked by: • We can also introduce the operator on matrices: Then the Hessian and the Airy operator are related by: Remark 4.2 (Formulas in different bases). We will often work with an additional orthonormal oriented basis (τ, ν) of V. Well-known formulas will be used without further ado. For instance, for a scalar field u, curl u = ∂ ν uτ − ∂ τ uν. For a column vector field one would write curl u = ∂ τ u · ν − ∂ ν u · τ . Notice however that we have chosen to use differential operators only row wise. For a row vector field u the formula thus becomes curl u = ∂ τ uν − ∂ ν uτ , given that τ and ν are column vectors.

Remark 4.3 (Differential operators act on rows)
. We let all differential operators act row wise. To account for differential operators acting column wise, we combine with transposition, denoted t : u → u t .
Remark 4.4. We let the curl operator act on matrices row wise. Then we have the formula: where J and K were introduced in (56) and (58). This can be used to prove invariance properties of the curl of matrices under orthogonal basis change.

Elasticity complexes and diagram chasing
The preceding operators may be arranged into sequences as follows: Proposition 4.1 (Elasticity stress complex). We have a sequence: Here: The kernel of the Airy operator consists of the affine functions. The divergence operator is surjective. If U is contractible the sequence is exact. Proposition 4.2 (Elasticity strain complex). We have a sequence: Here: The kernel of the deformation operator consists of the rigid motions. The sven operator is surjective. If U is contractible the sequence is exact.
These complexes can be deduced from vector valued de Rham sequences through a diagram chase of the following type: Proposition 4.3. Suppose we have two complexes linked in a commuting diagram Suppose furthermore that we have an index i such that: • for j < i, φ j is injective, • for j > i, φ j is surjective.
Proposition 4.5. We have a commuting diagram: Here: Moreover t denotes transposition and: The elasticity strain complex (63) comes from (69) by Proposition 4.3.

Remark 4.5 (a variant with modified regularity).
On the bottom row of diagram (67) the differential operators act row-wise and the Sobolev spaces are defined row-wise. The middle space in the bottom row is not stable under taking the transpose. One could achieve that by enforcing conditions on the vertical divergence throughout the complex: Here: The diagram chase applied to (72) yields the same elasticity complex as (67), namely (61).
Remark 4.6 (regularity couples rows). On the top-row of diagram (69) the differential operators act row-wise, but the Sobolev spaces have a regularity whose definition couples the rows. The middle space is stable under taking the transpose.
We are interested in finding the natural boundary conditions for this space. For u in this space the τ τ component on the boundary makes sense. Indeed, there is u ∈ H 1 (U, S) such that sven u = sven u and v ∈ H 1 (U, V) such that def v = u − u (on noncontractible domains use the regular decomposition). Now u has well defined traces on the boundary, in particular τ τ traces, say in L 2 (∂U ). Moreover: and the trace of (grad v)τ is well defined in H −1/2 (∂U ). But then, taking scalar product with τ is problematic. The components of τ are not in H 1/2 (∂U ) on non-smooth domains, due to possible discontinuities. On a triangle one can still conclude that the trace of ((def v)τ ) · τ is well defined in the dual of H 1/2 00 of each face. But it is problematic that this not enough to define τ τ traces in, say, the dual of Lipschitz functions on the whole boundary.
Characterizing the space of traces seems simpler for smooth boundaries. The difficulties with triangles can then be expressed as being related to the presence of infinite extrinsic curvature at vertices.
One should also prove that an element of H 0 sven (U, S) has zero τ τ traces iff it is the limit of elements of C ∞ c (U, S).
Remark 4.8 (Strain complex with lower regularity). As an alternative to the elasticity strain complex (63) we can consider: Here: This complex can be deduced from the following diagram: Here: One may check that the skew operator on H 1 (U, R) maps onto the antisymmetric elements of this space.
Remark 4.9 (regularity and partitions of unity). The spaces defined in (64) and (79) are stable under multiplication by smooth functions. Thus, partition of unity techniques apply to them, enabling one to glue together smooth enough fields. The regularity in (64) and (79) seems quite adapted to FES techniques, which is about gluing together finite element fields. Notice that, on the other hand, the L 2 based graph norm space: is not stable under multiplication by smooth functions. It seems less amenable to FES techniques. We see that we have several choices of Sobolev norms in the strain complex. In the stress complex (61), on the other hand, the situation is simpler. The last two spaces are dictated by the Hellinger-Reissner principle, and this fixes the first space.
In dimension 3, this suggests the following two choices of norms, in the elasticity complex. For higher regularity: with: For lower regularity: with: On the other hand, the complex: constructed on the L 2 based graph norm space: is avoided, as it is not stable under multiplication by smooth functions. For Regge calculus, the functional framework of [25] was based on: which corresponds to lowering the regularity by 1 throughout the complex (84). Exactness properties of all these sequences can be deduced from results in [11].

Poincaré operators
Poincaré and Koszul operators are central to the construction of the finite element de Rham complexes of [60][50] [6]. They were also essential to the construction of the smooth finite element de Rham complexes in [33]. In [34] we provided similar operators for elasticity complexes. We detail them here, for the two different complexes we have in dimension 2. They will also be used for the construction of finite element elasticity complexes.
Stress elasticity complex. We suppose that the domains U is starshaped with respect to the origin 0.
We define Poincaré operators for elasticity as follows.
Proposition 4.6. We have the following identities.
p 2 increases polynomial degree by 1 (at most).
With these we may define two koszul operators. We may remark that on constant fields these operators reduce to: and: Strain elasticity complex. For the strain elasticity complex we also have Poincaré operators (we use the same notation). We suppose that the domain U is starshaped with respect to the origin 0. In practice we will use these operators with respect to other base points, namely the central vertex of the Clough-Tocher split.
p 2 increases polynomial degree by 2 (at most).

Definition 4.3 (Koszul operators).
Taking the Poincaré operators on homogeneous polynomials of degree r, we obtain the Koszul operators. More generally, we let the operator K r 1 : Γ(U, S) → Γ(U, V) be given by: whereas the operator K r 2 : Γ(U, R) → Γ(U, S) is given by: The next lemma shows some regularity of these operators.
Lemma 4.8. For w piecewise smooth on R(T ), we have p 2 w ∈ H 0 curl (T, S).
Proof. By straightforward calculations, (x r x rt )x = 0 and t curl(x r x rt )x = 0. This implies that wx rt x r and w curl(x rt x r ) have continuous tangential components even if w is discontinuous across of the interior edges of the Clough-Tocher split.

Induced operators and discrete vectorbundle
We consider a vertex V , in an edge E, in a triangle T . The oriented unit tangent on E is denoted τ and the normal is denoted ν, so that (τ, ν) is an oriented orthonormal basis of V. Spaces and operators for the stress complex may be arranged in the following commuting diagram: Remark 5.1. With the preceding notations we have that: • K(T ) consists of affine functions.
• K(E) consists of pairs (u, v) of realvalued functions on E such that v is constant and u is affine.
To complete the picture we need to define a discrete vector bundle. We take the dual point of view developed in §3.3 page 20.
• We let M (T ) denote the space of rigid motions. The rigid motions appear as the kernel of the formal adjoint of div : Γ(T, S) → Γ(T, V), acting on symmetric matrices, namely the deformation operator def.
• We let M (E) be the kernel of the (formal) adjoint of d 0 Therefore M (E) consist of pairs (u, v) of functions on E, where u is constant and v is affine.
We define bijective restriction operators M (T ) → M (E) and M (E) → M (V ) as follows.
• We define M (T ) → M (E) as the map sending a rigid motion φ to the pair (u · τ, u · ν), where restriciton to E is implied. Indeed the tangent component of φ on E is constant and the normal component on E is affine.
which is independent of E.
We need bilinear pairings.

Discrete spaces
We now define finite element spaces for the stress complex, following [53]. Our main novely is that by highlighting degrees of freedom involving the M (T ) spaces, we are led to a natural reduction of this space (where stresses have dimension 9 on an element, reduced from 15), see Remark 5.4. Notice that in [53], the displacement is chosen in the finite element space P 1 (T, V) rather than P 0 (R(T ), V). Thus their finite element pair is not part of a complex. A point of view emphasizing discrete complexes is developed in [3]. For a triangle T we denote by R(T ) the Clough-Tocher split of T .
Definition 5.1 (FE for the stress complex).
• We define: A 0 (T ) is thus the Clough-Tocher space. The degrees of freedom are: values at vertices (3), values of the gradient at vertices (3 × 2), integral of normal derivative on edges (3).
• We define, following [53]: by which we mean symmetric matrix fields that are piecewise polynomial of degree at most 1, with continuous normal components on interior edges. The degrees of freedom are: on edges E, E uν · v for v ∈ P 1 (E, V) (3 × 4), and on T , integration against P 0 (T, S) (3).
• We define: that is the space of piecewise constant vector fields. The degrees of freedom are integration against P 1 (T, V) (6).
These finite element spaces are represented in Figure 1. • for A 0 (T ), pairings of restrictions with M (V ) at vertices V can be recovered from the DoFs.
• for A 1 (T ), pairings of restrictions with M (E) at edges E can be recovered from the DoFs.
• for A 2 (T ), pairings with M (T ) can be recovered from the DoFs.
Theorem 5.2 (FE for the stress complex).
• A 0 (T ) has dimension 12 and the provided DoFs are unisolvent.
• A 1 (T ) has dimension 15 and the provided DoFs are unisolvent.
• A 2 (T ) has dimension 6 and the provided DoFs are unisolvent.
• The complex A • (T ) is a resolution of the affine functions.
Proof. (i) For A 0 (T ) this is well known.
(ii) For A 2 (T ) one can consider the scalar analogue, and go via the dual result, that any u ∈ P 1 (T ) is uniquely determined by the integrals on the three small triangles in R(T ), since their isobarycentres are not colinear.
(iii) For A 1 (T ) this is proved in [53], also via the Clough-Tocher element. We provide a slight modification of that proof, exploiting the M DoFs (see Remark 5.2). The space P 1 (R(T ), S) has dimension 3 × 3 × 3 = 27. Imposing continuity of the normal component on interior edges can be expressed with 3 × 2 × 2 = 12 constraints, so dim A 1 (T ) ≥ 27 − 12 = 15. Now let u ∈ A 1 (T ) and suppose that its DoFs are 0. It then follows that div u ∈ A 2 (T ) has DoFs 0, by integration by parts. So div u = 0 and hence u = airy v for some v ∈ H 2 (T ). The second order derivatives of v are in P 1 (R(T )), hence v ∈ A 0 (T ). Now the M (E)-DoFs of airy v are 0, so there exists w ∈ P 1 (T ) such that v and w have the same M (V )-DoFs (see Lemma 5.1). We have u = airy(v − w). On a given edge ∂ τ ∂ ν (v − w) = 0 from the DoFs of u, so ∂ ν (v − w) is affine, in fact constant. Therefore v − w = 0 so u = 0.
(iv) Exactness of A • (T ) at index 1, was just proved. At index 0 the kernel is the space of affine functions. Exactness at index 2 then follows by dimension count.
Remark 5.3. For each edge E, let χ E be a nonzero affine map E → R such that χ E = 0.
In A 1 (T ) we may think of the provided DoFs attached to an edge E as • Integrals of uν ·ν against P 1 (E) and of uν ·τ against P 0 (E), which together constitute pairings with M (E).
From this point of view it seems natural to replace the edge DoF in A 0 (T ) (namely u → E ∂ ν u) by u → (∂ τ ∂ ν u)χ E . In particular (∂ τ ∂ ν u)χ E = 0 iff ∂ ν u is affine, which appeared as a step in the proof of Theorem 5.2.
Definition 5.2 (FES for the stress complex). We get a finite element system by appending the following spaces: A system of degrees of freedom is defined by: Proposition 5.3. The above finite element system is compatible and the system of degrees of freedom unisolvent.
Proof. The essential points were already proved in Theorem 5.2. The main addition is the unisolvence of the degrees of freedom on A 1 (E), which is straigthforward.
Remark 5.4 (Minimal spaces). The preceding FE complex may be reduced as follows: •Ã 0 (T ) is the reduced Clough-Tocher element where ∂ ν u is affine on edges. The DoFs are now only pairings with M (V ) at vertices. The dimension is 9.
•Ã 2 (T ) consists of divergence free elements of P 0 (R(T )), i.e. the elements that are continuous in the normal direction on interior edges. It can also be characterized as curl C 0 P 1 (R(T ), R). The dimension is 3 and the degrees of freedom are now pairings with M (T ) only. Figure 2: Diagram chase for the Johnson Mercier element.
•Ã 1 (T ) is the subspace of H 0 div P 1 (R(T ), S), consisting of elements u such that div u ∈Ã 2 (T ) and, on any edge, uν · τ ∈ P 0 (E). The dimension is 9 and the degrees of freedom are pairings with M (E) on edges, only.
A compatible FES is obtained by appending As already asserted, natural degrees of freedom are: Remark 5.5. Discrete diagram chase for the Johnson Mercier stress complex. See Figure 2.
One could also exhibit an alternative diagram chase, based on Remark 4.5.
Remark 5.6. High order finite element stress complexes are provided in [3]. They seem to fit in our framework too. Their lowest order complex starts with C 1 P 4 (R(T ), R) and ends in P 1 (T, V). The middle space of stresses augments P 2 (T, S) with a 3-dimensional space, so has dimension 21.
6 FES for the strain complex 6

.1 Induced operators and discrete vector bundle
We here consider now the complex (63). We identify induced spaces and operators on edges and vertices. They are summarized in the commuting diagram depicted in Figure 3. Again we consider a vertex V , in an edge E, in a triangle T . The oriented unit tangent on E is denoted τ and the normal is denoted ν, so that (τ, ν) is an oriented orthonormal basis of V. If we consider instead the complex of lower regularity, as in (77), the corresponding diagram is depicted in Figure 4.
The discrete vector bundle is defined as follows: • We let M (T ) be the kernel of the Airy operator, namely the space of affine functions.
• We let M (E) the space of pairs (u, v) where u is an affine function on E and v is constant.
We define restriction operators: which is independent of E. The pairings are: The Stokes-like identities are: • For u ∈ Γ(T, S) and φ ∈ M (T ) we have: • For (u, v, u ) ∈ Γ(E, R 2 × R) and (φ, ψ) ∈ M (E): Here the brackets denote differences between values at the two vertices of E. Proof. From Corollary 1.5. At index 0 the kernel of the discrete covariant exterior derivative has dimension 3. On the other hand the rigid motions naturally inject into this kernel. Surjectivity then follows from dimension equality.

Discrete spaces: higher regularity
For a triangle T we denote by R(T ) the Clough Tocher split of T .
Proposition 6.2. The provided DoFs give interpolators commuting with the differential operators.
The preceding finite element spaces can regarded as a finite element system, as follows. Definition 6.2 (FES for the strain complex).
• Differential operators and restrictions are defined according to Figure 3.
• The spaces A 0 (T ), A 1 (T ) and A 2 (T ) are defined as in Definition 6.1.
• Spaces on edges E and vertices V are defined by: • A system of degrees of freedom F on A is defined by: (other spaces are set to 0). Theorem 6.3 (FE for the strain complex).
• A 0 (T ) has dimension 24 and the provided DoFs are unisolvent.
• A 1 (T ) has dimension 24 and the provided DoFs are unisolvent.
• A 2 (T ) has dimension 3 and the provided DoFs are unisolvent.
• The sequence A • (T ) is a resolution of the rigid motions.
Proof. (i) For A 0 (T ) the dimension count, and the unisolvence of the degrees of freedom are standard. We also remark that the degrees of freedom corresponding to an edge E are unisolvent on A 0 (E). At vertices the corresponding result is trivial.
(ii) For A 2 (T ) the dimension count is trivial and the unisolvence of the DoFs is straightforward.
(iii) The space C 0 P 2 (R(T ), S) has dimension 10×3 = 30. For u ∈ C 0 P 2 (R(T ), S), in order to impose that sven u is integrable we impose that ∂ ν uτ ·τ is continuous on the three interior edges. Since these fields are linear, this can be expressed as 6 constraints. This shows that dim A 1 (T ) ≥ 24.
(iv) We now prove unisolvence of the degrees of freedom for A 1 (T ). Let u ∈ A 1 (T ) and suppose that the DoFs are all 0. Then sven u ∈ A 2 (T ) has 0 degrees of freedom, so it is 0. So we can choose v ∈ H 2 (T, V) so that def v = u. The second order derivatives of v can be recovered from the first order derivatives of u, and therefore turn out to be linear. Therefore v ∈ A 0 (T ).
Since the M (E)-dofs of def v are 0, there exists a rigid motion w, which has the same M (V )-dofs as w (ie value and value of curl), by Lemma 6.1. We notice that def(v − w) = u and proceed to show that v − w = 0, by showing that its degrees of freedom, as defined in A 0 (T ), are 0.
We have that v − w is 0 at vertices and that grad(v − w) is 0 at vertices (the symmetric part is def v = u and antisymmetric parts is essentially curl(v − w)).
-We have that E ∂ ν (v − w) · ν = 0 since this is one of the DOFs of u.
Remark 6.2. The preceding proof of unisolvence for A 1 (T ) was written from the point of view of Definition 6.1. From the point of view of FES, as in the extended Definition 6.2, one would go via Proposition 3.6, with similar arguments. That way yields the additional important information that degrees of freedom attached to edges are unisolvent on the space of restrictions to the edge, which guarantees the appropriate global continuity of the finite element fields defined piecewise by their DoFs. Remark 6.3. Discrete diagram chase for the finite element strain complex with high regularity. See Figure 6.
In the next section we will use the following consequence of Theorem 6.3: Proposition 6.4. For any v ∈ A 2 (T ) there is a unique u ∈ A 1 (T ) such that sven u = v, the restriction of u to ∂T is zero and ∂ ν uτ · τ is constant on each edge.
Proof. Indeed if we choose a constant c E for each edge E, the data (0, 0, 0, c E ) ∈ A 1 (E) is compatible at vertices, so can we extended to an element u of A 1 (T ), which is unique since there are no interior degrees of freedom. The Stokes identity then takes the form, for any affine φ on T : For any desired v = sven u ∈ A 2 (T ), this uniquely determines the coefficients c E of u.
Remark 6.4 (Minimal finite element strain complex with high regularity). One can get a minimal complex as follows.
We start with the modified reduced Clough Tocher space for A 0 (T ). Recall that one usually requires the normal derivative on edges to be affine. Instead we take the subspace of vectorfields u such that def u applied to the normal vector on edges is affine. The degrees of freedom are just vertex values and vertex values of the gradient.
For A 1 (T ) one takes the sum of the space def A 0 (T ) and the space defined in the preceding proposition, so that normal components on edges are affine. The degrees of freedom consisting of vertex values and pairing with M (E) for each edge E are then unisolvent. Then A 1 (T ) = P 2 (E) × P 1 (E) × P 1 (E) × P 1 (E).
The space A 2 (T ) is unchanged. The canonical DoFs give interpolators that commute with the differential operators.

Discrete spaces: lower regularity
In the following, for each edge E we let χ E be a nonzero affine map E → R such that χ E = 0. Notice that if u : E → R is affine and E uχ E = 0 then u is constant. Definition 6.3 (FE for the strain complex: low regularity).
• A 2 (T ) = P 0 (R(T ), R), the space of piecewise constants. The degrees of freedom are pairings with M (T ), namely integration against affine functions (3).
The finite element spaces of Definition 6.3 are represented in Figure 7.
• A 0 (T ) has dimension 15, and the provided degrees of freedom are unisolvent.
• A 1 (T ) has dimension 15, and the provided degrees of freedom are unisolvent.
• A 2 (T ) has dimension 3, and the provided degrees of freedom are unisolvent.
• The sequence A • (T ) resolves the rigid motions.
Proof. (i) We introduced this space in [33] (Proposition 3), but there we had different edge degrees of freedom (in particular we had u → E u · τ ). In any case this gives the dimension. We now check unisolvence. Let u ∈ A 0 (T ) have 0 DoFs. In particular u is 0 at vertices as well as curl u. Given preceding results it remains to be proved that u is 0 on each edge.
-We have that def uτ · τ = ∂ τ (u · τ ) is affine and orthogonal to χ E hence constant. So u · τ is linear on each edge, hence 0.
-We have that def uτ · ν = (1/2)(∂ τ (u · ν) + ∂ τ (u · ν)) affine and orthogonal to χ E , hence constant. Also ∂ τ (u · ν) − ∂ τ (u · ν) = curl t u = 0 on E, so ∂ τ (u · ν) is constant so u · ν = 0 on each edge. (iii) Unisolvence on A 1 (T ). Choose u ∈ A 1 (T ) with 0 degrees of freedom. We then get sven u = 0 from the Stokes identity, since the M (E) DoFs are zero. So we may choose v ∈ A 0 (T ) such that u = def v. We may find a rigid motion such that v and w have the same M (V ) degrees of freedom. Then def(v − w) = u moreover the vertex degrees of freedom of v − w are 0. The remaining edge degrees of freedom of v − w, are edge degrees of freedom of u, hence 0, so v − w = 0. Hence u = 0. Proposition 6.6. We get a compatible FES by appending the spaces: The degrees of freedom are now described as: Proof. Unisolvence of the edge DoFs on the edge spaces can be checked by similar arguments.
Remark 6.5 (Alternative definition of W (T )). One can replace the space W (T ) in Definition 6.3 by a construction with the Poincaré -Koszul operators. Indeed let x r be the the identity vector field, rotated by π/2, with respect to an origin located at the central vertex of R(T ). Consider the matrix field ω = x r (x r ) t . For any internal edge E, connecting the central vertex of R(T ) with one of the vertices of T , with tangent vector τ and normal vectror ν, the matrix field ω has the property that, on the edge, ωτ = 0 and ∂ ν ωτ · τ = 0. See Lemma 4.8. Furthermore on any (external) edge E of T , ωτ · τ ∈ P 0 (E), ωτ · ν ∈ P 1 (E) and ∂ ν ωτ · τ ∈ P 0 (E). Finally sven ω is constant.
We can define W (T ) = {uω : u ∈ P 0 (R(T ), R)}. Indeed such uω will have restriction 0 to internal edges, so the sven can be computed classically and is proportional to u. Remark 6.6. Discrete diagram chase for the finite element strain complex with lower regularity. See Figure 8. Remark 6.7 (minimal finite element strain complex with low regularity). We get a minimal complex by imposing all DoFs of the type, on edges E, integral against χ E to be zero. The dimensions are then 9 for the vector fields, 9 for the strain tensors and 3 for the scalars.