- Email: [email protected]
A uniform rationale for Whitney forms on various supporting shapes
Available online at www.sciencedirect.com
Mathematics and Computers in Simulation 80 (2010) 1567–1577
A uniform rationale for Whitney forms on various supporting shapes A. Bossavit LGEP (CNRS), 11 Rue Joliot-Curie, 91192 Gif-sur-Yvette CEDEX, France Received 5 November 2008; accepted 6 November 2008 Available online 25 November 2008
Abstract A method for constructing Whitney forms is proposed, which applies to tetrahedra, hexahedra, triangular prisms, and pyramids in a similar way, and proceeds from a unique generating principle, thus unifying their presentation. The principle automatically enforces conformity (i.e., “tangential” or “normal continuity” of the elementary proxy fields) at element interfaces, and generates a complex of forms with the “exact sequence property” built in. © 2008 IMACS. Published by Elsevier B.V. All rights reserved. PACS: 65N30 Keywords: Finite elements; Edge elements; Whitney forms
1. Introduction In the finite element practice for electromagnetism, vector-valued “edge elements” on tetrahedra are now universally accepted, and the advantages of thinking of them as belonging to a larger family of scalar- or vector-valued simplexrelated fields – the Whitney forms – are well known. But flexibility in solid modelling requires such elements on other shapes of solids than simplices: prisms, hexahedra, pyramids. (It can easily be seen, for instance, that a transition layer of elements between a simplicial mesh and a brick mesh must include pyramids [7].) There is a rationale about “why Whitney forms are what they are” on simplices [2](which one may think Whitney himself knew about, although this is not explicitly stated in [9]). But this rationale does not apply directly to hexahedra and prisms, for which the known constructions [8] seem to be of a different flavor. Worse, for the pyramid, the constructions described in the few published papers (e.g., [3,4]) seem to stem from specific heuristic ideas, not obviously related with those used for tetrahedra. It is therefore desirable to derive all Whitney forms – edge-based, face-based, etc., for tetrahedra, prisms, etc. – from a unique generating principle. To this effect, I present a recursive procedure, consisting in the iteration of two basic geometric transforms, “conation” and “extrusion”, both able to “lift up” to dimension n + 1 a given n-dimensional Whitney complex. Hence (by the recursive formulas (8)–(12) below), Whitney forms on all kinds of n-topes that can be built from a point by iterated conation or extrusion. In dimension 3, this gives Whitney forms on the tetrahedron, hexahedron, triangular prism, and pyramid, allowing the practitioner to freely mix such shapes, each cell (node, edge, face,. . .) of this mixed mesh being assigned its own Whitney form, whose integral is 1 over this cell and 0 over all others. Their compatibility at interfaces, the property also known as “conformity” (or “tangential continuity” in the
E-mail address: [email protected]. 0378-4754/$36.00 © 2008 IMACS. Published by Elsevier B.V. All rights reserved. doi:10.1016/j.matcom.2008.11.005
1568
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
case of vector-field representation of 1-forms), is automatically enforced. Moreover, the “exact sequence property” of the whole complex of forms, so useful in the avoidance of spurious modes in resonant cavity problems [1], also comes naturally. 2. Basic notions Some notation (italics, when not used for emphasis, warn about implied definitions): we work in A∞ , the ∞dimensional real affine space, whose associated vector space is denoted V∞ , and we deal with subspaces and submanifolds (with boundary, most often) of finite dimension p ≥ 0. Symbol x and similar ones (y, z, etc.), refer to points, symbols v, w, etc., denote vectors. If v translates x to y, we may write y = x + v. Barycenters are denoted by tx + (1 − t)y, where x and y are points and t a real number, or by x + tv. We use the Grassmann algebra of multivectors (symbol ∨ for the product) and multicovectors (symbol ∧). Vectors (and by extension, all Grassmann objects) can be free (i.e., elements of V∞ ) or bound (i.e., anchored at some point x of A∞ ). Singular p-chains are formal weighted sums, with real coefficients, of p-submanifolds. Cochains, or differential forms, are the dual objects. We shall denote by v; ω the pairing between vector v and covector ω, andmore generally between p-vector and p-covector, or chain and cochain. In particular, when referring to the integral M ω of a form ω over a manifold M, we shall prefer the alternative notation M; ω. The boundary operator ∂ maps p-chains to (p − 1)-chains. The exterior derivative d is defined as its dual, acting on cochains: ∂M; ω = M; dω. Let be given a cell-complex Kn of finite dimension n, with compact carrier |Kn |, borne by an n-dimensional affine subspace. (Precise definitions will be found in Topology courses such as [6], but the notion is familiar: finite element meshes, provided they have no “hanging” nodes or edges, are cell complexes, whose carrier is the meshed region.) Since our purpose is to lift Kn , somehow, to dimension n + 1, we shall refer to it as the horizontal complex, and use underlined symbols for its elements. Similar symbols, not underlined, will serve for elements of the enlarged, vertical complex Kn+1 , of dimension n + 1, that will result from our construction. The complex Kn is made of separately and arbitrarily oriented p-cells, each the homeomorphic image of an open p-ball, if p > 0, and 0-cells are points (whose orientation is the number 1 or −1). The generic cell is denoted by s or σ, with the convention that dim(σ) = dim(s) − 1. Cells are related by incidence numbers, forming incidence matrices p , with 0 < p ≤ n: thus, sσ p = ±1 if the (p − 1)-cell σ belongs to ∂s, the sign depending on whether orientations of σ and s match or not. We denote by Sp the set of oriented p-cells, with p ≤ n. No generality will be lost by assuming there is only one n-cell, carried by |Kn |, as in Fig. 1. To Kn is naturally associated a chain complex {C, ∂}, by taking the restriction ∂ of ∂ to cell-based chains (i.e., linearcombinations of cells with real sweights). Note the relation between and ∂: the boundary of the p-chain c = s ∈ S cs s is σ ∈ S , s ∈ S sσ p c σ. So if c (boldface) stands for the array of weights {cs : s ∈ Sp }, one has p−1 p p ∂c = σ ∈ S (p c)σ σ. More concisely said, the family of maps κp defined, for each p, by c → s ∈ S cs s, is a p−1 p chain map: ∂κp = κp−1 p . We now describe the generation of Kn+1 by conation or extrusion of Kn . How these operations affect cochains, especially Whitney forms, will come later.
Fig. 1. Notations, and conventions about orientation. On the left, the simplex s is the face in the foreground.
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1569
First, conation. Let us select a point a (for “apex”) outside the subspace supporting Kn , and consider the cone |Kn+1 | = {λa + (1 − λ)x : x ∈ |Kn |, λ ∈ [0, 1]} that will carry the conation of the complex. More generally, if M is a submanifold, not necessarily in |Kn |, one sets cone(a, M) = {λa + (1 − λ)x : x ∈ M, λ ∈ [0, 1]}. (To avoid irrelevant difficulties, we assume that no half-line {a + tv : t ≥ 0}, where v is a nonzero vector, will intersect M more than once. The apex a will be understood and omitted from now on.) If M is oriented, this induces an orientation of cone(M), the one such that (Fig. 1) the following equality between singular chains holds: M = ∂[cone(M)] + cone(∂M).
(1)
(By convention, cone(∅) = {a}. When M reduces to a point x, with orientation +, this implies that cone(x) is the segment ax, oriented from a to x.) If x = / a, the point where ax intersects |Kn | is the projection πx of x. This induces a projection operation on other objets, for instance bound vectors anchored at point x, the projection of v being π∗ v. The complex Kn lifts up as follows: each of its p-cells σ has an associated so-called vertical cell s = cone(σ), carried sσ by the set {λa + (1 − λ)x : x ∈ σ, λ ∈ ]0, 1[}, and one may define new incidence numbers by sσ p = 1 and p = 0 for other cells σ in Sp . If cells σ and s of Kn are incident, their cones are incident too, and orientations match or not the same way, so one will have cone(s)cone(σ) = sσ p p . The horizontal cells of Kn+1 are just those of K n , with the sσ sσ same incidence numbers, p = p , so p is embedded in p as a submatrix. Last, a itself forms a new 0-cell, and sa 0 = −1 for all vertical 1-cells s. The set of p-cells in Kn+1 is denoted Sp . Extrusion, now. This time, a fixed vector X, not parallel to the subspace supporting |Kn |, is given. Setting ϕ(t, x) = x + tX, one defines the extrusion extr(|Kn |) as the set |Kn+1 | = {ϕ(λ, x) : x ∈ |Kn |, λ ∈ [0, 1]}. (More generally, one might have a smooth, singularity-free vector field X, with flow ϕ.) A point x ∈ |Kn+1 | is ϕ(λ, x) ¯ = x¯ = ϕ(1, x) besides its bottom projection for some λ ∈ [0, 1] and some x ∈ |Kn |, so it has now a top projection πx πx = x. A p-vector v at x projects, similarly, to π∗ v and π¯ ∗ v. We shall denote by rear(x) the segment x¯x and define rear(M) as the set-union {rear(x) : x ∈ M} in the case of a submanifold M in |Kn+1 | to which X is nowhere tangent. If ¯ by continuity, and rear(M) as well, in such a way M is oriented, this orients its bottom and top projections M and M, that ¯ M = ∂[rear(M)] + rear(∂M) + M.
(2)
To a cell σ in Kn correspond now two horizontal cells in Kn+1 , namely σ itself and its image ϕ(1, σ), which we denote ¯ (Thus horizontal cells in Kn+1 form two disjoint complexes, Kn on the one hand, and its translate by X, that one by σ. ¯ n , on the other hand.) Vertical cells are extr(σ) for σ ∈ Sp . Again, p is embedded in p , with minor may denote by K extr(σ)σ
variations: incidence numbers between horizontal cells are kept unchanged. We set p
s σ and p = sp σ¯
=0
for vertical cells s
extr(σ)σ¯
= 1 and p
= −1
other than extr(σ). Last, if cells σ and s of Kn are incident, their extrusions are extr(s)extr(σ)
sσ
incident too, and orientations match or not the same way, so one will have p = p . Let us now introduce J. Harrison’s chainlets, which are elements of the completion of a space of chains with respect to an appropriate norm, the natural norm [5]. Some passages to the limit, involving Cauchy sequences in the sense of this norm, can then profitably be associated with an element of the Banach space of chainlets. We must refer to [5] for details, and just mention results of this theory that will play a role here. For instance, a vector at point x can be considered as the chainlet limit of the sequence of chains cn = n [x, x + v/n] (i.e., the segment from point x to point x + v/n, with weight n), which allows us to define cone(v) as the chainlet limit of cone(cn ). Hence a hybrid object, neither bivector nor surface, but yet a bona fide chainlet, that we shall call a needle of dimension 2, or 2-needle. A similar treatment applies to a p-vector, seen as the chainlet limit of a sequence of smaller and smaller p-dimensional patches with higher and higher weights, the cone of which is a (p + 1)-needle. In the same spirit, a p-manifold M can be approached, in the sense of p-chainlets, by a Cauchy sequence of finite sums of (bound) p-vectors, each of which approximates a p-dimensional chunk of M. The corresponding free p-vectors can be added, and the limit of this sum is a (free) p-vector vect(M), that we shall call the vectorial p-volume of M. This operation extends to chains, and to chainlets. (One may see vect(c) as the map ω → c; ω, where ω is on the left a p-covector, on the right the constant p-form equal to it, since such a map can be identified, by duality, with a
1570
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
p-vector.) One checks that vect(v), if the bound vector v is considered as a chainlet, is the free vector v. Note also that vect(∂c) = 0, whatever c. We shall have use for the following relation, where v is a p-vector anchored at a point x ∈ |Kn+1 |, at height λ < 1 (Fig. 1), with horizontal projection πx relative to the apex a: vect(cone(v)) = (1 − λ)p+1 vect(cone(π∗ v)).
(3)
This is easily proven, for p = 1, by working on the n th term of a Cauchy sequence of 1-chains with chainlet limit v, for instance the segments [x, x + v/n] with weight n, as above. One then compares the areas of the two triangles cone([x, x + v/n]) and cone([x, x + π∗ v/n]), to see that their ratio tends to (1 − λ)2 . Generalizing to p > 1 is easy. In the case of extrusion, there is a similar relation: vect(rear(v)) = (1 − λ)vect(extr(π∗ v)), but nothing that simple will hold if the extrusion field X is not constant. 3. Heuristic moves Consider the known context of Whitney forms on a simplicial mesh: to each p-simplex s is associated a p-form well
ws such that s ws ≡ s; ws = δss (the Kronecker delta). These forms, traditionally, are understood as interpolants, like finite element shape functions. Starting from a smooth p-form a (the mathematical representation of some physical field, such as a vector potential if p = 1), one may build an array{as : s ∈ Sp } of numbers, one per p-cell, by setting Whitney forms can then be used to generate a p-form s as ws , which is an approximation of a in the sense as = s; a. that c; a − s as ws = 0 for, not all p-chains c of course, but for all simplicial p-chains. (To check this, set c = s and
use s; ws = δss .) The finer the mesh, the better this approximation, and one may base error analysis on that. We shall look at Whitney forms, here, equivalent, dual viewpoint: as approximants of chains, rather than of cochains. from an For this, note that c; s as ws = s ∈ Sp c; ws s; a. So if one sets whit(c) = s ∈ Sp c; ws , this simplicial chain is an approximation of c, in the sense that c − whit(c); a = 0 for, not all p-forms a, but all those in the space W p
generated by the ws . (To check this, set a = ws and use s; ws = δss .) Whitney forms are thus a device that, starting from a singular p-chain (or chainlet) c embedded in the carrier of a cell complex, delivers a cellular p-chain whit(c), meant to approximate c. Generating Whitney forms (should we not know them) thus amounts to build the map whit. As this section will show, this map is determined by the following requirements:
(1) That whit(s) = s (which is the same as s; ws = δss ). s (2) That it be a chain map: whit(∂c) = ∂(whit(c)) (which entails, by duality, dwσ = s ∈ Sp sσ p w , for σ in Sp−1 , the very property that helped eliminate spurious modes [1]). (3) That whit preserve vectorial volume: vect(whit(c)) = vect(c) (by duality, this means that W p contains all constant p-forms, also a useful feature for error analysis). (4) That it behave naturally under conation, as explained below. (Conditions (1) and (2) would not suffice, as the following counter-example shows: In dimension 1, with 0-cells {0} and {1} and 1-cell ]0, 1[, the 0-forms (1 − x) + αx(1 − x) and x − αx(1 − x) would make, together with the 1-form [1 − α(1 − 2x)] dx, a Whitney complex according to these criteria, whatever α, whereas α = 1/2 for the standard one.) Back to the situation of the previous section, let us assume that a complex of Whitney forms is already known on Kn , and that spaces W p they span have the exact sequence property: dW p ⊂ W p+1 and ker(d; W p ) = dW p−1 for p ≥ 1. So we know how to map a horizontal singular chain c (i.e., one embedded in |Kn |) to a horizontal cell-chain whit(c). Our purpose is to build whit when c is embedded in the larger region |Kn+1 |. This should yield Whitney spaces W p for Kn+1 , with a natural embedding W p ⊂ W p . In two particular cases, whit is forced on us. It will then be seen that no choice is left about its general definition. The first case is c = cone(σ) ≡ s, for a horizontal cell σ. Condition (1) then imposes whit(cone(σ)) ≡ whit(s) = s ≡ cone(σ) = cone(whit(σ)), so whit and cone commute when applied to horizontal (i.e., embedded in Kn ) cellular chains. This suggests to define whit(cone(c)) as cone(whit(c)) for all horizontal c. This is the “naturality” alluded to above, point (4). The second case is c = cone(v), a (p + 1)-needle based on the p-vector v anchored at a point x ∈ |Kn+1 |, at height λ < 1, with horizontal projection πx. The projection π∗ v of v, a p-vector anchored at πx, is a chainlet with respect to the horizontal complex Kn . So whit(π∗ v) is a definite horizontal p-chain, namely s ∈ S π∗ v; ws (πx) s, where ws p
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1571
denotes the form in W p associated with s. According to the third condition above, and taking relation (3) into account, we require that whit(cone(v)) = (1 − λ)p+1 cone(whit(π∗ v))
(4)
in order to preserve vectorial volume. Let us unfold the implications. For a p-vector v, (4) entails whit(cone(v)) = (1 − λ)p+1 π∗ v; ws (πx) cone(s) s ∈ Sp
=
v; (1 − λ)p+1 (π∗ ws )(πx) cone(s)
(5)
s ∈ Sp
where π∗ denotes the pullback operation. (By definition, π∗ v; ω = v; π∗ ω for all v, ω.) Here, λ had its value at the anchor point of v, and the bracket paired a p-vector with a p-covector. But (1 − λ)p+1 π∗ ws makes sense as a p-form living on |Kn+1 |, with now λ a function of x. So we can replace v by a p-chain c and write whit(cone(c)) = c; (1 − λ)p+1 π∗ ws cone(s), s ∈ Sp
with the bracket now an integral over c (that is to say, the limit of a Riemann sum of terms similar to the bracket in (5)). The same formula applied to ∂c yields (the dummy variable changes from s to σ, according to our general convention) whit(cone(∂c))
=
∂c; (1 − λ)p π∗ wσ cone(σ)
σ ∈ Sp−1
=
c; d[(1 − λ)p π∗ wσ ]cone(σ)
σ ∈ Sp−1
so finally, using (1), whit(c) = whit(∂[cone(c)]) + whit(cone(∂c)) = ∂[whit(cone(c))] + whit(cone(∂c)) = c; (1 − λ)p+1 π∗ ws ∂(cone(s)) + c; d[(1 − λ)p π∗ wσ ]cone(σ). s ∈ Sp
σ ∈ Sp−1
(6)
On the other hand, whit(c) = s ∈ Sp c; ws s, where the p-cells s belong to Kn+1 , so we shall find what ws must be by the tedious process of comparing this with (6). Cells s in Sp are of two kinds: Horizontal ones, which already belonged to Sp , and vertical ones, of the form s = cone(σ), where σ spans Sp−1 , so we deal with them separately. (Cf. Fig. 2 for mnemonic aid.) For a horizontal cell s, the only contribution of the right-hand side of (6) to c; ws is c; (1 − λ)p+1 π∗ ws , since ∂(cone(s)) = s + (. . .), where the dots stand for a chain of other, vertical, cells. So we find ws = (1 − λ)p+1 π∗ ws for horizontal p-cells of Kn+1 . (Note that the“horizontal” Whitney form ws lives on |Kn+1 |, and that ws is its restriction
Fig. 2. Current notation.
1572
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
to |Kn |.) This helps simplify (6) a bit: c; ws ∂(cone(s)) + c; dwσ cone(σ). whit(c) = s ∈ Sp
(6 )
σ ∈ Sp−1
We also observe that, for a horizontal (p − 1)-cell σ, dwσ
= =
d[(1 − λ)p π∗ wσ ] = −p(1 − λ)p−1 dλ ∧ π∗ wσ + (1 − λ)p π∗ dwσ p ∗ s −p(1 − λ)p−1 dλ ∧ π∗ wσ + sσ p (1 − λ) π w , s ∈ Sp
since the exact sequence property for Kn implies dwσ = we have s (1 − λ)dwσ = −p dλ ∧ wσ + sσ p w ,
sσ s s ∈ Sp p w .
So, for horizontal (p − 1)-forms on Kn+1 , (7)
s ∈ Sp
a formula which will prove useful immediately.
For a vertical cell s = cone(σ), factors of cone(s) in (6 ), in addition to the term c; dwσ that comes from the second
sum on the right of (6 ), will also come from the first sum, on behalf of cells s ∈ Sp incident on σ. Specifically, since
∂s = σ ∈ S sp σ σ , one has p−1
∂(cone(s )) = s − cone(∂s ) = s − sp σ cone(σ ), σ ∈ Sp−1
from which one gets, using (6 ) then (7),
s σ σ σ c; ws = c; dwσ − σs p c; w = c; dw − (1 − λ)dw − p dλ ∧ w . s ∈ Sp
So, finally, ws = λ dwσ − p dλ ∧ wσ for a vertical p-cell s = cone(σ). We could now repeat all these moves in the case of extrusion. This is straightforward. One will notice however, in so doing, that the extrusion vector field must be constant for vectorial volume to be preserved. 4. Whitney forms on simplices and bricks Our heuristic explorations have thus led to generating formulas. In the case of conation, (a) For 0 ≤ p ≤ n, and s ∈ Sp , the p-form ws , living on Kn , is lifted to Kn+1 by the formula ws = (1 − λ)p+1 π∗ ws ,
(8)
hence the horizontal p-form ws . (b) For 1 ≤ p ≤ n + 1, and σ ∈ Sp−1 , the vertical p-form associated with s = cone(σ) is ws = λ dwσ − p dλ ∧ wσ .
(9)
(c) The apical 0-form, or height of point x, denoted λ(x), in W 0 , completes the inventory. In the case of extrusion, λ being the function such that x = ϕ(λ(x), πx), (a) For 0 ≤ p ≤ n, and s ∈ Sp , each p-form ws , living on Kn , spawns two horizontal p-forms, the bottom one: ws = (1 − λ) π∗ ws ,
(10)
and the top one: wϕ(1,s) = λ π∗ ws .
(11)
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1573
(b) For 1 ≤ p ≤ n + 1, and σ ∈ Sp−1 , the vertical p-form associated with s = extr(σ) is ws = −dλ ∧ π∗ wσ ≡ −dλ ∧ (wσ + wϕ(1,σ) ),
(12)
as seen by adding (10) and (11). Iterated conation under Rules (a)(b)(c) generates the classical Whitney forms for any degree and spatial dimension, starting from the virtual nothing that K0 is. (When n = 0, the complex K0 reduces to a single 0-cell σ and its 0-form wσ = 1.) Let us formally prove this. First, any p-form, 0 ≤ p ≤ n, in the complex Kn+1 over the (n + 1)-simplex |Kn+1 |, does come, by rule (a), from a p-form of the complex Kn over some n-dimensional face of |Kn+1 |. To prove this, start from the known generic expression of Whitney forms [9], wς = p!
p
(−1)i λς(i) dλς(0) ∧ dλς(1) ∧ · · ·(i)· · · ∧ dλς(p) ,
(13)
i=0
where ς is any injective map from the integer segment [0, p] into the set of nodes of Kn+1 , and (i) means “omit the term indexed by i”. (Thus, ς labels an oriented p-face, and there is a Whitney form for each of these.) Suppose p ≤ n. There is then at least one node distinct from any of the ς(i)’s, that we take as cone’s tip and call a, denoting by λ its nodal function and by π the projection from a to the opposite n-face, the n-simplex |Kn |. Call λς(i) the restriction of λς(i) to the latter. Then λς(i) = (1 − λ)π∗ λς(i) . Since d commutes with pullback, d(π∗ λς(i) ) = π∗ (dλς(i) ). Then dλς(i) = (1 − λ)π∗ (dλς(i) ) − π∗ (λς(i) dλ). Substituting in (13), one gets two parts, one with dλ factored out where terms cancel two by two, and another part that is p! (1 − λ)p+1
p
π∗ ((−1)i λς(i) dλς(0) ∧ · · ·(i)· · · ∧ dλς(p) ).
i=0
(1 − λ)p+1 π∗ (wς ),
= confirming (8). It remains to be proven that the top form (p = n + 1) in Kn+1 comes from the top form (p = n) of some Kn by rule (b). So here, s is the (n + 1)-simplex |Kn+1 | and σ is one of its n-faces. Let’s label the nodes of σ by 0, 1, . . . , n and the corresponding barycentrics by λi , while λ is the barycentric of the node of s thus left out. The restriction of λi to σ is λi , and λi = (1 − λ)π∗ λi . The top form wσ in Kn is (up to sign) n! dλ1 ∧ · · · ∧ dλn and lifts to wσ = (1 − λ)n+1 π∗ wσ . Since the degree of wσ is n, one has dwσ = 0, hence d(π∗ wσ ) ≡ π∗ (dwσ ) = 0, so dwσ = −(n + 1)(1 − λ)n dλ ∧ π∗ wσ . Using (9), with p = n + 1, we obtain So
wς
ws
= −(n + 1)(1 − λ)n [λ dλ ∧ π∗ wσ + (1 − λ) dλ ∧ π∗ wσ ] = −(n + 1)(1 − λ)n dλ ∧ π∗ wσ
−(n + 1)! (1 − λ)n dλ ∧ π∗ (dλ1 ∧ · · · ∧ dλn ) = −(n + 1)! dλ ∧ i=1,...,n (1 − λ) d(π∗ λi ), i ∗ i where the expression headed by the symbol is the wedge product ofi its terms. Since dλ = d[(1 − λ) π λ ] = ∗ i ∗ i s (1 − λ)π (dλ ) − π λ dλ, we see that w = −(n + 1)! dλ ∧ i=1,...,n dλ , which is indeed the top form in Kn+1 , up to sign. So iterated conation, with rules (8) and (9), does generate the simplicial Whitney forms. We may expect iterated extrusion to generate Whitney forms on a hypercube. Verifying this fact is a lot easier than what precedes, as we presently see. Let v1 , . . . , vn be independent vectors (in the infinite-dimensional real vector space V∞ ), and consider the parallelotope built (in the associated affine space A ∞ ) on v1 , . . . , vn , anchored at some origin point o. Call this |Qn |. A natural representation of its points is x = o + i=1,...,n xi vi , where 0 ≤ xi ≤ 1 ∀i. (Here one ought to distinguish x, a point in A∞ , from the real n-tuple x = {xi : i = 1, . . . , n} of its coordinates. Denoting both entities by the same symbol x should not be too confusing in what follows.) The Whitney complex Qn over |Qn | is known [8] to consist of the following elements: =
1574
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
Fig. 3. Labelling scheme and orientations, for n = 3. Under each label, the Whitney form associated with this p-face. For better readability, coordinates x1 , x2 , x3 are here x, y, z. The top form is dx ∧ dy ∧ dz, associated with the cube, the label of which is ••• . 1−k
k
• The 0-forms wk (x) = i=1,...,n (xi ) i (1 − xi ) i ≡ pk (x) where ki = 0 or 1. Here, k is the integer i=1,...,n 2i−1 ki , whose binary representation is obtained by concatenating k1 , k2 , . . . , kn , which makes a convenient label for the 2n nodes of the n-cube and for the 0-forms wk . k 1−k • The 1-forms wk,j (x) = i =/ j (xi ) i (1 − xi ) i dxj , associated with edges parallel to vj . ki 1−ki k,j i i • The 2-forms w (x) = i∈J dxj1 ∧ dxj2 / (x ) (1 − x ) where J is the subset {j1 , j2 } of [1, n], with j1 < j2 , which determines a set of faces, those parallel to the 2-vector vj1 ∨ vj2 . • And so forth, according to a pattern that Fig. 3 should help understand: let us denote by σ an increasing injection from [1, p] into [1, n], if p ≥ 1, and by Jσ its image. (If p = 0, set Jσ = ∅.) This way, each σ labels the set of p-faces parallel to the p-vector vσ(1) ∨ · · · ∨ vσ(p) . Now, each of these 2n−p faces is labelled by an integer k, as above (with some harmless redundancy, since k and k point to the same p-face when kj = kj ∀j ∈ Jσ ). Thus the pair {k, σ} labels a definite p-face. We shall set dxσ = dxσ(1) ∧ dxσ(2) ∧ · · · ∧ dxσ(p) , that is to say, dxσ = dxj1 ∧ dxj2 ∧ · · · ∧ dxjp where j1 < j2 < · · · < jp compose the set Jσ , and ki
1−ki
i i pk,σ (x) = i∈J / σ (x ) (1 − x )
.
(14)
This polynomial, of order n − p, is linear with respect to each component xi taken separately (“n-linear” polynomial). Note that it coincides with the above-mentioned pk , i.e., the nodal function, when Jσ is empty. With such notation, and the convention that dxσ is the scalar constant 1 when Jσ is empty, all p-forms of the complex, 0 ≤ p ≤ n, can be written wk,σ (x) = pk,σ (x) dxσ .
(15)
Fig. 3 illustrates all that in dimension 3. Note the use of dots to place the set Jσ in [1, n]. The exact sequence property is now easy to establish: Applying d to (15) gives ∂i pk,σ dxi ∧ dxσ , dwk,σ = i∈J / σ
where ∂i is differentiation with respect to xi . After (14), ∂i pk,σ does depend on xi any longer, and dxi ∧ dxσ = ±dxs , not k,s k,σ where s maps onto the set Js = Jσ ∪ {i}. Therefore dw = s ± p dxs , where s spans the set of injections from / Jσ . [1, p + 1] to [1, n + 1] whose images are the sets Js . There is one such s for each i ∈ Now, let us add another independent vector vn+1 , call xn+1 what was λ in formulas (10)–(12), and extrude the k 1−k previous complex via vn+1 . A factor (xn+1 ) n (1 − xn+1 ) n appears when using (11)–(13) throws in a new term dxn+1 . The pattern (14) and (15) is thus respected, so extrusion of Qn does yield the Whitney complex Qn+1 on the (n + 1)-cube |Qn+1 |.
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1575
Parallelotopes thus obtained are not convenient, however, for finite element practice, which may require distorted (but still, stackable) hexahedra. For this, one may use the well known “isoparametric” trick, which consists in mapping the cube vertex labelled by k to some point uk in A∞ , with interpolation by the nodal functions. Point x then goes to u(x) = pk (x)uk . (16) k
Under some mild restrictions about how to place the points uk , this map will be injective, with a smooth inverse u−1 , that makes a convenient chart from u(|Qn |) to Rn . To get the Whitney complex on u(|Qn |), we just collect the forms ∗ (u−1 ) wk,σ obtained by pullback, using this chart. There is one map v u like u in (16) for each elementary volume, with −1
shared subsets of vertices uk . Maps (v u)−1 and (v u) thus coincide on p-faces common to volumes v and v in the actual mesh, and so do the restrictions of their pullbacks, which guarantees stackability and conformity. A well-known drawback of this isoparametric approach is that constant p-forms do not belong to the Whitney space thus obtained, for p ≥ 2 : this results from losing vectorial volume conservation. 5. Exact sequence property Let us now verify from scratch the fact that both processes, whatever the kind of complex they act on, provided it has the exactness property, generate a complex that does inherit this property. All we have to show in fact (thanks to the 1–1 correspondence between p-forms thus built and p-cells) is that each dwσ , where σ is a (p − 1)-cell, is a linear combination of p-forms ws of the conation or extrusion, as the case may be. First, conation. Call W p and W p the spaces generated by the p-forms of the vertical and the horizontal complexes. Assuming that dW p−1 ⊂ W p , we need to show that dwσ ∈ W p when wσ belongs to W p−1 . There are two cases, depending on whether σ is a horizontal or vertical cell. (Cf. Fig. 2 again.) For a horizontal cell σ ∈ Sp−1 , use (8), i.e., wσ = (1 − λ)p π∗ wσ . Then dwσ = −p (1 − λ)p−1 dλ ∧ π∗ wσ + (1 − λ)p π∗ (dwσ ), so (1 − λ) dwσ = −p dλ ∧ wσ + (1 − λ)p+1 π∗ (dwσ ), and dwσ = λ dwσ − p dλ ∧ wσ + (1 − λ)p+1 π∗ (dwσ ). After (9), λ dwσ − p dλ ∧ wσ is ws for the cell s = cone(σ), so we have dwσ = ws + (1 − λ)p+1 π∗ (dwσ ).
By hypothesis, dwσ = s ∈ S sp σ ws . Therefore, carrying on, p
dwσ = ws +
s ∈ Sp
sp σ (1 − λ)p+1 π∗ ws = ws +
s ∈ Sp
sp σ ws
by using (8) again for p-cells s . So dwσ ∈ W p . For a vertical cell s, we have s = cone(σ) for some σ ∈ Sp−1 . First, a useful lemma: Lemma. For s vertical, i.e., s = cone(σ), one has dws = (p + 1)dλ ∧ dwσ ,
λ dws − (p + 1)dλ ∧ ws = 0.
Proof. By the recursive formula, ws = λ dwσ − p dλ ∧ wσ , so dws = dλ ∧ dwσ + p dλ ∧ dwσ = (p + 1)dλ ∧ dwσ . Next, dλ ∧ ws = dλ ∧ [λ dwσ − p dλ ∧ wσ ] = λ dλ ∧ dwσ , and hence, λ dws − (p + 1) dλ ∧ ws = 0.
1576
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
Fig. 4. Labelling scheme and orientations for the pyramid, as seen from above (•• 0 is the bottom square). Under or below each label, the Whitney form associated with this p-face. On the right, the base square’s 2-form (above •• 0) and the pyramid’s 3-form (above ••• ).
Now, since all (p + 1)-faces on which s is incident are vertical, we may suspect that dws =
s σ cone(s ) , s ∈ Sp p w
give or take a sign, so let’s evaluate this sum, using (9) to express wcone(s ) :
sp σ wcone(s ) = sp σ [λ dws − (p + 1) dλ ∧ ws ] s ∈ Sp
s ∈ Sp
⎛
=
λd⎝
s ∈ Sp
⎞
sp σ ws ⎠ − (p + 1) dλ ∧
s ∈ Sp
sp σ ws
=
λ d(dwσ − ws ) − (p + 1) dλ ∧ (dwσ − ws ) = −(p + 1) dλ ∧ dwσ − [λ dws − (p + 1) dλ ∧ ws ]
which is −dws indeed, thanks to the Lemma. So dws ∈ W p , and this ends the proof of exactness. σ Extrusion, now: start from wσ = (1 − λ)π∗ wσ for a (p − 1)-cell σ (cf. (10)). Then dw ∧ π∗ wσ + (1 − = −dλ
σ s
∗ σ s ∗ σ σ s λ) d(π w ), which is w + (1 − λ)π dw after (12), with s = extr(σ). Again, dw = s ∈ S p w , and ws = p
(1 − λ)π∗ ws by using (10), so we do have dwσ = ws + s ∈ S sp σ ws . p
6. Whitney forms for other shapes: pyramid, prism On a point o and two independent vectors vx and vy , let us build, by two successive extrusions, the complex K2 made of • The 0-forms (1 − x)(1 − y), (1 − x)y, x(1 − y), xy, on the 0-cells 00, 01, 10, 11, respectively; • The 1-forms (1 − y) dx, y dx, (1 − x) dy, x dy, on the 1-cells • 0, • 1, 0• , 1• , respectively; • The 2-form dx ∧ dy on the horizontal parallelogram, the 2-cell •• in this notational scheme. They generate a Whitney complex {W 0 , W 1 , W 2 }, as we have verified earlier. Now pick an out-of-plane point a, and make it the apex of a pyramid based on •• . For the sake of notational uniformity, let us denote by z what was formerly λ. Points p of this pyramid are parameterized as follows: p(x, y, z) = za + (1 − z)(x vx + y vy ), 0 ≤ x ≤ 1, 0 ≤ y ≤ 1, 0 ≤ z ≤ 1.
(17)
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1577
Of course, these are not Cartesian coordinates (Fig. 4). The projection π is given by π(p(x, y, z)) = x vx + y vy : note the invariance of coordinates x and y by projection, which makes this chart a convenient one for pullback. As an alternative to (17), we have p(x, y, z) = (1 − z){x[(1 − y)10 + y 11] + (1 − x)[(1 − y)00 + y 01]} + za, which better shows the weights of p with respect to each of the five nodes. The conified complex K3 has horizontal 0-cells 000, 010, 100, 110, 1-cells • 00, • 10, 0• 0, 1• 0, and 2-cell •• 0, according to an obvious nomenclature (Fig. 4), whose conations form the vertical 1-cells 00• , 01• , 10• , 11• , 2-cells • 0• , • 1• , 0•• , 1•• , and 3-cell ••• (the pyramid itself), to which one must add the apical 0-cell a. Next, we apply the recipe (8), with z substituted for λ, then find the vertical forms using (9). This tedious process, details of which are left out, results in Fig. 4. The top form (p = 3) comes from σ=•• 0, with wσ (p(x, y, z)) = (1 − z)3 dx ∧ dy. Then s ≡ cone(σ)is ••• , the pyramid itself, and ws (p(x, y, z)) =
z d[(1 − z)3 dx ∧ dy] − 3(1 − z)3 dz ∧ dx ∧ dy
=
−3z(1 − z)2 dz ∧ dx ∧ dy − 3(1 − z)3 dz ∧ dx ∧ dy
=
−3(1 − z)2 dx ∧ dy ∧ dz
which is the expected volume form for these coordinates. (The “wrong” sign is due to our general orientation rules.) The apical 0-form, wa (p(x, y, z)) = z, completes the list (Fig. 4). A word, finally, on the triangular prism. The 3D simplex was generated from naught by conation (or extrusion: the same thing in dimension 1), conation, conation, the brick by conation, extrusion, extrusion, the pyramid by conation, extrusion, conation. The prism stems from conation, conation, extrusion. This is left as an easy exercise. For both the pyramid and the prism, warped shapes are allowed if one proceeds by mapping and pullback, as in the case of bricks. The same caveat applies about the representation of constant fields. Let us finally remark that, whatever the volume, traces on faces only depend on the face itself (be it a triangle or a quadrilateral). Thus conformity is built in. References [1] A. Bossavit, Solving Maxwell’s Equations in a Closed Cavity, and the Question of Spurious Modes, IEEE Trans. MAG-26 (1990) 702–705. [2] A. Bossavit, A new rationale for edge-elements, Int. Compumag Soc. Newslett. 1 (1995) 3–6. [3] J.L. Coulomb, F.X. Zgainski, Y. Maréchal, A pyramidal element to link hexahedral, prismatic and tetrahedral edge finite elements, IEEE Trans. MAG-33 (1997) 1362–1365. [4] V. Gradinaru, R. Hiptmair, Whitney elements on pyramids, ETNA 8 (1999) 154–168. [5] J. Harrison, Lectures on Chainlet Geometry - New Topological Methods in Geometric Measure Theory, arXiv:math-ph/0505063, 153 pp., 19 May 2005. [6] A. Hatcher, Algebraic Topology, Cambridge U.P., Cambridge, 2002. [7] S.J. Owen, S. Saigal, Formation of pyramid elements for hexahedra to tetrahedra transitions, Comput. Methods Appl. Mech. Eng. 190 (2001) 4505–4518. [8] S. Van Welij, Calculation of Eddy currents in terms of H on hexahedra, IEEE Trans. MAG-21 (1985) 2239–2241. [9] H. Whitney, Geometric Integration Theory, Princeton U.P., Princeton, 1957.
Mathematics and Computers in Simulation 80 (2010) 1567–1577
A uniform rationale for Whitney forms on various supporting shapes A. Bossavit LGEP (CNRS), 11 Rue Joliot-Curie, 91192 Gif-sur-Yvette CEDEX, France Received 5 November 2008; accepted 6 November 2008 Available online 25 November 2008
Abstract A method for constructing Whitney forms is proposed, which applies to tetrahedra, hexahedra, triangular prisms, and pyramids in a similar way, and proceeds from a unique generating principle, thus unifying their presentation. The principle automatically enforces conformity (i.e., “tangential” or “normal continuity” of the elementary proxy fields) at element interfaces, and generates a complex of forms with the “exact sequence property” built in. © 2008 IMACS. Published by Elsevier B.V. All rights reserved. PACS: 65N30 Keywords: Finite elements; Edge elements; Whitney forms
1. Introduction In the finite element practice for electromagnetism, vector-valued “edge elements” on tetrahedra are now universally accepted, and the advantages of thinking of them as belonging to a larger family of scalar- or vector-valued simplexrelated fields – the Whitney forms – are well known. But flexibility in solid modelling requires such elements on other shapes of solids than simplices: prisms, hexahedra, pyramids. (It can easily be seen, for instance, that a transition layer of elements between a simplicial mesh and a brick mesh must include pyramids [7].) There is a rationale about “why Whitney forms are what they are” on simplices [2](which one may think Whitney himself knew about, although this is not explicitly stated in [9]). But this rationale does not apply directly to hexahedra and prisms, for which the known constructions [8] seem to be of a different flavor. Worse, for the pyramid, the constructions described in the few published papers (e.g., [3,4]) seem to stem from specific heuristic ideas, not obviously related with those used for tetrahedra. It is therefore desirable to derive all Whitney forms – edge-based, face-based, etc., for tetrahedra, prisms, etc. – from a unique generating principle. To this effect, I present a recursive procedure, consisting in the iteration of two basic geometric transforms, “conation” and “extrusion”, both able to “lift up” to dimension n + 1 a given n-dimensional Whitney complex. Hence (by the recursive formulas (8)–(12) below), Whitney forms on all kinds of n-topes that can be built from a point by iterated conation or extrusion. In dimension 3, this gives Whitney forms on the tetrahedron, hexahedron, triangular prism, and pyramid, allowing the practitioner to freely mix such shapes, each cell (node, edge, face,. . .) of this mixed mesh being assigned its own Whitney form, whose integral is 1 over this cell and 0 over all others. Their compatibility at interfaces, the property also known as “conformity” (or “tangential continuity” in the
E-mail address: [email protected]. 0378-4754/$36.00 © 2008 IMACS. Published by Elsevier B.V. All rights reserved. doi:10.1016/j.matcom.2008.11.005
1568
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
case of vector-field representation of 1-forms), is automatically enforced. Moreover, the “exact sequence property” of the whole complex of forms, so useful in the avoidance of spurious modes in resonant cavity problems [1], also comes naturally. 2. Basic notions Some notation (italics, when not used for emphasis, warn about implied definitions): we work in A∞ , the ∞dimensional real affine space, whose associated vector space is denoted V∞ , and we deal with subspaces and submanifolds (with boundary, most often) of finite dimension p ≥ 0. Symbol x and similar ones (y, z, etc.), refer to points, symbols v, w, etc., denote vectors. If v translates x to y, we may write y = x + v. Barycenters are denoted by tx + (1 − t)y, where x and y are points and t a real number, or by x + tv. We use the Grassmann algebra of multivectors (symbol ∨ for the product) and multicovectors (symbol ∧). Vectors (and by extension, all Grassmann objects) can be free (i.e., elements of V∞ ) or bound (i.e., anchored at some point x of A∞ ). Singular p-chains are formal weighted sums, with real coefficients, of p-submanifolds. Cochains, or differential forms, are the dual objects. We shall denote by v; ω the pairing between vector v and covector ω, andmore generally between p-vector and p-covector, or chain and cochain. In particular, when referring to the integral M ω of a form ω over a manifold M, we shall prefer the alternative notation M; ω. The boundary operator ∂ maps p-chains to (p − 1)-chains. The exterior derivative d is defined as its dual, acting on cochains: ∂M; ω = M; dω. Let be given a cell-complex Kn of finite dimension n, with compact carrier |Kn |, borne by an n-dimensional affine subspace. (Precise definitions will be found in Topology courses such as [6], but the notion is familiar: finite element meshes, provided they have no “hanging” nodes or edges, are cell complexes, whose carrier is the meshed region.) Since our purpose is to lift Kn , somehow, to dimension n + 1, we shall refer to it as the horizontal complex, and use underlined symbols for its elements. Similar symbols, not underlined, will serve for elements of the enlarged, vertical complex Kn+1 , of dimension n + 1, that will result from our construction. The complex Kn is made of separately and arbitrarily oriented p-cells, each the homeomorphic image of an open p-ball, if p > 0, and 0-cells are points (whose orientation is the number 1 or −1). The generic cell is denoted by s or σ, with the convention that dim(σ) = dim(s) − 1. Cells are related by incidence numbers, forming incidence matrices p , with 0 < p ≤ n: thus, sσ p = ±1 if the (p − 1)-cell σ belongs to ∂s, the sign depending on whether orientations of σ and s match or not. We denote by Sp the set of oriented p-cells, with p ≤ n. No generality will be lost by assuming there is only one n-cell, carried by |Kn |, as in Fig. 1. To Kn is naturally associated a chain complex {C, ∂}, by taking the restriction ∂ of ∂ to cell-based chains (i.e., linearcombinations of cells with real sweights). Note the relation between and ∂: the boundary of the p-chain c = s ∈ S cs s is σ ∈ S , s ∈ S sσ p c σ. So if c (boldface) stands for the array of weights {cs : s ∈ Sp }, one has p−1 p p ∂c = σ ∈ S (p c)σ σ. More concisely said, the family of maps κp defined, for each p, by c → s ∈ S cs s, is a p−1 p chain map: ∂κp = κp−1 p . We now describe the generation of Kn+1 by conation or extrusion of Kn . How these operations affect cochains, especially Whitney forms, will come later.
Fig. 1. Notations, and conventions about orientation. On the left, the simplex s is the face in the foreground.
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1569
First, conation. Let us select a point a (for “apex”) outside the subspace supporting Kn , and consider the cone |Kn+1 | = {λa + (1 − λ)x : x ∈ |Kn |, λ ∈ [0, 1]} that will carry the conation of the complex. More generally, if M is a submanifold, not necessarily in |Kn |, one sets cone(a, M) = {λa + (1 − λ)x : x ∈ M, λ ∈ [0, 1]}. (To avoid irrelevant difficulties, we assume that no half-line {a + tv : t ≥ 0}, where v is a nonzero vector, will intersect M more than once. The apex a will be understood and omitted from now on.) If M is oriented, this induces an orientation of cone(M), the one such that (Fig. 1) the following equality between singular chains holds: M = ∂[cone(M)] + cone(∂M).
(1)
(By convention, cone(∅) = {a}. When M reduces to a point x, with orientation +, this implies that cone(x) is the segment ax, oriented from a to x.) If x = / a, the point where ax intersects |Kn | is the projection πx of x. This induces a projection operation on other objets, for instance bound vectors anchored at point x, the projection of v being π∗ v. The complex Kn lifts up as follows: each of its p-cells σ has an associated so-called vertical cell s = cone(σ), carried sσ by the set {λa + (1 − λ)x : x ∈ σ, λ ∈ ]0, 1[}, and one may define new incidence numbers by sσ p = 1 and p = 0 for other cells σ in Sp . If cells σ and s of Kn are incident, their cones are incident too, and orientations match or not the same way, so one will have cone(s)cone(σ) = sσ p p . The horizontal cells of Kn+1 are just those of K n , with the sσ sσ same incidence numbers, p = p , so p is embedded in p as a submatrix. Last, a itself forms a new 0-cell, and sa 0 = −1 for all vertical 1-cells s. The set of p-cells in Kn+1 is denoted Sp . Extrusion, now. This time, a fixed vector X, not parallel to the subspace supporting |Kn |, is given. Setting ϕ(t, x) = x + tX, one defines the extrusion extr(|Kn |) as the set |Kn+1 | = {ϕ(λ, x) : x ∈ |Kn |, λ ∈ [0, 1]}. (More generally, one might have a smooth, singularity-free vector field X, with flow ϕ.) A point x ∈ |Kn+1 | is ϕ(λ, x) ¯ = x¯ = ϕ(1, x) besides its bottom projection for some λ ∈ [0, 1] and some x ∈ |Kn |, so it has now a top projection πx πx = x. A p-vector v at x projects, similarly, to π∗ v and π¯ ∗ v. We shall denote by rear(x) the segment x¯x and define rear(M) as the set-union {rear(x) : x ∈ M} in the case of a submanifold M in |Kn+1 | to which X is nowhere tangent. If ¯ by continuity, and rear(M) as well, in such a way M is oriented, this orients its bottom and top projections M and M, that ¯ M = ∂[rear(M)] + rear(∂M) + M.
(2)
To a cell σ in Kn correspond now two horizontal cells in Kn+1 , namely σ itself and its image ϕ(1, σ), which we denote ¯ (Thus horizontal cells in Kn+1 form two disjoint complexes, Kn on the one hand, and its translate by X, that one by σ. ¯ n , on the other hand.) Vertical cells are extr(σ) for σ ∈ Sp . Again, p is embedded in p , with minor may denote by K extr(σ)σ
variations: incidence numbers between horizontal cells are kept unchanged. We set p
s σ and p = sp σ¯
=0
for vertical cells s
extr(σ)σ¯
= 1 and p
= −1
other than extr(σ). Last, if cells σ and s of Kn are incident, their extrusions are extr(s)extr(σ)
sσ
incident too, and orientations match or not the same way, so one will have p = p . Let us now introduce J. Harrison’s chainlets, which are elements of the completion of a space of chains with respect to an appropriate norm, the natural norm [5]. Some passages to the limit, involving Cauchy sequences in the sense of this norm, can then profitably be associated with an element of the Banach space of chainlets. We must refer to [5] for details, and just mention results of this theory that will play a role here. For instance, a vector at point x can be considered as the chainlet limit of the sequence of chains cn = n [x, x + v/n] (i.e., the segment from point x to point x + v/n, with weight n), which allows us to define cone(v) as the chainlet limit of cone(cn ). Hence a hybrid object, neither bivector nor surface, but yet a bona fide chainlet, that we shall call a needle of dimension 2, or 2-needle. A similar treatment applies to a p-vector, seen as the chainlet limit of a sequence of smaller and smaller p-dimensional patches with higher and higher weights, the cone of which is a (p + 1)-needle. In the same spirit, a p-manifold M can be approached, in the sense of p-chainlets, by a Cauchy sequence of finite sums of (bound) p-vectors, each of which approximates a p-dimensional chunk of M. The corresponding free p-vectors can be added, and the limit of this sum is a (free) p-vector vect(M), that we shall call the vectorial p-volume of M. This operation extends to chains, and to chainlets. (One may see vect(c) as the map ω → c; ω, where ω is on the left a p-covector, on the right the constant p-form equal to it, since such a map can be identified, by duality, with a
1570
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
p-vector.) One checks that vect(v), if the bound vector v is considered as a chainlet, is the free vector v. Note also that vect(∂c) = 0, whatever c. We shall have use for the following relation, where v is a p-vector anchored at a point x ∈ |Kn+1 |, at height λ < 1 (Fig. 1), with horizontal projection πx relative to the apex a: vect(cone(v)) = (1 − λ)p+1 vect(cone(π∗ v)).
(3)
This is easily proven, for p = 1, by working on the n th term of a Cauchy sequence of 1-chains with chainlet limit v, for instance the segments [x, x + v/n] with weight n, as above. One then compares the areas of the two triangles cone([x, x + v/n]) and cone([x, x + π∗ v/n]), to see that their ratio tends to (1 − λ)2 . Generalizing to p > 1 is easy. In the case of extrusion, there is a similar relation: vect(rear(v)) = (1 − λ)vect(extr(π∗ v)), but nothing that simple will hold if the extrusion field X is not constant. 3. Heuristic moves Consider the known context of Whitney forms on a simplicial mesh: to each p-simplex s is associated a p-form well
ws such that s ws ≡ s; ws = δss (the Kronecker delta). These forms, traditionally, are understood as interpolants, like finite element shape functions. Starting from a smooth p-form a (the mathematical representation of some physical field, such as a vector potential if p = 1), one may build an array{as : s ∈ Sp } of numbers, one per p-cell, by setting Whitney forms can then be used to generate a p-form s as ws , which is an approximation of a in the sense as = s; a. that c; a − s as ws = 0 for, not all p-chains c of course, but for all simplicial p-chains. (To check this, set c = s and
use s; ws = δss .) The finer the mesh, the better this approximation, and one may base error analysis on that. We shall look at Whitney forms, here, equivalent, dual viewpoint: as approximants of chains, rather than of cochains. from an For this, note that c; s as ws = s ∈ Sp c; ws s; a. So if one sets whit(c) = s ∈ Sp c; ws , this simplicial chain is an approximation of c, in the sense that c − whit(c); a = 0 for, not all p-forms a, but all those in the space W p
generated by the ws . (To check this, set a = ws and use s; ws = δss .) Whitney forms are thus a device that, starting from a singular p-chain (or chainlet) c embedded in the carrier of a cell complex, delivers a cellular p-chain whit(c), meant to approximate c. Generating Whitney forms (should we not know them) thus amounts to build the map whit. As this section will show, this map is determined by the following requirements:
(1) That whit(s) = s (which is the same as s; ws = δss ). s (2) That it be a chain map: whit(∂c) = ∂(whit(c)) (which entails, by duality, dwσ = s ∈ Sp sσ p w , for σ in Sp−1 , the very property that helped eliminate spurious modes [1]). (3) That whit preserve vectorial volume: vect(whit(c)) = vect(c) (by duality, this means that W p contains all constant p-forms, also a useful feature for error analysis). (4) That it behave naturally under conation, as explained below. (Conditions (1) and (2) would not suffice, as the following counter-example shows: In dimension 1, with 0-cells {0} and {1} and 1-cell ]0, 1[, the 0-forms (1 − x) + αx(1 − x) and x − αx(1 − x) would make, together with the 1-form [1 − α(1 − 2x)] dx, a Whitney complex according to these criteria, whatever α, whereas α = 1/2 for the standard one.) Back to the situation of the previous section, let us assume that a complex of Whitney forms is already known on Kn , and that spaces W p they span have the exact sequence property: dW p ⊂ W p+1 and ker(d; W p ) = dW p−1 for p ≥ 1. So we know how to map a horizontal singular chain c (i.e., one embedded in |Kn |) to a horizontal cell-chain whit(c). Our purpose is to build whit when c is embedded in the larger region |Kn+1 |. This should yield Whitney spaces W p for Kn+1 , with a natural embedding W p ⊂ W p . In two particular cases, whit is forced on us. It will then be seen that no choice is left about its general definition. The first case is c = cone(σ) ≡ s, for a horizontal cell σ. Condition (1) then imposes whit(cone(σ)) ≡ whit(s) = s ≡ cone(σ) = cone(whit(σ)), so whit and cone commute when applied to horizontal (i.e., embedded in Kn ) cellular chains. This suggests to define whit(cone(c)) as cone(whit(c)) for all horizontal c. This is the “naturality” alluded to above, point (4). The second case is c = cone(v), a (p + 1)-needle based on the p-vector v anchored at a point x ∈ |Kn+1 |, at height λ < 1, with horizontal projection πx. The projection π∗ v of v, a p-vector anchored at πx, is a chainlet with respect to the horizontal complex Kn . So whit(π∗ v) is a definite horizontal p-chain, namely s ∈ S π∗ v; ws (πx) s, where ws p
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1571
denotes the form in W p associated with s. According to the third condition above, and taking relation (3) into account, we require that whit(cone(v)) = (1 − λ)p+1 cone(whit(π∗ v))
(4)
in order to preserve vectorial volume. Let us unfold the implications. For a p-vector v, (4) entails whit(cone(v)) = (1 − λ)p+1 π∗ v; ws (πx) cone(s) s ∈ Sp
=
v; (1 − λ)p+1 (π∗ ws )(πx) cone(s)
(5)
s ∈ Sp
where π∗ denotes the pullback operation. (By definition, π∗ v; ω = v; π∗ ω for all v, ω.) Here, λ had its value at the anchor point of v, and the bracket paired a p-vector with a p-covector. But (1 − λ)p+1 π∗ ws makes sense as a p-form living on |Kn+1 |, with now λ a function of x. So we can replace v by a p-chain c and write whit(cone(c)) = c; (1 − λ)p+1 π∗ ws cone(s), s ∈ Sp
with the bracket now an integral over c (that is to say, the limit of a Riemann sum of terms similar to the bracket in (5)). The same formula applied to ∂c yields (the dummy variable changes from s to σ, according to our general convention) whit(cone(∂c))
=
∂c; (1 − λ)p π∗ wσ cone(σ)
σ ∈ Sp−1
=
c; d[(1 − λ)p π∗ wσ ]cone(σ)
σ ∈ Sp−1
so finally, using (1), whit(c) = whit(∂[cone(c)]) + whit(cone(∂c)) = ∂[whit(cone(c))] + whit(cone(∂c)) = c; (1 − λ)p+1 π∗ ws ∂(cone(s)) + c; d[(1 − λ)p π∗ wσ ]cone(σ). s ∈ Sp
σ ∈ Sp−1
(6)
On the other hand, whit(c) = s ∈ Sp c; ws s, where the p-cells s belong to Kn+1 , so we shall find what ws must be by the tedious process of comparing this with (6). Cells s in Sp are of two kinds: Horizontal ones, which already belonged to Sp , and vertical ones, of the form s = cone(σ), where σ spans Sp−1 , so we deal with them separately. (Cf. Fig. 2 for mnemonic aid.) For a horizontal cell s, the only contribution of the right-hand side of (6) to c; ws is c; (1 − λ)p+1 π∗ ws , since ∂(cone(s)) = s + (. . .), where the dots stand for a chain of other, vertical, cells. So we find ws = (1 − λ)p+1 π∗ ws for horizontal p-cells of Kn+1 . (Note that the“horizontal” Whitney form ws lives on |Kn+1 |, and that ws is its restriction
Fig. 2. Current notation.
1572
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
to |Kn |.) This helps simplify (6) a bit: c; ws ∂(cone(s)) + c; dwσ cone(σ). whit(c) = s ∈ Sp
(6 )
σ ∈ Sp−1
We also observe that, for a horizontal (p − 1)-cell σ, dwσ
= =
d[(1 − λ)p π∗ wσ ] = −p(1 − λ)p−1 dλ ∧ π∗ wσ + (1 − λ)p π∗ dwσ p ∗ s −p(1 − λ)p−1 dλ ∧ π∗ wσ + sσ p (1 − λ) π w , s ∈ Sp
since the exact sequence property for Kn implies dwσ = we have s (1 − λ)dwσ = −p dλ ∧ wσ + sσ p w ,
sσ s s ∈ Sp p w .
So, for horizontal (p − 1)-forms on Kn+1 , (7)
s ∈ Sp
a formula which will prove useful immediately.
For a vertical cell s = cone(σ), factors of cone(s) in (6 ), in addition to the term c; dwσ that comes from the second
sum on the right of (6 ), will also come from the first sum, on behalf of cells s ∈ Sp incident on σ. Specifically, since
∂s = σ ∈ S sp σ σ , one has p−1
∂(cone(s )) = s − cone(∂s ) = s − sp σ cone(σ ), σ ∈ Sp−1
from which one gets, using (6 ) then (7),
s σ σ σ c; ws = c; dwσ − σs p c; w = c; dw − (1 − λ)dw − p dλ ∧ w . s ∈ Sp
So, finally, ws = λ dwσ − p dλ ∧ wσ for a vertical p-cell s = cone(σ). We could now repeat all these moves in the case of extrusion. This is straightforward. One will notice however, in so doing, that the extrusion vector field must be constant for vectorial volume to be preserved. 4. Whitney forms on simplices and bricks Our heuristic explorations have thus led to generating formulas. In the case of conation, (a) For 0 ≤ p ≤ n, and s ∈ Sp , the p-form ws , living on Kn , is lifted to Kn+1 by the formula ws = (1 − λ)p+1 π∗ ws ,
(8)
hence the horizontal p-form ws . (b) For 1 ≤ p ≤ n + 1, and σ ∈ Sp−1 , the vertical p-form associated with s = cone(σ) is ws = λ dwσ − p dλ ∧ wσ .
(9)
(c) The apical 0-form, or height of point x, denoted λ(x), in W 0 , completes the inventory. In the case of extrusion, λ being the function such that x = ϕ(λ(x), πx), (a) For 0 ≤ p ≤ n, and s ∈ Sp , each p-form ws , living on Kn , spawns two horizontal p-forms, the bottom one: ws = (1 − λ) π∗ ws ,
(10)
and the top one: wϕ(1,s) = λ π∗ ws .
(11)
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1573
(b) For 1 ≤ p ≤ n + 1, and σ ∈ Sp−1 , the vertical p-form associated with s = extr(σ) is ws = −dλ ∧ π∗ wσ ≡ −dλ ∧ (wσ + wϕ(1,σ) ),
(12)
as seen by adding (10) and (11). Iterated conation under Rules (a)(b)(c) generates the classical Whitney forms for any degree and spatial dimension, starting from the virtual nothing that K0 is. (When n = 0, the complex K0 reduces to a single 0-cell σ and its 0-form wσ = 1.) Let us formally prove this. First, any p-form, 0 ≤ p ≤ n, in the complex Kn+1 over the (n + 1)-simplex |Kn+1 |, does come, by rule (a), from a p-form of the complex Kn over some n-dimensional face of |Kn+1 |. To prove this, start from the known generic expression of Whitney forms [9], wς = p!
p
(−1)i λς(i) dλς(0) ∧ dλς(1) ∧ · · ·(i)· · · ∧ dλς(p) ,
(13)
i=0
where ς is any injective map from the integer segment [0, p] into the set of nodes of Kn+1 , and (i) means “omit the term indexed by i”. (Thus, ς labels an oriented p-face, and there is a Whitney form for each of these.) Suppose p ≤ n. There is then at least one node distinct from any of the ς(i)’s, that we take as cone’s tip and call a, denoting by λ its nodal function and by π the projection from a to the opposite n-face, the n-simplex |Kn |. Call λς(i) the restriction of λς(i) to the latter. Then λς(i) = (1 − λ)π∗ λς(i) . Since d commutes with pullback, d(π∗ λς(i) ) = π∗ (dλς(i) ). Then dλς(i) = (1 − λ)π∗ (dλς(i) ) − π∗ (λς(i) dλ). Substituting in (13), one gets two parts, one with dλ factored out where terms cancel two by two, and another part that is p! (1 − λ)p+1
p
π∗ ((−1)i λς(i) dλς(0) ∧ · · ·(i)· · · ∧ dλς(p) ).
i=0
(1 − λ)p+1 π∗ (wς ),
= confirming (8). It remains to be proven that the top form (p = n + 1) in Kn+1 comes from the top form (p = n) of some Kn by rule (b). So here, s is the (n + 1)-simplex |Kn+1 | and σ is one of its n-faces. Let’s label the nodes of σ by 0, 1, . . . , n and the corresponding barycentrics by λi , while λ is the barycentric of the node of s thus left out. The restriction of λi to σ is λi , and λi = (1 − λ)π∗ λi . The top form wσ in Kn is (up to sign) n! dλ1 ∧ · · · ∧ dλn and lifts to wσ = (1 − λ)n+1 π∗ wσ . Since the degree of wσ is n, one has dwσ = 0, hence d(π∗ wσ ) ≡ π∗ (dwσ ) = 0, so dwσ = −(n + 1)(1 − λ)n dλ ∧ π∗ wσ . Using (9), with p = n + 1, we obtain So
wς
ws
= −(n + 1)(1 − λ)n [λ dλ ∧ π∗ wσ + (1 − λ) dλ ∧ π∗ wσ ] = −(n + 1)(1 − λ)n dλ ∧ π∗ wσ
−(n + 1)! (1 − λ)n dλ ∧ π∗ (dλ1 ∧ · · · ∧ dλn ) = −(n + 1)! dλ ∧ i=1,...,n (1 − λ) d(π∗ λi ), i ∗ i where the expression headed by the symbol is the wedge product ofi its terms. Since dλ = d[(1 − λ) π λ ] = ∗ i ∗ i s (1 − λ)π (dλ ) − π λ dλ, we see that w = −(n + 1)! dλ ∧ i=1,...,n dλ , which is indeed the top form in Kn+1 , up to sign. So iterated conation, with rules (8) and (9), does generate the simplicial Whitney forms. We may expect iterated extrusion to generate Whitney forms on a hypercube. Verifying this fact is a lot easier than what precedes, as we presently see. Let v1 , . . . , vn be independent vectors (in the infinite-dimensional real vector space V∞ ), and consider the parallelotope built (in the associated affine space A ∞ ) on v1 , . . . , vn , anchored at some origin point o. Call this |Qn |. A natural representation of its points is x = o + i=1,...,n xi vi , where 0 ≤ xi ≤ 1 ∀i. (Here one ought to distinguish x, a point in A∞ , from the real n-tuple x = {xi : i = 1, . . . , n} of its coordinates. Denoting both entities by the same symbol x should not be too confusing in what follows.) The Whitney complex Qn over |Qn | is known [8] to consist of the following elements: =
1574
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
Fig. 3. Labelling scheme and orientations, for n = 3. Under each label, the Whitney form associated with this p-face. For better readability, coordinates x1 , x2 , x3 are here x, y, z. The top form is dx ∧ dy ∧ dz, associated with the cube, the label of which is ••• . 1−k
k
• The 0-forms wk (x) = i=1,...,n (xi ) i (1 − xi ) i ≡ pk (x) where ki = 0 or 1. Here, k is the integer i=1,...,n 2i−1 ki , whose binary representation is obtained by concatenating k1 , k2 , . . . , kn , which makes a convenient label for the 2n nodes of the n-cube and for the 0-forms wk . k 1−k • The 1-forms wk,j (x) = i =/ j (xi ) i (1 − xi ) i dxj , associated with edges parallel to vj . ki 1−ki k,j i i • The 2-forms w (x) = i∈J dxj1 ∧ dxj2 / (x ) (1 − x ) where J is the subset {j1 , j2 } of [1, n], with j1 < j2 , which determines a set of faces, those parallel to the 2-vector vj1 ∨ vj2 . • And so forth, according to a pattern that Fig. 3 should help understand: let us denote by σ an increasing injection from [1, p] into [1, n], if p ≥ 1, and by Jσ its image. (If p = 0, set Jσ = ∅.) This way, each σ labels the set of p-faces parallel to the p-vector vσ(1) ∨ · · · ∨ vσ(p) . Now, each of these 2n−p faces is labelled by an integer k, as above (with some harmless redundancy, since k and k point to the same p-face when kj = kj ∀j ∈ Jσ ). Thus the pair {k, σ} labels a definite p-face. We shall set dxσ = dxσ(1) ∧ dxσ(2) ∧ · · · ∧ dxσ(p) , that is to say, dxσ = dxj1 ∧ dxj2 ∧ · · · ∧ dxjp where j1 < j2 < · · · < jp compose the set Jσ , and ki
1−ki
i i pk,σ (x) = i∈J / σ (x ) (1 − x )
.
(14)
This polynomial, of order n − p, is linear with respect to each component xi taken separately (“n-linear” polynomial). Note that it coincides with the above-mentioned pk , i.e., the nodal function, when Jσ is empty. With such notation, and the convention that dxσ is the scalar constant 1 when Jσ is empty, all p-forms of the complex, 0 ≤ p ≤ n, can be written wk,σ (x) = pk,σ (x) dxσ .
(15)
Fig. 3 illustrates all that in dimension 3. Note the use of dots to place the set Jσ in [1, n]. The exact sequence property is now easy to establish: Applying d to (15) gives ∂i pk,σ dxi ∧ dxσ , dwk,σ = i∈J / σ
where ∂i is differentiation with respect to xi . After (14), ∂i pk,σ does depend on xi any longer, and dxi ∧ dxσ = ±dxs , not k,s k,σ where s maps onto the set Js = Jσ ∪ {i}. Therefore dw = s ± p dxs , where s spans the set of injections from / Jσ . [1, p + 1] to [1, n + 1] whose images are the sets Js . There is one such s for each i ∈ Now, let us add another independent vector vn+1 , call xn+1 what was λ in formulas (10)–(12), and extrude the k 1−k previous complex via vn+1 . A factor (xn+1 ) n (1 − xn+1 ) n appears when using (11)–(13) throws in a new term dxn+1 . The pattern (14) and (15) is thus respected, so extrusion of Qn does yield the Whitney complex Qn+1 on the (n + 1)-cube |Qn+1 |.
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1575
Parallelotopes thus obtained are not convenient, however, for finite element practice, which may require distorted (but still, stackable) hexahedra. For this, one may use the well known “isoparametric” trick, which consists in mapping the cube vertex labelled by k to some point uk in A∞ , with interpolation by the nodal functions. Point x then goes to u(x) = pk (x)uk . (16) k
Under some mild restrictions about how to place the points uk , this map will be injective, with a smooth inverse u−1 , that makes a convenient chart from u(|Qn |) to Rn . To get the Whitney complex on u(|Qn |), we just collect the forms ∗ (u−1 ) wk,σ obtained by pullback, using this chart. There is one map v u like u in (16) for each elementary volume, with −1
shared subsets of vertices uk . Maps (v u)−1 and (v u) thus coincide on p-faces common to volumes v and v in the actual mesh, and so do the restrictions of their pullbacks, which guarantees stackability and conformity. A well-known drawback of this isoparametric approach is that constant p-forms do not belong to the Whitney space thus obtained, for p ≥ 2 : this results from losing vectorial volume conservation. 5. Exact sequence property Let us now verify from scratch the fact that both processes, whatever the kind of complex they act on, provided it has the exactness property, generate a complex that does inherit this property. All we have to show in fact (thanks to the 1–1 correspondence between p-forms thus built and p-cells) is that each dwσ , where σ is a (p − 1)-cell, is a linear combination of p-forms ws of the conation or extrusion, as the case may be. First, conation. Call W p and W p the spaces generated by the p-forms of the vertical and the horizontal complexes. Assuming that dW p−1 ⊂ W p , we need to show that dwσ ∈ W p when wσ belongs to W p−1 . There are two cases, depending on whether σ is a horizontal or vertical cell. (Cf. Fig. 2 again.) For a horizontal cell σ ∈ Sp−1 , use (8), i.e., wσ = (1 − λ)p π∗ wσ . Then dwσ = −p (1 − λ)p−1 dλ ∧ π∗ wσ + (1 − λ)p π∗ (dwσ ), so (1 − λ) dwσ = −p dλ ∧ wσ + (1 − λ)p+1 π∗ (dwσ ), and dwσ = λ dwσ − p dλ ∧ wσ + (1 − λ)p+1 π∗ (dwσ ). After (9), λ dwσ − p dλ ∧ wσ is ws for the cell s = cone(σ), so we have dwσ = ws + (1 − λ)p+1 π∗ (dwσ ).
By hypothesis, dwσ = s ∈ S sp σ ws . Therefore, carrying on, p
dwσ = ws +
s ∈ Sp
sp σ (1 − λ)p+1 π∗ ws = ws +
s ∈ Sp
sp σ ws
by using (8) again for p-cells s . So dwσ ∈ W p . For a vertical cell s, we have s = cone(σ) for some σ ∈ Sp−1 . First, a useful lemma: Lemma. For s vertical, i.e., s = cone(σ), one has dws = (p + 1)dλ ∧ dwσ ,
λ dws − (p + 1)dλ ∧ ws = 0.
Proof. By the recursive formula, ws = λ dwσ − p dλ ∧ wσ , so dws = dλ ∧ dwσ + p dλ ∧ dwσ = (p + 1)dλ ∧ dwσ . Next, dλ ∧ ws = dλ ∧ [λ dwσ − p dλ ∧ wσ ] = λ dλ ∧ dwσ , and hence, λ dws − (p + 1) dλ ∧ ws = 0.
1576
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
Fig. 4. Labelling scheme and orientations for the pyramid, as seen from above (•• 0 is the bottom square). Under or below each label, the Whitney form associated with this p-face. On the right, the base square’s 2-form (above •• 0) and the pyramid’s 3-form (above ••• ).
Now, since all (p + 1)-faces on which s is incident are vertical, we may suspect that dws =
s σ cone(s ) , s ∈ Sp p w
give or take a sign, so let’s evaluate this sum, using (9) to express wcone(s ) :
sp σ wcone(s ) = sp σ [λ dws − (p + 1) dλ ∧ ws ] s ∈ Sp
s ∈ Sp
⎛
=
λd⎝
s ∈ Sp
⎞
sp σ ws ⎠ − (p + 1) dλ ∧
s ∈ Sp
sp σ ws
=
λ d(dwσ − ws ) − (p + 1) dλ ∧ (dwσ − ws ) = −(p + 1) dλ ∧ dwσ − [λ dws − (p + 1) dλ ∧ ws ]
which is −dws indeed, thanks to the Lemma. So dws ∈ W p , and this ends the proof of exactness. σ Extrusion, now: start from wσ = (1 − λ)π∗ wσ for a (p − 1)-cell σ (cf. (10)). Then dw ∧ π∗ wσ + (1 − = −dλ
σ s
∗ σ s ∗ σ σ s λ) d(π w ), which is w + (1 − λ)π dw after (12), with s = extr(σ). Again, dw = s ∈ S p w , and ws = p
(1 − λ)π∗ ws by using (10), so we do have dwσ = ws + s ∈ S sp σ ws . p
6. Whitney forms for other shapes: pyramid, prism On a point o and two independent vectors vx and vy , let us build, by two successive extrusions, the complex K2 made of • The 0-forms (1 − x)(1 − y), (1 − x)y, x(1 − y), xy, on the 0-cells 00, 01, 10, 11, respectively; • The 1-forms (1 − y) dx, y dx, (1 − x) dy, x dy, on the 1-cells • 0, • 1, 0• , 1• , respectively; • The 2-form dx ∧ dy on the horizontal parallelogram, the 2-cell •• in this notational scheme. They generate a Whitney complex {W 0 , W 1 , W 2 }, as we have verified earlier. Now pick an out-of-plane point a, and make it the apex of a pyramid based on •• . For the sake of notational uniformity, let us denote by z what was formerly λ. Points p of this pyramid are parameterized as follows: p(x, y, z) = za + (1 − z)(x vx + y vy ), 0 ≤ x ≤ 1, 0 ≤ y ≤ 1, 0 ≤ z ≤ 1.
(17)
A. Bossavit / Mathematics and Computers in Simulation 80 (2010) 1567–1577
1577
Of course, these are not Cartesian coordinates (Fig. 4). The projection π is given by π(p(x, y, z)) = x vx + y vy : note the invariance of coordinates x and y by projection, which makes this chart a convenient one for pullback. As an alternative to (17), we have p(x, y, z) = (1 − z){x[(1 − y)10 + y 11] + (1 − x)[(1 − y)00 + y 01]} + za, which better shows the weights of p with respect to each of the five nodes. The conified complex K3 has horizontal 0-cells 000, 010, 100, 110, 1-cells • 00, • 10, 0• 0, 1• 0, and 2-cell •• 0, according to an obvious nomenclature (Fig. 4), whose conations form the vertical 1-cells 00• , 01• , 10• , 11• , 2-cells • 0• , • 1• , 0•• , 1•• , and 3-cell ••• (the pyramid itself), to which one must add the apical 0-cell a. Next, we apply the recipe (8), with z substituted for λ, then find the vertical forms using (9). This tedious process, details of which are left out, results in Fig. 4. The top form (p = 3) comes from σ=•• 0, with wσ (p(x, y, z)) = (1 − z)3 dx ∧ dy. Then s ≡ cone(σ)is ••• , the pyramid itself, and ws (p(x, y, z)) =
z d[(1 − z)3 dx ∧ dy] − 3(1 − z)3 dz ∧ dx ∧ dy
=
−3z(1 − z)2 dz ∧ dx ∧ dy − 3(1 − z)3 dz ∧ dx ∧ dy
=
−3(1 − z)2 dx ∧ dy ∧ dz
which is the expected volume form for these coordinates. (The “wrong” sign is due to our general orientation rules.) The apical 0-form, wa (p(x, y, z)) = z, completes the list (Fig. 4). A word, finally, on the triangular prism. The 3D simplex was generated from naught by conation (or extrusion: the same thing in dimension 1), conation, conation, the brick by conation, extrusion, extrusion, the pyramid by conation, extrusion, conation. The prism stems from conation, conation, extrusion. This is left as an easy exercise. For both the pyramid and the prism, warped shapes are allowed if one proceeds by mapping and pullback, as in the case of bricks. The same caveat applies about the representation of constant fields. Let us finally remark that, whatever the volume, traces on faces only depend on the face itself (be it a triangle or a quadrilateral). Thus conformity is built in. References [1] A. Bossavit, Solving Maxwell’s Equations in a Closed Cavity, and the Question of Spurious Modes, IEEE Trans. MAG-26 (1990) 702–705. [2] A. Bossavit, A new rationale for edge-elements, Int. Compumag Soc. Newslett. 1 (1995) 3–6. [3] J.L. Coulomb, F.X. Zgainski, Y. Maréchal, A pyramidal element to link hexahedral, prismatic and tetrahedral edge finite elements, IEEE Trans. MAG-33 (1997) 1362–1365. [4] V. Gradinaru, R. Hiptmair, Whitney elements on pyramids, ETNA 8 (1999) 154–168. [5] J. Harrison, Lectures on Chainlet Geometry - New Topological Methods in Geometric Measure Theory, arXiv:math-ph/0505063, 153 pp., 19 May 2005. [6] A. Hatcher, Algebraic Topology, Cambridge U.P., Cambridge, 2002. [7] S.J. Owen, S. Saigal, Formation of pyramid elements for hexahedra to tetrahedra transitions, Comput. Methods Appl. Mech. Eng. 190 (2001) 4505–4518. [8] S. Van Welij, Calculation of Eddy currents in terms of H on hexahedra, IEEE Trans. MAG-21 (1985) 2239–2241. [9] H. Whitney, Geometric Integration Theory, Princeton U.P., Princeton, 1957.