Conventions

The finite element spaces H1, HCurl, HDiv and L2 of NGSolve form exact sequences \[H^1 \xrightarrow{\ \nabla\ } H(\operatorname{curl}) \xrightarrow{\ \operatorname{curl}\ } H(\operatorname{div}) \xrightarrow{\ \operatorname{div}\ } L^2\] on triangles, quadrilaterals, tetrahedra, prisms and hexahedra, and on meshes that mix these shapes. One order convention holds for all shapes: the order \(k\) is the polynomial degree of the trace, \(\mathcal P_k\) on edges and triangular faces, \(\mathcal Q_k\) on quadrilateral faces. Since the traces do not depend on the element shape or on the flags below, all variants are conforming to each other.

• H1(k): \(\mathcal P_k\) on simplices, \(\mathcal Q_k\) on quads and hexes, \(\mathcal P_k(x,y)\otimes\mathcal P_k(z)\) on prisms. L2(k) uses the same polynomials without continuity.

• HCurl(k) \(= \nabla\,\)H1(k+1) \(\oplus\, Z_k\). The tangential trace has degree \(k\). The non-gradient part \(Z_k\) has injective curl, and \(\operatorname{curl}Z_k\) lies in HDiv(k-1). With type1=True only the gradients of H1(k) are used; \(Z_k\) stays the same. Thus type1 removes the top-degree gradients, and the edge traces drop to degree \(k-1\).

• HDiv(k): the normal trace has degree \(k\). The space consists of facet functions, divergence-free cell bubbles \(B_k\), and a cell part \(D_q\) with \(\operatorname{div}D_q =\) L2(q)\(/\mathbb R\). With RT=True, \(q=k\) and \(\operatorname{div}\) maps onto L2(k); with RT=False, \(q=k-1\). Thus RT adds the top-degree divergences. If RT is not given, quads, prisms and hexes use RT=True, and trigs, tets and pyramids use RT=False.

Below, HCurl(k,type1), HDiv(k,RT) and HDiv(k,BDM) are short for type1=True, RT=True and RT=False. The spaces are nested, \[\texttt{HCurl(k,type1)} \subset \texttt{HCurl(k)} \subset \texttt{HCurl(k+1,type1)},\qquad \texttt{HDiv(k,BDM)} \subset \texttt{HDiv(k,RT)} \subset \texttt{HDiv(k+1,BDM)},\] and they contain the full polynomials \(\mathcal P_k^d \subset\) HCurl(k), \(\mathcal P_{k-1}^d\subset\) HCurl(k,type1), \(\mathcal P_k^d\subset\) HDiv(k,BDM). On simplices these are the classical elements: HCurl(k,type1) and HCurl(k) are the Nédélec elements of the first [4] and second kind [5], HDiv(k,RT) and HDiv(k,BDM) the Raviart–Thomas and Brezzi–Douglas–Marini elements.

The basis functions are hierarchical and follow the construction of [1, 2]: each space consists of lowest-order functions, of the images of the high-order functions of the previous space, and of functions that complete the space. This splitting holds separately on every edge, face and cell (local complete sequence property). It is the reason why gradients or divergences can be switched on and off by flags, and it is shown in the tables of Section 3.

Exact sequences

For every choice of the two flags the spaces fit together. In three dimensions the sequences \[\begin{aligned} \text{(A)}\quad & \texttt{H1(k)} \xrightarrow{\nabla} \texttt{HCurl(k,type1)} \xrightarrow{\operatorname{curl}} \texttt{HDiv(k-1,RT)} \xrightarrow{\operatorname{div}} \texttt{L2(k-1)}\\ \text{(B)}\quad & \texttt{H1(k)} \xrightarrow{\nabla} \texttt{HCurl(k,type1)} \xrightarrow{\operatorname{curl}} \texttt{HDiv(k-1,BDM)} \xrightarrow{\operatorname{div}} \texttt{L2(k-2)}\\ \text{(C)}\quad & \texttt{H1(k+1)} \xrightarrow{\nabla} \texttt{HCurl(k)} \xrightarrow{\operatorname{curl}} \texttt{HDiv(k-1,RT)} \xrightarrow{\operatorname{div}} \texttt{L2(k-1)}\\ \text{(D)}\quad & \texttt{H1(k+1)} \xrightarrow{\nabla} \texttt{HCurl(k)} \xrightarrow{\operatorname{curl}} \texttt{HDiv(k-1,BDM)} \xrightarrow{\operatorname{div}} \texttt{L2(k-2)} \end{aligned}\] are exact on contractible domains: the range of each operator is the kernel of the next one, and \(\operatorname{div}\) is onto. On tetrahedra these are the four sequences of Arnold, Falk and Winther [3]; (A) is the sequence of the first kind, (D) of the second kind. The same four sequences exist on prisms, hexes and mixed meshes. In two dimensions the sequences are H1(k) \(\to\) HCurl(k,type1) \(\to\) L2(k-1) and H1(k+1) \(\to\) HCurl(k) \(\to\) L2(k-1). The 2D HDiv spaces are the rotated HCurl spaces, HDiv(k,BDM) \(=\) rot HCurl(k) and HDiv(k,RT) \(=\) rot HCurl(k+1,type1).

