diff --git a/python/dune/perftool/options.py b/python/dune/perftool/options.py index b294fbb15cb73b3e07f981fbce33d1b964009fc3..75d879949e6930313b69faec4022bd70a8287bee 100644 --- a/python/dune/perftool/options.py +++ b/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") diff --git a/python/dune/perftool/pdelab/driver.py b/python/dune/perftool/pdelab/driver.py index 6d47731038bb96652cf7f50fa27f3791a224423e..9b8efbea778fd2b510fa0534c996363dd242f173 100644 --- a/python/dune/perftool/pdelab/driver.py +++ b/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") diff --git a/python/dune/perftool/pdelab/geometry.py b/python/dune/perftool/pdelab/geometry.py index 63ffcd08ad8362c276b54034f9bc4a4e7a74ca47..5b2049b6cafb821a85bf99d1dc5fa68d1e9f1290 100644 --- a/python/dune/perftool/pdelab/geometry.py +++ b/python/dune/perftool/pdelab/geometry.py @@ -1,14 +1,18 @@ 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