Alternative names | Whitney (triangle,tetrahedron), Nédélec, Q H(curl) (quadrilateral,hexahedron), Raviart–Thomas cubical H(curl) (quadrilateral), Nédélec cubical H(curl) (hexahedron) |
De Rham complex families | \(\left[S_{2,k}^\unicode{0x25FA}\right]_{1}\) / \(\mathcal{P}^-_{k}\Lambda^{1}(\Delta_d)\), \(\left[S_{4,k}^\square\right]_{1}\) / \(\mathcal{Q}^-_{k}\Lambda^{1}(\square_d)\) |
Abbreviated names | N1curl, NC, RTce (quadrilateral), Nce (hexahedron) |
Variants | Legendre: Integral moments are taken against orthonormal polynomials Lagrange: Integral moments are taken against (Lagrange)[element:lagrange] basis functions |
Orders | \(1\leqslant k\) |
Reference elements | triangle, tetrahedron, quadrilateral, hexahedron, prism |
Polynomial set | \(\mathcal{P}_{k-1}^d \oplus \mathcal{Z}^{(20)}_{k}\) (triangle, tetrahedron)
\(\mathcal{Q}_{k-1}^d \oplus \mathcal{Z}^{(21)}_{k}\) (quadrilateral, hexahedron)
↓ Show polynomial set definitions ↓↑ Hide polynomial set definitions ↑\(\mathcal{P}_k=\operatorname{span}\left\{\prod_{i=1}^dx_i^{p_i}\middle|\sum_{i=1}^dp_i\leqslant k\right\}\)
\(\mathcal{Z}^{(20)}_k=\left\{\boldsymbol{p}\in\hat{\mathcal{P}}_{k}^d\middle|\boldsymbol{p}(\boldsymbol{x})\cdot \boldsymbol{x}=0\right\}\)
\(\hat{\mathcal{P}}_k=\operatorname{span}\left\{\prod_{i=1}^dx_i^{p_i}\middle|\sum_{i=1}^dp_i=k\right\}=\mathcal{P}_k\setminus\mathcal{P}_{k-1}\)
\(\mathcal{Q}_k=\operatorname{span}\left\{\prod_{i=1}^dx_i^{p_i}\middle|\max_i(p_i)\leqslant k\right\}\)
\(\mathcal{Z}^{(21)}_k=\left\{\boldsymbol{q}\in\hat{\mathcal{Q}}_{k}\middle|\boldsymbol{q}(\boldsymbol{x})\cdot x_i\boldsymbol{e}_i\in\mathcal{Q}_{k}\text{ for }i=1,\dots,d\right\}\)
\(\hat{\mathcal{Q}}_k=\operatorname{span}\left\{\prod_{i=1}^dx_i^{p_i}\middle|\max_i(p_i)=k\right\}=\mathcal{Q}_k\setminus\mathcal{Q}_{k-1}\) |
DOFs | On each edge: tangent integral moments with an order \(k-1\) Lagrange space
On each face (triangle): integral moments with an order \(k-2\) vector Lagrange space
On each face (quadrilateral): integral moments with an order \(k-1\) Q H(div) space
On each volume (tetrahedron): integral moments with an order \(k-3\) vector Lagrange space
On each volume (hexahedron): integral moments with an order \(k-1\) Q H(div) space |
Number of DOFs | triangle: \(k(k+2)\) (A005563) tetrahedron: \(k(k+2)(k+3)/2\) (A005564) quadrilateral: \(2k(k+1)\) (A046092) hexahedron: \(3k(k+1)^2\) (A059986) prism: \(3k(k+2)(k+1)/2\) |
Mapping | covariant Piola |
continuity | Components tangential to facets are continuous |
Categories | Vector-valued elements, H(curl) conforming elements |
Basix | basix.ElementFamily.N1E ↓ Show Basix examples ↓↑ Hide Basix examples ↑Before trying this example, you must install Basix: pip3 install git+https://github.com/FEniCS/basix.git This element can then be created with the following lines of Python: import basix
# Create Nedelec (first kind) (Lagrange variant) order 1 on a triangle element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.triangle, 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a triangle element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.triangle, 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a quadrilateral element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a quadrilateral element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a tetrahedron element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.tetrahedron, 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a tetrahedron element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.tetrahedron, 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a hexahedron element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.hexahedron, 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a hexahedron element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.hexahedron, 2)
# Create Nedelec (first kind) (Legendre variant) order 1 on a triangle element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.triangle, 1, lagrange_variant=basix.LagrangeVariant.legendre)
# Create Nedelec (first kind) (Legendre variant) order 2 on a triangle element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.triangle, 2, lagrange_variant=basix.LagrangeVariant.legendre)
# Create Nedelec (first kind) (Legendre variant) order 1 on a quadrilateral element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 1, lagrange_variant=basix.LagrangeVariant.legendre)
# Create Nedelec (first kind) (Legendre variant) order 2 on a quadrilateral element = basix.create_element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 2, lagrange_variant=basix.LagrangeVariant.legendre) This implementation is correct for all the examples below that it supports. Correct: triangle,1,lagrange; triangle,2,lagrange; quadrilateral,1,lagrange; quadrilateral,2,lagrange; tetrahedron,1,lagrange; tetrahedron,2,lagrange; hexahedron,1,lagrange; hexahedron,2,lagrange; triangle,1,legendre; triangle,2,legendre; quadrilateral,1,legendre; quadrilateral,2,legendre Not implemented: prism,1,lagrange; prism,2,lagrange |
Basix.UFL | basix.ElementFamily.N1E ↓ Show Basix.UFL examples ↓↑ Hide Basix.UFL examples ↑Before trying this example, you must install Basix.UFL: pip3 install git+https://github.com/FEniCS/basix.git git+https://github.com/FEniCS/ufl.git This element can then be created with the following lines of Python: import basix import basix.ufl
# Create Nedelec (first kind) (Lagrange variant) order 1 on a triangle element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.triangle, 1, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a triangle element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.triangle, 2, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a quadrilateral element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 1, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a quadrilateral element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 2, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a tetrahedron element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.tetrahedron, 1, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a tetrahedron element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.tetrahedron, 2, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a hexahedron element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.hexahedron, 1, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a hexahedron element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.hexahedron, 2, lagrange_variant=basix.LagrangeVariant.equispaced)
# Create Nedelec (first kind) (Legendre variant) order 1 on a triangle element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.triangle, 1)
# Create Nedelec (first kind) (Legendre variant) order 2 on a triangle element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.triangle, 2)
# Create Nedelec (first kind) (Legendre variant) order 1 on a quadrilateral element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 1)
# Create Nedelec (first kind) (Legendre variant) order 2 on a quadrilateral element = basix.ufl.element(basix.ElementFamily.N1E, basix.CellType.quadrilateral, 2) This implementation is correct for all the examples below that it supports. Correct: triangle,1,lagrange; triangle,2,lagrange; quadrilateral,1,lagrange; quadrilateral,2,lagrange; tetrahedron,1,lagrange; tetrahedron,2,lagrange; hexahedron,1,lagrange; hexahedron,2,lagrange; triangle,1,legendre; triangle,2,legendre; quadrilateral,1,legendre; quadrilateral,2,legendre Not implemented: prism,1,lagrange; prism,2,lagrange |
Bempp | "SNC" ↓ Show Bempp examples ↓↑ Hide Bempp examples ↑Before trying this example, you must install Bempp: pip3 install bempp-cl This element can then be created with the following lines of Python: import bempp.api grid = bempp.api.shapes.regular_sphere(1) |
FIAT | FIAT.Nedelec(..., variant="integral") ↓ Show FIAT examples ↓↑ Hide FIAT examples ↑Before trying this example, you must install FIAT: pip3 install git+https://github.com/firedrakeproject/fiat.git This element can then be created with the following lines of Python: import FIAT
# Create Nedelec (first kind) (Legendre variant) order 1 element = FIAT.Nedelec(FIAT.ufc_cell("triangle"), 1, variant="integral")
# Create Nedelec (first kind) (Legendre variant) order 2 element = FIAT.Nedelec(FIAT.ufc_cell("triangle"), 2, variant="integral") This implementation is incorrect for this element. |
Symfem | "N1curl", variant="legendre" (triangle, Legendre; tetrahedron, Legendre)
"Qcurl", variant="legendre" (quadrilateral, Legendre; hexahedron, Legendre)
"Ncurl", variant="legendre" (prism, Legendre)
"N1curl" (triangle, Lagrange; tetrahedron, Lagrange)
"Qcurl" (quadrilateral, Lagrange; hexahedron, Lagrange)
"Ncurl" (prism, Lagrange) ↓ Show Symfem examples ↓↑ Hide Symfem examples ↑Before trying this example, you must install Symfem: pip3 install symfem This element can then be created with the following lines of Python: import symfem
# Create Nedelec (first kind) (Lagrange variant) order 1 on a triangle element = symfem.create_element("triangle", "N1curl", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a triangle element = symfem.create_element("triangle", "N1curl", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a quadrilateral element = symfem.create_element("quadrilateral", "Qcurl", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a quadrilateral element = symfem.create_element("quadrilateral", "Qcurl", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a tetrahedron element = symfem.create_element("tetrahedron", "N1curl", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a tetrahedron element = symfem.create_element("tetrahedron", "N1curl", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a hexahedron element = symfem.create_element("hexahedron", "Qcurl", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a hexahedron element = symfem.create_element("hexahedron", "Qcurl", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a prism element = symfem.create_element("prism", "Ncurl", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a prism element = symfem.create_element("prism", "Ncurl", 2)
# Create Nedelec (first kind) (Legendre variant) order 1 on a triangle element = symfem.create_element("triangle", "N1curl", 1, variant="legendre")
# Create Nedelec (first kind) (Legendre variant) order 2 on a triangle element = symfem.create_element("triangle", "N1curl", 2, variant="legendre")
# Create Nedelec (first kind) (Legendre variant) order 1 on a quadrilateral element = symfem.create_element("quadrilateral", "Qcurl", 1, variant="legendre")
# Create Nedelec (first kind) (Legendre variant) order 2 on a quadrilateral element = symfem.create_element("quadrilateral", "Qcurl", 2, variant="legendre") This implementation is used to compute the examples below and verify other implementations. |
(legacy) UFL | "N1curl" (triangle, Lagrange; tetrahedron, Lagrange)
"RTCE" (quadrilateral, Lagrange)
"NCE" (hexahedron, Lagrange) ↓ Show (legacy) UFL examples ↓↑ Hide (legacy) UFL examples ↑Before trying this example, you must install (legacy) UFL: pip3 install git+https://github.com/FEniCS/ufl.git@ufl_legacy This element can then be created with the following lines of Python: import ufl_legacy
# Create Nedelec (first kind) (Lagrange variant) order 1 on a triangle element = ufl_legacy.FiniteElement("N1curl", "triangle", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a triangle element = ufl_legacy.FiniteElement("N1curl", "triangle", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a quadrilateral element = ufl_legacy.FiniteElement("RTCE", "quadrilateral", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a quadrilateral element = ufl_legacy.FiniteElement("RTCE", "quadrilateral", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a tetrahedron element = ufl_legacy.FiniteElement("N1curl", "tetrahedron", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a tetrahedron element = ufl_legacy.FiniteElement("N1curl", "tetrahedron", 2)
# Create Nedelec (first kind) (Lagrange variant) order 1 on a hexahedron element = ufl_legacy.FiniteElement("NCE", "hexahedron", 1)
# Create Nedelec (first kind) (Lagrange variant) order 2 on a hexahedron element = ufl_legacy.FiniteElement("NCE", "hexahedron", 2) |