If RT is not given, the divergence maps onto L2(k-1) on quads, prisms and hexes, and onto L2(k-2) on trigs and tets. The space V.CreateDivergenceSpace() is the L2 space with these element-wise orders, so that \[\texttt{H1(k+1)} \xrightarrow{\nabla} \texttt{HCurl(k)} \xrightarrow{\operatorname{curl}} \texttt{V = HDiv(k-1)} \xrightarrow{\operatorname{div}} \texttt{V.CreateDivergenceSpace()}\] is exact on mixed meshes as well, also with HCurl(k,type1) and H1(k) in front.

Degrees of freedom

The tables show the element spaces of the sequences that start from H1(2), with the number of dofs per vertex (V), edge (E), triangular face (T), quadrilateral face (Q) and cell (C), and the element dimension. Black is the lowest-order dof (Whitney edge function, RT\(_0\) face function), green dofs span the kernel of the outgoing operator (gradients in HCurl, divergence-free functions in HDiv), blue dofs are mapped injectively into the next space. The green dofs of one space are the images of the blue dofs of the previous one, node by node.

On the quad, HCurl(k,type1) is the tensor-product Nédélec space \(N_{[k-1]}=\mathcal Q_{k-1,k}\times\mathcal Q_{k,k-1}\). The default space HCurl(k) \(=\nabla\mathcal Q_{k+1}+N_{[k-1]}\) has dimension \((k{+}2)^2-1+k^2\); its curl is \(\mathcal Q_{k-1}\), as \(\mathcal P_{k-1}\) on the triangle. The cell dofs of order \(k\) are

The alternating sum of the element dimensions is 1 for each sequence, for example (A) on the prism: \(18-36+25-6=1\), and (C) on the hex: \(27-31+6-1=1\). The face dofs of the 3D elements are those of the 2D elements: a quadrilateral face of HCurl(2,type1) carries 1+3 dofs like the cell of the quad, a triangular face 2. This makes meshes of tets, prisms and hexes conforming. The cell dofs show the bubble spaces: the blue cell dofs of HCurl and the green cell dofs of HDiv both count the divergence-free bubbles \(B_1\) (tet 0, prism 2, hex 5), and the blue cell dofs of HDiv count \(\dim\)L2\(\,-\,1\).

Element spaces.

The following table lists the spaces for general order. HCurl(k) is HCurl(k,type1) plus the gradients of H1(k+1) on every shape.

Here \(\tilde{\mathcal P}_k\) are the homogeneous polynomials, and \(N^\triangle_k\), RT\(^\triangle_k\) are the triangle spaces HCurl(k,type1) and HDiv(k,RT). HDiv(k,BDM) is the subspace of HDiv(k,RT) with divergence in L2(k-1); on tets it is \(\mathcal P_k^3\). Element dimensions for \(k=1,2,3\):

Flags and helper functions

The flags change only cell dofs of HDiv and only gradient dofs of HCurl, so every combination is conforming. HDivSurface takes the same RT flags.

V.CreateDivergenceSpace() returns the L2 space onto which \(\operatorname{div}\) maps V: order \(k\) on elements with RT, order \(k-1\) on the others, and order \(0\) with hodivfree. On a mixed mesh with default flags neither L2(mesh, order=k) nor L2(mesh, order=k-1) is the right pressure space for HDiv(mesh, order=k); the first is too large on simplices, the second too small on tensor-product elements. G, W0 = W1.CreateGradient() returns the matching space W0, which is H1(k+1) or with type1 H1(k), and the sparse matrix G of the embedding \(\nabla:\) W0 \(\to\) W1.

Example.

The script builds the default sequence on an hp-refined mesh of the Fichera corner, which consists of 501 tets, 209 prisms and 37 hexes.

from ngsolve import *
from netgen.occ import *

