Implement constant_transformation_matrix option

python/dune/perftool/options.py

@@ -36,7 +36,8 @@ def get_form_compiler_arguments():
     parser.add_argument("--print-transformations", action="store_true", help="print out dot files after ufl tree transformations")
     parser.add_argument("--print-transformations-dir", type=str, help="place where to put dot files (can be omitted)")
     parser.add_argument("--quadrature-order", type=int, help="Quadrature order used for all integrals.")
-    parser.add_argument("--diagonal-transformation-matrix", action="store_true", help="set option if the jacoby of the transformation is diagonal (axiparallel grids)")
+    parser.add_argument("--diagonal-transformation-matrix", action="store_true", help="set option if the jacobian of the transformation is diagonal (axiparallel grids)")
+    parser.add_argument("--constant-transformation-matrix", action="store_true", help="set option if the jacobian of the transformation is constant on a cell")
     parser.add_argument("--ini-file", type=str, help="An inifile to use. A generated driver will be hard-coded to it, a [formcompiler] section will be used as default values to form compiler arguments (use snake case)")
     parser.add_argument("--timer", action="store_true", help="measure times")
     parser.add_argument("--project-basedir", type=str, help="The base (build) directory of the dune-perftool project")

python/dune/perftool/pdelab/driver.py

@@ -176,9 +176,10 @@ def name_dimension():
 def typedef_grid(name):
     dim = name_dimension()
     if any(_driver_data['form'].ufl_cell().cellname() in x for x in ["vertex", "interval", "quadrilateral", "hexahedron"]):
-        # For Yasp Grids the jacobi of the transformation is diagonal
+        # For Yasp Grids the jacobi of the transformation is diagonal and constant on each cell
         from dune.perftool.options import set_option
         set_option('diagonal_transformation_matrix', True)
+        set_option('constant_transformation_matrix', True)
 
         gridt = "Dune::YaspGrid<{}>".format(dim)
         include_file("dune/grid/yaspgrid.hh", filetag="driver")

python/dune/perftool/pdelab/geometry.py

 from dune.perftool.ufl.modified_terminals import Restriction
 from dune.perftool.pdelab.restriction import restricted_name
-from dune.perftool.generation import (cached,
+from dune.perftool.generation import (backend,
+                                      cached,
                                       domain,
                                       get_backend,
                                       get_global_context_value,
+                                      globalarg,
                                       iname,
                                       include_file,
                                       preamble,
                                       temporary_variable,
+                                      valuearg,
                                       )
+from dune.perftool.options import option_switch
 from dune.perftool.pdelab.quadrature import (pymbolic_quadrature_position_in_cell,
                                              quadrature_preamble,
                                              )
@@ -17,6 +21,8 @@ from ufl.algorithms import MultiFunction
 from pymbolic.primitives import Variable
 from pymbolic.primitives import Expression as PymbolicExpression
 
+import numpy as np
+
 
 @preamble
 def define_reference_element(name):
@@ -281,6 +287,22 @@ def define_jacobian_inverse_transposed_temporary(restriction):
     return _define_jacobian_inverse_transposed_temporary
 
 
+@preamble
+@backend(interface="define_jit", name="constant_transformation_matrix")
+def define_constant_jacobian_inveser_transposed(name, restriction):
+    geo = name_cell_geometry(restriction)
+    pos = name_localcenter()
+    dim = name_dimension()
+
+    globalarg(name, dtype=np.float64, shape=(dim, dim), managed=False)
+
+    return 'auto {} = {}.jacobianInverseTransposed({});'.format(name,
+                                                                geo,
+                                                                pos,
+                                                                )
+
+
+@backend(interface="define_jit", name="default")
 def define_jacobian_inverse_transposed(name, restriction):
     dim = name_dimension()
     temporary_variable(name, decl_method=define_jacobian_inverse_transposed_temporary(restriction), shape=(dim, dim))
@@ -297,10 +319,25 @@ def define_jacobian_inverse_transposed(name, restriction):
 
 def name_jacobian_inverse_transposed(restriction):
     name = restricted_name("jit", restriction)
-    define_jacobian_inverse_transposed(name, restriction)
+    get_backend(interface="define_jit", selector=option_switch("constant_transformation_matrix"))(name, restriction)
     return name
 
 
+@backend(interface="detjac", name="constant_transformation_matrix")
+@preamble
+def define_constant_jacobian_determinant(name):
+    geo = name_geometry()
+    pos = name_localcenter()
+
+    valuearg(name, dtype=np.float64)
+
+    return "auto {} = {}.integrationElement({});".format(name,
+                                                         geo,
+                                                         pos,
+                                                         )
+
+
+@backend(interface="detjac", name="default")
 def define_jacobian_determinant(name):
     temporary_variable(name, shape=())
     geo = name_geometry()
@@ -317,11 +354,11 @@ def define_jacobian_determinant(name):
 
 def name_jacobian_determinant():
     name = 'detjac'
-    define_jacobian_determinant(name)
+    get_backend(interface="detjac", selector=option_switch("constant_transformation_matrix"))(name)
     return name
 
 
 def name_facet_jacobian_determinant():
     name = 'fdetjac'
-    define_jacobian_determinant(name)
+    get_backend(interface="detjac", selector=option_switch("constant_transformation_matrix"))(name)
     return name