From aa6ffb8c18ad75be6535e6efd5b9b3c0216e2206 Mon Sep 17 00:00:00 2001
From: Dominic Kempf <dominic.kempf@iwr.uni-heidelberg.de>
Date: Fri, 9 Mar 2018 13:42:28 +0100
Subject: [PATCH] Implement possibility of functions that provide their own
visitor
This will be needed to replace nonlinearities in Richards code
with a cubic spline evaluation without bloating dune-perftool with
the implementation details of such procedure.
Also provides a test, that defines a custom square function.
---
python/dune/perftool/ufl/visitor.py | 16 ++++++++++++
test/poisson/CMakeLists.txt | 6 +++++
test/poisson/poisson_customfunction.mini | 18 +++++++++++++
test/poisson/poisson_customfunction.ufl | 33 ++++++++++++++++++++++++
4 files changed, 73 insertions(+)
create mode 100644 test/poisson/poisson_customfunction.mini
create mode 100644 test/poisson/poisson_customfunction.ufl
diff --git a/python/dune/perftool/ufl/visitor.py b/python/dune/perftool/ufl/visitor.py
index 3c4562c9..cbede6c3 100644
--- a/python/dune/perftool/ufl/visitor.py
+++ b/python/dune/perftool/ufl/visitor.py
@@ -343,6 +343,22 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker):
def min_value(self, o):
return self._minmax_impl(min, "min", tuple(self.call(op) for op in o.ufl_operands))
+ def math_function(self, o):
+ # MathFunction is a base class for unary functions. We use this to provide
+ # custom functions. Such a custom functions inherits from it and defines the
+ # following methods:
+ # * visit: This function is called from here to delegate the visiting process
+ # to the user code. The only argument is this visitor instance.
+ # * derivative: It is called from UFL AD code to determine the derivative.
+ # Upstream documentation indicates that FEniCS allows the same
+ # (ab)use of the MathFunction node.
+ # Note that if the __init__ method of your function differs from MathFunction,
+ # you also need to implement the method _ufl_expr_reconstruct_
+ if hasattr(o, "visit"):
+ return o.visit(self)
+ else:
+ raise NotImplementedError("Function {} is not known to dune-perftool.".format(o._name))
+
#
# Handler for conditionals, use pymbolic base implementation
#
diff --git a/test/poisson/CMakeLists.txt b/test/poisson/CMakeLists.txt
index 68c635bc..a05cbe8c 100644
--- a/test/poisson/CMakeLists.txt
+++ b/test/poisson/CMakeLists.txt
@@ -67,6 +67,12 @@ dune_add_formcompiler_system_test(UFLFILE poisson_dg_tensor.ufl
INIFILE poisson_dg_tensor.mini
)
+# 12. Poisson Test Case with a custom function
+dune_add_formcompiler_system_test(UFLFILE poisson_customfunction.ufl
+ BASENAME poisson_customfunction
+ INIFILE poisson_customfunction.mini
+ )
+
# the reference vtk file
add_executable(poisson_dg_ref reference_main.cc)
set_target_properties(poisson_dg_ref PROPERTIES EXCLUDE_FROM_ALL 1)
diff --git a/test/poisson/poisson_customfunction.mini b/test/poisson/poisson_customfunction.mini
new file mode 100644
index 00000000..c45e50fa
--- /dev/null
+++ b/test/poisson/poisson_customfunction.mini
@@ -0,0 +1,18 @@
+__name = poisson_customfunction_{__exec_suffix}
+__exec_suffix = numdiff, symdiff | expand num
+
+lowerleft = 0.0 0.0
+upperright = 1.0 1.0
+elements = 32 32
+elementType = simplical
+
+[wrapper.vtkcompare]
+name = {__name}
+reference = poisson_ref
+extension = vtu
+
+[formcompiler]
+compare_l2errorsquared = 1e-7
+
+[formcompiler.r]
+numerical_jacobian = 1, 0 | expand num
diff --git a/test/poisson/poisson_customfunction.ufl b/test/poisson/poisson_customfunction.ufl
new file mode 100644
index 00000000..026e069d
--- /dev/null
+++ b/test/poisson/poisson_customfunction.ufl
@@ -0,0 +1,33 @@
+import ufl
+
+cell = triangle
+
+x = SpatialCoordinate(cell)
+
+class SquareFct(ufl.classes.MathFunction):
+ def __init__(self, arg):
+ ufl.classes.MathFunction.__init__(self, 'square', arg)
+
+ def _ufl_expr_reconstruct_(self, *operands):
+ return SquareFct(*operands)
+
+ def derivative(self):
+ return 2 * self.ufl_operands[0]
+
+ def visit(self, visitor):
+ op = visitor.call(self.ufl_operands[0])
+ return op * op
+
+
+c = SquareFct(0.5-x[0]) + SquareFct(0.5-x[1])
+g = exp(-1.*c)
+f = 4*(1.-c)*g
+
+V = FiniteElement("CG", cell, 1)
+u = TrialFunction(V)
+v = TestFunction(V)
+
+r = (inner(grad(u), grad(v)) - f*v)*dx
+exact_solution = g
+interpolate_expression = g
+is_dirichlet = 1
\ No newline at end of file
--
GitLab