From a0380133509f62fa8d87c54d592b44ae0a468ba7 Mon Sep 17 00:00:00 2001
From: Marcel Koch <marcel.koch@uni-muenster.de>
Date: Thu, 13 Jul 2017 09:36:17 +0200
Subject: [PATCH] Moves application of transformation from micro to reference
 element from visitor to PDELab interface

---
 .../dune/perftool/blockstructured/__init__.py |  1 +
 .../dune/perftool/blockstructured/geometry.py | 17 +++++++++++
 python/dune/perftool/pdelab/__init__.py       | 19 +++++-------
 python/dune/perftool/pdelab/basis.py          | 30 +++++--------------
 python/dune/perftool/pdelab/geometry.py       | 10 +++++++
 python/dune/perftool/ufl/preprocess.py        |  4 ---
 python/dune/perftool/ufl/visitor.py           |  6 ++--
 7 files changed, 46 insertions(+), 41 deletions(-)
 create mode 100644 python/dune/perftool/blockstructured/geometry.py

diff --git a/python/dune/perftool/blockstructured/__init__.py b/python/dune/perftool/blockstructured/__init__.py
index 50214de6..b82070c9 100644
--- a/python/dune/perftool/blockstructured/__init__.py
+++ b/python/dune/perftool/blockstructured/__init__.py
@@ -1,5 +1,6 @@
 import dune.perftool.blockstructured.accumulation
 import dune.perftool.blockstructured.quadrature
 import dune.perftool.blockstructured.argument
+import dune.perftool.blockstructured.geometry
 import dune.perftool.blockstructured.spaces
 import dune.perftool.blockstructured.basis
diff --git a/python/dune/perftool/blockstructured/geometry.py b/python/dune/perftool/blockstructured/geometry.py
new file mode 100644
index 00000000..565f74d3
--- /dev/null
+++ b/python/dune/perftool/blockstructured/geometry.py
@@ -0,0 +1,17 @@
+from dune.perftool.generation import backend
+from dune.perftool.options import get_option
+from dune.perftool.pdelab.geometry import (name_jacobian_determinant,
+                                           name_jacobian_inverse_transposed)
+import pymbolic.primitives as prim
+
+
+@backend(interface="pymbolic_jacobian_determinant", name="blockstructured")
+def pymbolic_jacobian_determinant():
+    return prim.Quotient(prim.Variable(name_jacobian_determinant()),
+                         prim.Power(get_option("number_of_blocks"),2))
+
+
+@backend(interface="pymbolic_jacobian_inverse_transposed", name="blockstructured")
+def pymbolic_jacobian_inverse_transposed(i,j,restriction):
+    return prim.Product((get_option("number_of_blocks"),
+                         prim.Subscript(prim.Variable(name_jacobian_inverse_transposed(restriction)), (j, i))))
diff --git a/python/dune/perftool/pdelab/__init__.py b/python/dune/perftool/pdelab/__init__.py
index 32cec784..b232b899 100644
--- a/python/dune/perftool/pdelab/__init__.py
+++ b/python/dune/perftool/pdelab/__init__.py
@@ -3,7 +3,7 @@
 # Trigger some imports that are needed to have all backend implementations visible
 # to the selection mechanisms
 from dune.perftool.generation import (get_backend)
