Update Basix interpolation matrix shapes

ffcx/codegeneration/basix_custom_element_template.py

@@ -28,6 +28,7 @@
   .map_type = {map_type},
   .discontinuous = {discontinuous},
   .highest_complete_degree = {highest_complete_degree},
+  .interpolation_nderivs = {interpolation_nderivs},
   .highest_degree = {highest_degree}
 }};
 

ffcx/codegeneration/finite_element.py

@@ -123,6 +123,7 @@ def generate_custom_element(name, ir):
     d["highest_complete_degree"] = ir.highest_complete_degree
     d["highest_degree"] = ir.highest_degree
     d["discontinuous"] = "true" if ir.discontinuous else "false"
+    d["interpolation_nderivs"] = ir.interpolation_nderivs
 
     import ffcx.codegeneration.C.cnodes as L
 
@@ -158,12 +159,14 @@ def generate_custom_element(name, ir):
     ndofs = []
     M = []
     for entity in ir.M:
-        for matrices in entity:
-            ndofs.append(matrices.shape[0])
-            for mat in matrices:
-                for row in mat:
-                    for i in row:
-                        M.append(i)
+        for mat4d in entity:
+            ndofs.append(mat4d.shape[0])
+            for mat3d in mat4d:
+                for mat2d in mat3d:
+                    for row in mat2d:
+                        for i in row:
+                            M.append(i)
+
     d["ndofs"] = f"ndofs_{name}"
     d["ndofs_init"] = f"int ndofs_{name}[{len(ndofs)}] = "
     d["ndofs_init"] += "{" + ",".join([f" {i}" for i in ndofs]) + "};"

ffcx/codegeneration/ufcx.h

@@ -191,6 +191,9 @@ extern "C"
     /// The highest degree full polynomial space contained in this element
     int highest_complete_degree;
 
+    /// The number of derivatives needed when interpolating
+    int interpolation_nderivs;
+
     /// The highest degree of a polynomial in the element
     int highest_degree;
   } ufcx_basix_custom_finite_element;

ffcx/element_interface.py

@@ -281,6 +281,12 @@ def discontinuous(self) -> bool:
         """True if the discontinuous version of the element is used."""
         pass
 
+    @property
+    @abstractmethod
+    def interpolation_nderivs(self) -> int:
+        """The number of derivatives needed when interpolating."""
+        pass
+
     @property
     def is_custom_element(self) -> bool:
         """True if the element is a custom Basix element."""
@@ -378,6 +384,10 @@ def cell_type(self):
     def discontinuous(self):
         return self.element.discontinuous
 
+    @property
+    def interpolation_nderivs(self) -> int:
+        return self.element.interpolation_nderivs
+
     @property
     def is_custom_element(self) -> bool:
         """True if the element is a custom Basix element."""
@@ -477,6 +487,10 @@ def cell_type(self):
     def discontinuous(self):
         return self.element.discontinuous
 
+    @property
+    def interpolation_nderivs(self) -> int:
+        return self.element.interpolation_nderivs
+
 
 class MixedElement(BaseElement):
     """A mixed element that combines two or more elements."""
@@ -603,6 +617,10 @@ def cell_type(self):
     def discontinuous(self):
         return False
 
+    @property
+    def interpolation_nderivs(self) -> int:
+        return max([e.interpolation_nderivs for e in self.sub_elements])
+
 
 class BlockedElement(BaseElement):
     """An element with a block size that contains multiple copies of a sub element."""
@@ -711,6 +729,10 @@ def cell_type(self):
     def discontinuous(self):
         return self.sub_element.discontinuous
 
+    @property
+    def interpolation_nderivs(self) -> int:
+        return self.sub_element.interpolation_nderivs
+
 
 class QuadratureElement(BaseElement):
     """A quadrature element."""
@@ -821,3 +843,7 @@ def cell_type(self) -> None:
     @property
     def discontinuous(self):
         return False
+
+    @property
+    def interpolation_nderivs(self) -> int:
+        return 0

ffcx/ir/representation.py

@@ -60,7 +60,7 @@
     'original_coefficient_positions', 'points', 'coefficient_names', 'constant_names', 'needs_facet_permutations',
     'function_spaces', 'name_from_uflfile'])
 ir_custom_element = namedtuple('ir_custom_element', [
-    'cell_type', 'value_shape', 'wcoeffs', 'x', 'M', 'map_type',
+    'cell_type', 'value_shape', 'wcoeffs', 'x', 'M', 'map_type', 'interpolation_nderivs',
     'discontinuous', 'highest_complete_degree', 'highest_degree'])
 
 ir_data = namedtuple('ir_data', ['elements', 'dofmaps', 'integrals', 'forms', 'expressions'])
@@ -177,6 +177,7 @@ def _compute_custom_element_ir(basix_element):
     ir["M"] = basix_element.M
     ir["map_type"] = basix_element.map_type
     ir["discontinuous"] = basix_element.discontinuous
+    ir["interpolation_nderivs"] = basix_element.interpolation_nderivs
     ir["highest_complete_degree"] = basix_element.highest_complete_degree
     ir["highest_degree"] = basix_element.highest_degree