# hp-refined Fichera corner: tets, prisms and hexes
cube = Box((0,0,0), (2,2,2))
cube.edges.hpref = 1
geo = OCCGeometry(Box((-1,-1,-1), (1,1,1)) - cube)
mesh = Mesh(geo.GenerateMesh(maxh=0.4))
mesh.RefineHP(levels=2, factor=0.2)

k = 2
W0 = H1(mesh, order=k+1)
W1 = HCurl(mesh, order=k)            # or HCurl(mesh, order=k+1, type1=True)
W2 = HDiv(mesh, order=k-1)           # RT on prisms/hexes, BDM on tets
W3 = W2.CreateDivergenceSpace()      # L2(k-1) there, L2(k-2) on tets
nd = [W.FreeDofs().NumSet() for W in (W0, W1, W2, W3)]   # used dofs
print(nd, nd[0] - nd[1] + nd[2] - nd[3])

# curl of an HCurl function is represented exactly in HDiv
u = GridFunction(W1)
u.Set(CF((sin(y*z), x*z, exp(x)*y)))
w = GridFunction(W2)
w.Set(curl(u))
print(sqrt(Integrate((w - curl(u))**2, mesh)))
print(sqrt(Integrate(div(w)**2, mesh)))

It prints the dimensions [7363, 13261, 7950, 2051] with alternating sum 1, the Euler characteristic of the domain, and the errors \(7\cdot10^{-15}\) and \(2\cdot10^{-13}\). The count uses FreeDofs instead of ndof, since RefineHP leaves unused edges and faces in the mesh, and their dofs are contained in ndof.

Verification and limits

Tests.

Exactness was tested algebraically, by comparing the ranks of the operators with the space dimensions and by checking that the \(L^2\)-projection of \(\nabla u\), \(\operatorname{curl}u\), \(\operatorname{div}u\) onto the next space reproduces it.

• Mesh of one hex, two prisms and two tets, with two vertex numberings: sequences (A)–(D) and the default sequences with CreateDivergenceSpace are exact for \(k\le 3\).

• Mesh of 12 trigs and 10 quads: both HCurl sequences, and H1(k+1) \(\to\) HDiv(k) \(\to\) divergence space for RT set, unset and False, are exact for \(k\le4\).

• The inclusions of the full polynomials and the nesting of the spaces of Section 1 hold on the mixed 3D mesh for \(k\le 3\) and \(k \le 2\), respectively. The dof counts of Section 3 are those reported by the spaces.

Changes against version 6.2.2604.

Triangles and tets are unchanged. On quads, prisms and hexes the default HCurl(k) is smaller: the curl now has degree \(k-1\) as on simplices (quad HCurl(1): 12 \(\to\) 9 dofs, prism 30 \(\to\) 21, hex 54 \(\to\) 31). On quads and hexes the former space is available as HCurl(k+1,type1). type1 on prisms and hexes was a mixture of both kinds whose curl kernel was larger than \(\nabla\,\)H1(k); HCurl(2,type1) went from 67 to 36 dofs on the prism and from 132 to 54 on the hex. The HDiv prism is now the tensor-product Raviart–Thomas space (HDiv(1): 27 \(\to\) 25 dofs). RT=False was ignored on quads, prisms and hexes and now selects the smaller space.

Limits.

• Pyramids are conforming to all neighbours, and \(\ker\operatorname{curl}=\nabla H^1\) holds. For HCurl of order \(k\ge2\) the curl is not contained in HDiv(k-1) on pyramids, so sequences on meshes with pyramids are exact only at lowest order.

• On non-affine quads and hexes the stages \(\nabla\) and, in 3D, \(\operatorname{curl}\) stay exact. The last operator maps into \(J^{-1}\,\)L2 instead of L2, where \(J\) is the Jacobian determinant; the ranks of the operators are unchanged.

• Anisotropic polynomial orders on prisms and hexes are not tested with the new spaces.

References

1. S. Zaglmayr. High Order Finite Element Methods for Electromagnetic Field Computation. PhD thesis, Johannes Kepler University Linz, 2006.

2. J. Schöberl and S. Zaglmayr. High order Nédélec elements with local complete sequence properties. COMPEL 24(2):374–384, 2005.

3. D. N. Arnold, R. S. Falk and R. Winther. Finite element exterior calculus, homological techniques, and applications. Acta Numer. 15:1–155, 2006.

4. J.-C. Nédélec. Mixed finite elements in \(\mathbb R^3\). Numer. Math. 35:315–341, 1980.

5. J.-C. Nédélec. A new family of mixed finite elements in \(\mathbb R^3\). Numer. Math. 50:57–81, 1986.