Adds micro to macro conversion for 2d boundaries, not quite working

python/dune/perftool/blockstructured/__init__.py

@@ -4,15 +4,15 @@ import dune.perftool.blockstructured.geometry
 import dune.perftool.blockstructured.spaces
 import dune.perftool.blockstructured.basis
 
+from dune.perftool.pdelab.quadrature import pymbolic_quadrature_position
+# from dune.perftool.pdelab.geometry import to_global
 from dune.perftool.blockstructured.spaces import lfs_inames
-from dune.perftool.blockstructured.quadrature import pymbolic_quadrature_position
 from dune.perftool.blockstructured.geometry import (pymbolic_jacobian_inverse_transposed,
                                                     pymbolic_jacobian_determinant,
                                                     pymbolic_facet_jacobian_determinant,
                                                     to_global)
 
 from dune.perftool.pdelab import PDELabInterface
-from dune.perftool.pdelab.quadrature import pymbolic_quadrature_position
 
 
 class BlockStructuredInterface(PDELabInterface):

python/dune/perftool/blockstructured/tools.py

+from dune.perftool.ufl.modified_terminals import Restriction
 from dune.perftool.generation import (iname,
-                                      domain)
-from dune.perftool.pdelab.geometry import local_dimension
+                                      domain,
+                                      get_global_context_value,
+                                      temporary_variable,
+                                      instruction)
+from dune.perftool.pdelab.geometry import (local_dimension,
+                                           world_dimension,
+                                           name_in_cell_coordinates,
+                                           name_localcenter,
+                                           name_in_cell_geometry)
+
+from dune.perftool.pdelab.quadrature import (pymbolic_quadrature_position_in_cell,
+                                             quadrature_preamble,
+                                             quadrature_inames)
 from dune.perftool.generation.counter import get_counted_variable
 from dune.perftool.options import get_option
 import pymbolic.primitives as prim
@@ -17,6 +29,17 @@ def sub_element_inames():
         domain(name_counted, get_option("number_of_blocks"))
     return inames
 
+
+def do_stupid_stuff_until_i_find_a_better_solution():
+    localcenter = name_localcenter()
+    name = localcenter + "_in_cell"
+    temporary_variable(name, shape_impl=('fv',), shape=(world_dimension(),))
+    geo = name_in_cell_geometry(Restriction.NEGATIVE)
+    instruction(code="{} = {}.global({});".format(name, geo, localcenter),
+                assignees=frozenset({name}))
+    return name
+
+
 # TODO 3d
 def micro_index_to_macro_index(index):
     if not isinstance(index, prim.Variable):
@@ -26,12 +49,49 @@ def micro_index_to_macro_index(index):
 
     k = get_option("number_of_blocks")
 
-    modified_index = prim.Sum((prim.Variable(subelem_inames[0]), index))
+    it = get_global_context_value("integral_type")
+    if it == "cell":
+        x = subelem_inames[0]
+        if local_dimension()>1:
+            y = subelem_inames[1]
+            if world_dimension()>2:
+                z = subelem_inames[2]
+    elif it == "exterior_facet" or it == "interior_facet":
+        center = do_stupid_stuff_until_i_find_a_better_solution()
+        predicate = lambda x: frozenset([prim.Comparison(prim.Subscript(prim.Variable(center),0),x, 0.5)])
+        x = "x"
+        temporary_variable(x)
+        instruction(assignee=prim.Variable(x),
+                    expression=prim.Variable(subelem_inames[0]),
+                    within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
+                    predicates=predicate("==")
+                    )
+        instruction(assignee=prim.Variable(x),
+                    expression=prim.Product(((k-1),prim.Subscript(prim.Variable(center), (0,)))),
+                    within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
+                    predicates=predicate("!=")
+                    )
+        y = "y"
+        temporary_variable(y)
+        instruction(assignee=prim.Variable(y),
+                    expression=prim.Product(((k-1),prim.Subscript(prim.Variable(center), (1,)))),
+                    within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
+                    predicates=predicate("==")
+                    )
+        instruction(assignee=prim.Variable(y),
+                    expression=prim.Variable(subelem_inames[0]),
+                    within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
+                    predicates=predicate("!=")
+                    )
+
+
+
+    modified_index = prim.Sum((prim.Variable(x), index))
     if local_dimension()>1:
-        modified_index = prim.Sum((modified_index, prim.Product(((k+1), prim.Variable(subelem_inames[1])))))
-        index_div_2 = prim.FloorDiv(index, 2)
-        index_div_2 = prim.Product((index_div_2, k-1))
-        modified_index = prim.Sum((modified_index, index_div_2))
+        modified_index = prim.Sum((modified_index, prim.Product(((k+1), prim.Variable(y)))))
+    index_div_2 = prim.FloorDiv(index, 2)
+    index_div_2 = prim.Product((index_div_2, k-1))
+    modified_index = prim.Sum((modified_index, index_div_2))
 
     return modified_index