diff --git a/python/dune/perftool/compile.py b/python/dune/perftool/compile.py index ec456ec942f5477f6feb7b6d7cf52cce935d4d4a..cfb4b9f3d533db42371e39cc55f3b4695b37dd84 100644 --- a/python/dune/perftool/compile.py +++ b/python/dune/perftool/compile.py @@ -20,7 +20,7 @@ from dune.perftool.pdelab.localoperator import (generate_localoperator_basefile, name_localoperator_file) from dune.perftool.ufl.preprocess import preprocess_form -import os.path +from os.path import splitext, basename # Disable loopy caching before we do anything else! diff --git a/python/dune/perftool/pdelab/driver.py b/python/dune/perftool/pdelab/driver.py index 39ef2e247a8206c78339d3e6c6dad0ec242d8c50..8303cf8a0ea9d6afd8636c468bcdc90658e59066 100644 --- a/python/dune/perftool/pdelab/driver.py +++ b/python/dune/perftool/pdelab/driver.py @@ -8,6 +8,7 @@ gained there. """ from dune.perftool.error import PerftoolCodegenError from dune.perftool.generation import (generator_factory, + global_context, include_file, cached, preamble, @@ -748,13 +749,24 @@ def define_vector(name, formdata): @preamble def define_boundary_lambda(boundary, name): + from dune.perftool.ufl.execution import Expression + from ufl.classes import Expr if boundary is None: return "auto {} = [&](const auto& x){{ return 0.0; }};".format(name) - if boundary.is_global: - return "auto {} = [&](const auto& x){{ {} }};".format(name, boundary.c_expr[0]) - else: - return "auto {} = [&](const auto& e, const auto& x){{ {} }};".format(name, boundary.c_expr[0]) - + elif isinstance(boundary, Expression): + if boundary.is_global: + return "auto {} = [&](const auto& x){{ {} }};".format(name, boundary.c_expr[0]) + else: + return "auto {} = [&](const auto& e, const auto& x){{ {} }};".format(name, boundary.c_expr[0]) + elif isinstance(boundary, Expr): + # Set up a visitor + with global_context(integral_type="exterior_facet", formdata=_driver_data["formdata"], driver=True): + from dune.perftool.ufl.visitor import UFL2LoopyVisitor + from dune.perftool.pdelab import PDELabInterface + visitor = UFL2LoopyVisitor(PDELabInterface(), "exterior_facet", {}) + return "auto {} = [&](const auto& x){{ return {}; }};".format(name, visitor(boundary)) + + raise ValueError("Expression not understood") def name_boundary_lambda(boundary, name): define_boundary_lambda(boundary, name + "lambda") diff --git a/python/dune/perftool/pdelab/geometry.py b/python/dune/perftool/pdelab/geometry.py index 1be633ff09c7c5bd8ab10554f0e3331ebcde20a0..bcb17d07e46928f3d9ae4cff02777be50b5d3a37 100644 --- a/python/dune/perftool/pdelab/geometry.py +++ b/python/dune/perftool/pdelab/geometry.py @@ -387,7 +387,7 @@ def name_facet_jacobian_determinant(): def apply_to_global_transformation(name, local): - temporary_variable(name, shape=(name_dimension(),), shape_impl=("fv",)) + temporary_variable(name, shape=(world_dimension(),), shape_impl=("fv",)) geo = name_geometry() code = "{} = {}.global({});".format(name, geo, diff --git a/python/dune/perftool/ufl/execution.py b/python/dune/perftool/ufl/execution.py index 57d337fbe541beb5de6afe429c5a3930dccf962e..112aa5ec2ed30a1fce5e2689fa93ef42a7f56759 100644 --- a/python/dune/perftool/ufl/execution.py +++ b/python/dune/perftool/ufl/execution.py @@ -145,13 +145,14 @@ class FiniteElement(ufl.FiniteElement): self.dirichlet_constraints = kwargs.pop('dirichlet_constraints', 'return true;') if isinstance(self.dirichlet_constraints, str): self.dirichlet_constraints = Expression(self.dirichlet_constraints) - assert isinstance(self.dirichlet_constraints, Expression) + from ufl.classes import Expr + assert isinstance(self.dirichlet_constraints, (Expression, Expr)) # Get dirichlet_constraints and convert it to Expression if necessary! self.dirichlet_expression = kwargs.pop('dirichlet_expression', 'return 0.0;') if isinstance(self.dirichlet_expression, str): self.dirichlet_expression = Expression(self.dirichlet_expression) - assert isinstance(self.dirichlet_expression, Expression) + assert isinstance(self.dirichlet_expression, (Expression, Expr)) # Initialize the original finite element from ufl ufl.FiniteElement.__init__(self, *args, **kwargs) diff --git a/python/dune/perftool/ufl/visitor.py b/python/dune/perftool/ufl/visitor.py index a476668f2c6dc52d1b861c7394f4cced3a1fc742..0ad6a9d33463c70ee58e04f9d1acc240f886b8b8 100644 --- a/python/dune/perftool/ufl/visitor.py +++ b/python/dune/perftool/ufl/visitor.py @@ -305,7 +305,11 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker): # def spatial_coordinate(self, o): - return self.interface.pymbolic_spatial_coordinate(self.restriction) + # If this is called from the driver, we just want to return x + if get_global_context_value("driver", False): + return prim.Variable("x") + else: + return self.interface.pymbolic_spatial_coordinate(self.restriction) def facet_normal(self, o): # The normal must be restricted to be well-defined diff --git a/test/poisson/poisson.ufl b/test/poisson/poisson.ufl index 27395230801bb9fa14dd29ba83749a591a2af7e5..73515675141be1bd95dfd4d3d65d7c7195bb5bb3 100644 --- a/test/poisson/poisson.ufl +++ b/test/poisson/poisson.ufl @@ -1,13 +1,14 @@ -g = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return std::exp(-1.*c.two_norm2());") cell = triangle +x = SpatialCoordinate(cell) + +c = (0.5-x[0])**2 + (0.5-x[1])**2 +g = exp(-1.*c) +f = 4*(1.-c)*g + V = FiniteElement("CG", cell, 1, dirichlet_expression=g) u = TrialFunction(V) v = TestFunction(V) -x = SpatialCoordinate(cell) - -c = (0.5-x[0])**2 + (0.5-x[1])**2 -f = 4*(1.-c)*exp(-1.*c) forms = [(inner(grad(u), grad(v)) - f*v)*dx]