Add transition element (#103)

* check that the order is within allowable range

* pass kwargs rather that every element taking a variant

* default variant

* remove variant test

* remove variant from demo

* Add transition element

* add transition element to docstring

Author: mscroggs

demo/custom_element.py

@@ -8,7 +8,7 @@
 class CustomElement(CiarletElement):
     """Custom element on a quadrilateral."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order):
         zero = sympy.Integer(0)
         one = sympy.Integer(1)
         half = sympy.Rational(1, 2)

symfem/create.py

@@ -111,7 +111,7 @@ def create_reference(cell_type, vertices=None):
         raise ValueError(f"Unknown cell type: {cell_type}")
 
 
-def create_element(cell_type, element_type, order, variant="equispaced"):
+def create_element(cell_type, element_type, order, **kwargs):
     """Make a finite element.
 
     Parameters
@@ -172,19 +172,37 @@ def create_element(cell_type, element_type, order, variant="equispaced"):
         Bogner-Fox-Schmit, BFS,
         Fortin-Soulie, FS,
         Bernardi-Raugel,
-        Wu-Xu
+        Wu-Xu,
+        transition
     order : int
         The order of the element.
-    variant : str
-        The arrangement type of the points used the define the space.
-
-        Supported values:
-        equispaced, lobatto, radau, legendre
     """
     reference = create_reference(cell_type)
 
     if element_type in _elementmap:
-        assert reference.name in _elementmap[element_type]
-        return _elementmap[element_type][reference.name](reference, order, variant=variant)
+        if reference.name not in _elementmap[element_type]:
+            raise ValueError(f"{element_type} element cannot be created on a {reference.name}.")
+        element_class = _elementmap[element_type][reference.name]
+        if not _order_is_allowed(element_class, reference.name, order):
+            raise ValueError(f"Order {order} {element_type} element cannot be created.")
+        return element_class(reference, order, **kwargs)
 
     raise ValueError(f"Unsupported element type: {element_type}")
+
+
+def _order_is_allowed(element_class, ref, order):
+    if hasattr(element_class, "min_order"):
+        if isinstance(element_class.min_order, dict):
+            if ref in element_class.min_order:
+                if order < element_class.min_order[ref]:
+                    return False
+        elif order < element_class.min_order:
+            return False
+    if hasattr(element_class, "max_order"):
+        if isinstance(element_class.max_order, dict):
+            if ref in element_class.max_order:
+                if order > element_class.max_order[ref]:
+                    return False
+        elif order > element_class.max_order:
+            return False
+    return True

symfem/elements/argyris.py

@@ -15,7 +15,7 @@
 class Argyris(CiarletElement):
     """Argyris finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order):
         from symfem import create_reference
         assert order == 5
         assert reference.name == "triangle"

symfem/elements/aw.py

@@ -17,7 +17,7 @@
 class ArnoldWinther(CiarletElement):
     """An Arnold-Winther element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order, variant="equispaced"):
         from symfem import create_reference
         assert reference.name == "triangle"
         poly = [(p[0], p[1], p[1], p[2])

symfem/elements/bddm.py

@@ -33,14 +33,13 @@ def bddf_polyset(reference, order):
 class BDDF(CiarletElement):
     """Brezzi-Douglas-Duran-Fortin Hdiv finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order, variant="equispaced"):
         poly = bddf_polyset(reference, order)
 
         dofs = make_integral_moment_dofs(
             reference,
-            facets=(NormalIntegralMoment, DiscontinuousLagrange, order),
-            cells=(IntegralMoment, VectorDiscontinuousLagrange, order - 2),
-            variant=variant
+            facets=(NormalIntegralMoment, DiscontinuousLagrange, order, {"variant": variant}),
+            cells=(IntegralMoment, VectorDiscontinuousLagrange, order - 2, {"variant": variant})
         )
 
         super().__init__(reference, order, poly, dofs, reference.tdim, reference.tdim)

symfem/elements/bdfm.py

@@ -47,14 +47,13 @@ def bdfm_polyset(reference, order):
 class BDFM(CiarletElement):
     """Brezzi-Douglas-Fortin-Marini Hdiv finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order, variant="equispaced"):
         poly = bdfm_polyset(reference, order)
 
         dofs = make_integral_moment_dofs(
             reference,
-            facets=(NormalIntegralMoment, DiscontinuousLagrange, order - 1),
-            cells=(IntegralMoment, VectorDiscontinuousLagrange, order - 2),
-            variant=variant
+            facets=(NormalIntegralMoment, DiscontinuousLagrange, order - 1, {"variant": variant}),
+            cells=(IntegralMoment, VectorDiscontinuousLagrange, order - 2, {"variant": variant}),
         )
 
         super().__init__(reference, order, poly, dofs, reference.tdim, reference.tdim)

symfem/elements/bdm.py

@@ -15,14 +15,13 @@
 class BDM(CiarletElement):
     """Brezzi-Douglas-Marini Hdiv finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order, variant="equispaced"):
         poly = polynomial_set(reference.tdim, reference.tdim, order)
 
         dofs = make_integral_moment_dofs(
             reference,
-            facets=(NormalIntegralMoment, DiscontinuousLagrange, order),
-            cells=(IntegralMoment, NedelecFirstKind, order - 1),
-            variant=variant
+            facets=(NormalIntegralMoment, DiscontinuousLagrange, order, {"variant": variant}),
+            cells=(IntegralMoment, NedelecFirstKind, order - 1, {"variant": variant}),
         )
 
         super().__init__(reference, order, poly, dofs, reference.tdim, reference.tdim)

symfem/elements/bell.py

@@ -14,7 +14,7 @@
 class Bell(CiarletElement):
     """Bell finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order, variant="equispaced"):
         assert reference.name == "triangle"
         assert order == 5
         dofs = []
@@ -27,8 +27,7 @@ def __init__(self, reference, order, variant):
             dofs.append(DerivativePointEvaluation(v, (0, 2), entity=(0, v_n)))
         dofs += make_integral_moment_dofs(
             reference,
-            edges=(NormalDerivativeIntegralMoment, DiscontinuousLagrange, 0),
-            variant=variant
+            edges=(NormalDerivativeIntegralMoment, DiscontinuousLagrange, 0, {"variant": variant}),
         )
 
         super().__init__(

symfem/elements/bernardi_raugel.py

@@ -17,7 +17,7 @@
 class BernardiRaugel(CiarletElement):
     """Bernardi-Raugel Hdiv finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order):
         poly = polynomial_set(reference.tdim, reference.tdim, order)
 
         p = Lagrange(reference, 1, variant="equispaced")
@@ -41,8 +41,8 @@ def __init__(self, reference, order, variant):
 
         dofs += make_integral_moment_dofs(
             reference,
-            facets=(NormalIntegralMoment, DiscontinuousLagrange, 0, "contravariant"),
-            variant=variant
+            facets=(NormalIntegralMoment, DiscontinuousLagrange, 0, "contravariant",
+                    {"variant": "equispaced"}),
         )
 
         if order > 1:

symfem/elements/bernstein.py

@@ -68,7 +68,7 @@ def bernstein_polynomials(n, d):
 class Bernstein(CiarletElement):
     """Bernstein finite element."""
 
-    def __init__(self, reference, order, variant):
+    def __init__(self, reference, order):
 
         dofs = [
             IntegralAgainst(reference, p, entity=(reference.tdim, 0))