-from dune.perftool.options import get_option
+from dune.perftool.options import option_switch
 from dune.perftool.pdelab.argument import (pymbolic_apply_function,
                                            pymbolic_apply_function_gradient,
                                            pymbolic_trialfunction,
@@ -17,8 +17,8 @@ from dune.perftool.pdelab.geometry import (component_iname,
                                            name_cell_volume,
                                            name_facet_area,
                                            name_facet_jacobian_determinant,
-                                           name_jacobian_determinant,
-                                           name_jacobian_inverse_transposed,
+                                           pymbolic_jacobian_determinant,
+                                           pymbolic_jacobian_inverse_transposed,
                                            name_unit_inner_normal,
                                            name_unit_outer_normal,
                                            to_global,
@@ -58,10 +58,7 @@ class PDELabInterface(object):
     #
 
     def lfs_inames(self, element, restriction, number=None, context=''):
-        if get_option("blockstructured"):
-            return get_backend("lfs_inames", selector=lambda: "blockstructured")(element, restriction, number, context)
-        else:
-            return get_backend("lfs_inames")(element, restriction, number, context)
+        return get_backend("lfs_inames", selector=option_switch("blockstructured"))(element, restriction, number, context)
 
     #
     # Test and trial function related generator functions
@@ -115,11 +112,11 @@ class PDELabInterface(object):
     def name_facet_jacobian_determinant(self):
         return name_facet_jacobian_determinant()
 
-    def name_jacobian_determinant(self):
-        return name_jacobian_determinant()
+    def pymbolic_jacobian_determinant(self):
+        return get_backend("pymbolic_jacobian_determinant", selector=option_switch("blockstructured"))()
 
-    def name_jacobian_inverse_transposed(self, restriction):
-        return name_jacobian_inverse_transposed(restriction)
+    def pymbolic_jacobian_inverse_transposed(self, i, j, restriction):
+        return get_backend("pymbolic_jacobian_inverse_transposed",selector=option_switch("blockstructured"))(i,j,restriction)
 
     def name_unit_inner_normal(self):
         return name_unit_inner_normal()
diff --git a/python/dune/perftool/pdelab/basis.py b/python/dune/perftool/pdelab/basis.py
index 2475b6ed..92113832 100644
--- a/python/dune/perftool/pdelab/basis.py
+++ b/python/dune/perftool/pdelab/basis.py
@@ -8,7 +8,7 @@ from dune.perftool.generation import (backend,
                                       kernel_cached,
                                       temporary_variable,
                                       )
-from dune.perftool.options import get_option
+from dune.perftool.options import option_switch
 from dune.perftool.pdelab.spaces import (lfs_child,
                                          name_leaf_lfs,
                                          name_lfs,
@@ -79,14 +79,10 @@ def pymbolic_basis(leaf_element, restriction, number, context=''):
     assert leaf_element.num_sub_elements() == 0
     name = "phi_{}".format(FEM_name_mangling(leaf_element))
     name = restricted_name(name, restriction)
-    if get_option("blockstructured"):
-        get_backend("evaluate_basis", selector=lambda: "blockstructured")(leaf_element, name, restriction)
-        iname, = get_backend("lfs_inames", selector=lambda: "blockstructured")(leaf_element, restriction, number, context=context)
-    else:
-        get_backend("evaluate_basis")(leaf_element, name, restriction)
-        iname, = get_backend("lfs_inames")(leaf_element, restriction, number, context=context)
+    get_backend("evaluate_basis", selector=option_switch("blockstructured"))(leaf_element, name, restriction)
+    iname, = get_backend("lfs_inames", selector=option_switch("blockstructured"))(leaf_element, restriction, number, context=context)
 
-    return Subscript(Variable(name), (Variable(iname),))
+    return Subscript(Variable(name), (Variable(iname), ))
 
 
 @backend(interface="evaluate_grad")
@@ -111,12 +107,8 @@ def pymbolic_reference_gradient(leaf_element, restriction, number, context=''):
     assert leaf_element.num_sub_elements() == 0
     name = "js_{}".format(FEM_name_mangling(leaf_element))
     name = restricted_name(name, restriction)
-    if get_option("blockstructured"):
-        get_backend("evaluate_grad", selector=lambda: "blockstructured")(leaf_element, name, restriction)
-        iname, = get_backend("lfs_inames", selector=lambda: "blockstructured")(leaf_element, restriction, number, context=context)
-    else:
-        get_backend("evaluate_grad")(leaf_element, name, restriction)
-        iname, = get_backend("lfs_inames")(leaf_element, restriction, number, context=context)
+    get_backend("evaluate_grad", selector=option_switch("blockstructured"))(leaf_element, name, restriction)
+    iname, = get_backend("lfs_inames", selector=option_switch("blockstructured"))(leaf_element, restriction, number, context=context)
 
     return Subscript(Variable(name), (Variable(iname), 0))
 
@@ -154,10 +146,7 @@ def evaluate_coefficient(element, name, container, restriction, tree_path):
     if isinstance(sub_element, (VectorElement, TensorElement)):
         lfs = lfs_child(lfs, idims, shape=shape_as_pymbolic(shape), symmetry=element.symmetry())
 
-    if get_option('blockstructured'):
-        coeff = get_backend("pymbolic_coefficient", selector=lambda: "blockstructured")(container, lfs, index)
-    else:
-        coeff = get_backend("pymbolic_coefficient")(container, lfs, index)
+    coeff = get_backend("pymbolic_coefficient", selector=option_switch("blockstructured"))(container, lfs, index)
     assignee = Subscript(Variable(name), tuple(Variable(i) for i in idims))
 
     reduction_expr = Product((coeff, basis))
@@ -194,10 +183,7 @@ def evaluate_coefficient_gradient(element, name, container, restriction, tree_pa
     if isinstance(sub_element, (VectorElement, TensorElement)):
         lfs = lfs_child(lfs, idims[:-1], shape=shape_as_pymbolic(shape[:-1]), symmetry=element.symmetry())
 
-    if get_option('blockstructured'):
-        coeff = get_backend("pymbolic_coefficient", selector=lambda: "blockstructured")(container, lfs, index)
-    else:
-        coeff = get_backend("pymbolic_coefficient")(container, lfs, index)
+    coeff = get_backend("pymbolic_coefficient", selector=option_switch("blockstructured"))(container, lfs, index)
     assignee = Subscript(Variable(name), tuple(Variable(i) for i in idims))
 
     reduction_expr = Product((coeff, Subscript(Variable(get_pymbolic_basename(basis)), basis.index + (Variable(idims[-1]),))))
diff --git a/python/dune/perftool/pdelab/geometry.py b/python/dune/perftool/pdelab/geometry.py
index 9e4c7482..90b2ac90 100644
--- a/python/dune/perftool/pdelab/geometry.py
+++ b/python/dune/perftool/pdelab/geometry.py
@@ -352,6 +352,12 @@ def name_jacobian_inverse_transposed(restriction):
     get_backend(interface="define_jit", selector=option_switch("constant_transformation_matrix"))(name, restriction)
     return name
 
+@backend(interface="pymbolic_jacobian_inverse_transposed", name="default")
+def pymbolic_jacobian_inverse_transposed(i,j,restriction):
+    # Dune only has JacobianInverseTransposed as a first class citizen,
+    # so we need to switch the indices around!
+    return prim.Subscript(prim.Variable(name_jacobian_inverse_transposed(restriction)), (j, i))
+
 
 @backend(interface="detjac", name="constant_transformation_matrix")
 @preamble
@@ -392,6 +398,10 @@ def name_jacobian_determinant():
     return name
 
 
+@backend(interface="pymbolic_jacobian_determinant", name="default")
+def pymbolic_jacobian_determinant():
+    return prim.Variable(name_jacobian_determinant())
+
 def name_facet_jacobian_determinant():
     name = 'fdetjac'
     get_backend(interface="detjac", selector=option_switch("constant_transformation_matrix"))(name)
diff --git a/python/dune/perftool/ufl/preprocess.py b/python/dune/perftool/ufl/preprocess.py
index 61337edb..fba501de 100644
--- a/python/dune/perftool/ufl/preprocess.py
+++ b/python/dune/perftool/ufl/preprocess.py
@@ -30,10 +30,6 @@ def apply_default_transformations(form):
     from dune.perftool.ufl.transformations.reindexing import reindexing
     from dune.perftool.ufl.transformations.unroll import unroll_dimension_loops
 
-    if get_option("blockstructured"):
-        from dune.perftool.ufl.transformations.blockstructured import blockstructured
-        form = transform_form(form, blockstructured)
-
     form = transform_form(form, unroll_dimension_loops)
     form = transform_form(form, pushdown_indexed)
     form = transform_form(form, reindexing)
diff --git a/python/dune/perftool/ufl/visitor.py b/python/dune/perftool/ufl/visitor.py
index e1769ae2..0bbb13f5 100644
--- a/python/dune/perftool/ufl/visitor.py
+++ b/python/dune/perftool/ufl/visitor.py
@@ -368,7 +368,7 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker):
         return self.interface.pymbolic_quadrature_weight()
 
     def jacobian_determinant(self, o):
-        return Variable(self.interface.name_jacobian_determinant())
+        return self.interface.pymbolic_jacobian_determinant()
 
     def jacobian_inverse(self, o):
         restriction = self.restriction
@@ -383,9 +383,7 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker):
         if isinstance(i, int) and isinstance(j, int) and i != j:
             return 0
 
-        # Dune only has JacobianInverseTransposed as a first class citizen,
-        # so we need to switch the indices around!
-        return Subscript(Variable(self.interface.name_jacobian_inverse_transposed(restriction)), (j, i))
+        return self.interface.pymbolic_jacobian_inverse_transposed(i, j, restriction)
 
     def jacobian(self, o):
         raise NotImplementedError("How did you get Jacobian into your form? We only support JacobianInverse right now. Report!")
-- 
GitLab