From 1886ef250d71e39c80c84b150db8654a17aa9ffb Mon Sep 17 00:00:00 2001
From: Dominic Kempf <dominic.kempf@iwr.uni-heidelberg.de>
Date: Fri, 6 Oct 2017 16:15:04 +0200
Subject: [PATCH] Implement facet jacobian determinants
This time only for sumfact, as we need direction dependency.
---
python/dune/perftool/sumfact/__init__.py | 8 +++++
python/dune/perftool/sumfact/geometry.py | 45 ++++++++++++++++++++++++
2 files changed, 53 insertions(+)
diff --git a/python/dune/perftool/sumfact/__init__.py b/python/dune/perftool/sumfact/__init__.py
index a4fce883..3a96f275 100644
--- a/python/dune/perftool/sumfact/__init__.py
+++ b/python/dune/perftool/sumfact/__init__.py
@@ -67,6 +67,10 @@ class SumFactInterface(PDELabInterface):
def pymbolic_quadrature_weight(self):
return quadrature_weight(self.visitor)
+ #
+ # Geometry related generator functions
+ #
+
def pymbolic_spatial_coordinate(self):
import dune.perftool.sumfact.geometry
ret = get_backend(interface="spatial_coordinate", selector=option_switch("diagonal_transformation_matrix"))(self.visitor.do_predicates, self.visitor)
@@ -83,3 +87,7 @@ class SumFactInterface(PDELabInterface):
ret, indices = pymbolic_unit_inner_normal(self.visitor.indices)
self.visitor.indices = indices
return ret
+
+ def pymbolic_facet_jacobian_determinant(self):
+ from dune.perftool.sumfact.geometry import pymbolic_facet_jacobian_determinant
+ return pymbolic_facet_jacobian_determinant()
diff --git a/python/dune/perftool/sumfact/geometry.py b/python/dune/perftool/sumfact/geometry.py
index 7f94fceb..8df024fe 100644
--- a/python/dune/perftool/sumfact/geometry.py
+++ b/python/dune/perftool/sumfact/geometry.py
@@ -1,6 +1,7 @@
""" Sum factorized geometry evaluations """
from dune.perftool.generation import (backend,
+ class_member,
domain,
get_backend,
get_counted_variable,
@@ -16,6 +17,7 @@ from dune.perftool.pdelab.geometry import (local_dimension,
world_dimension,
name_geometry,
)
+from dune.perftool.sumfact.switch import get_facedir
from dune.perftool.sumfact.symbolic import SumfactKernelInputBase
from dune.perftool.sumfact.vectorization import attach_vectorization_info
from dune.perftool.options import get_option, option_switch
@@ -194,3 +196,46 @@ def pymbolic_unit_inner_normal(visitor_indices):
else:
from dune.perftool.pdelab.geometry import pymbolic_unit_inner_normal as _norm
return _norm(), visitor_indices
+
+
+def pymbolic_facet_jacobian_determinant():
+ if get_option("constant_transformation_matrix"):
+ return pymbolic_constant_facet_jacobian_determinant()
+ else:
+ from dune.perftool.pdelab.geometry import pymbolic_facet_jacobian_determinant as _norm
+ return _norm()
+
+
+def pymbolic_constant_facet_jacobian_determinant():
+ facedir = get_facedir(Restriction.NEGATIVE)
+ assert isinstance(facedir, int)
+
+ name = "fdetjac"
+ define_constant_facet_jacobian_determinant(name)
+ globalarg(name, dtype=np.float64, shape=(world_dimension(),))
+ return prim.Subscript(prim.Variable(name), (facedir,))
+
+
+@class_member(classtag="operator")
+def define_constant_facet_jacobian_determinant(name):
+ define_constant_facet_jacobian_determinant_eval(name)
+ from dune.perftool.pdelab.localoperator import lop_template_ansatz_gfs
+ gfst = lop_template_ansatz_gfs()
+ return "std::array<typename {}::Traits::GridView::template Codim<0>::Geometry::ctype, {}> {};".format(gfst, world_dimension(), name)
+
+
+@kernel_cached(kernel="operator")
+def define_constant_facet_jacobian_determinant_eval(name):
+ from dune.perftool.pdelab.localoperator import name_ansatz_gfs_constructor_param, lop_template_range_field
+ gfs = name_ansatz_gfs_constructor_param()
+ rft = lop_template_range_field()
+ code = ["{",
+ " auto e = *({}.gridView().template begin<0>());".format(gfs),
+ " int dir=0;",
+ " for(auto is = {0}.gridView().ibegin(e); is != {0}.gridView().iend(e); ++(++is), ++dir)".format(gfs),
+ " {}[dir] = is->geometry().integrationElement(Dune::FieldVector<{}, {}>());".format(name, rft, world_dimension() - 1),
+ "}",
+ ]
+ instruction(code="\n".join(code),
+ kernel="operator",
+ )
--
GitLab