Lagrange variants

FIAT/gauss_lobatto_legendre.py

@@ -8,27 +8,10 @@
 #
 # Modified by Pablo D. Brubeck ([email protected]), 2021
 
-from FIAT import finite_element, polynomial_set, lagrange
-from FIAT.reference_element import LINE, TRIANGLE, TETRAHEDRON
-from FIAT.barycentric_interpolation import LagrangePolynomialSet
+from FIAT import lagrange
 
 
-class GaussLobattoLegendre(finite_element.CiarletElement):
+class GaussLobattoLegendre(lagrange.Lagrange):
     """Simplicial continuous element with nodes at the (recursive) Gauss-Lobatto points."""
     def __init__(self, ref_el, degree):
-        if ref_el.shape not in {LINE, TRIANGLE, TETRAHEDRON}:
-            raise ValueError("Gauss-Lobatto-Legendre elements are only defined on simplices.")
-        dual = lagrange.LagrangeDualSet(ref_el, degree, variant="gll")
-        if ref_el.shape == LINE:
-            # In 1D we can use the primal basis as the expansion set,
-            # avoiding any round-off coming from a basis transformation
-            points = []
-            for node in dual.nodes:
-                # Assert singleton point for each node.
-                pt, = node.get_point_dict().keys()
-                points.append(pt)
-            poly_set = LagrangePolynomialSet(ref_el, points)
-        else:
-            poly_set = polynomial_set.ONPolynomialSet(ref_el, degree)
-        formdegree = 0  # 0-form
-        super(GaussLobattoLegendre, self).__init__(poly_set, dual, degree, formdegree)
+        super(GaussLobattoLegendre, self).__init__(ref_el, degree, variant="gll")

FIAT/lagrange.py

@@ -7,14 +7,16 @@
 
 from FIAT import finite_element, polynomial_set, dual_set, functional
 from FIAT.orientation_utils import make_entity_permutations_simplex
+from FIAT.barycentric_interpolation import LagrangePolynomialSet
+from FIAT.reference_element import LINE
 
 
 class LagrangeDualSet(dual_set.DualSet):
     """The dual basis for Lagrange elements.  This class works for
     simplices of any dimension.  Nodes are point evaluation at
     equispaced points."""
 
-    def __init__(self, ref_el, degree, variant=None):
+    def __init__(self, ref_el, degree, variant="equispaced"):
         entity_ids = {}
         nodes = []
         entity_permutations = {}
@@ -44,8 +46,18 @@ def __init__(self, ref_el, degree, variant=None):
 class Lagrange(finite_element.CiarletElement):
     """The Lagrange finite element.  It is what it is."""
 
-    def __init__(self, ref_el, degree):
-        poly_set = polynomial_set.ONPolynomialSet(ref_el, degree)
-        dual = LagrangeDualSet(ref_el, degree)
+    def __init__(self, ref_el, degree, variant="equispaced"):
+        dual = LagrangeDualSet(ref_el, degree, variant=variant)
+        if ref_el.shape == LINE and variant != "equispaced":
+            # In 1D we can use the primal basis as the expansion set,
+            # avoiding any round-off coming from a basis transformation
+            points = []
+            for node in dual.nodes:
+                # Assert singleton point for each node.
+                pt, = node.get_point_dict().keys()
+                points.append(pt)
+            poly_set = LagrangePolynomialSet(ref_el, points)
+        else:
+            poly_set = polynomial_set.ONPolynomialSet(ref_el, degree)
         formdegree = 0  # 0-form
         super(Lagrange, self).__init__(poly_set, dual, degree, formdegree)