diff --git a/python/dune/perftool/pdelab/driver/constraints.py b/python/dune/perftool/pdelab/driver/constraints.py index 5c9319726fefbebdaeb57fdcc28f90f05b9a806d..2a093a9cbe4ae723dd7e538ac6ed02fce1d4caac 100644 --- a/python/dune/perftool/pdelab/driver/constraints.py +++ b/python/dune/perftool/pdelab/driver/constraints.py @@ -96,17 +96,11 @@ def name_bctype_lambda(name, func): @preamble def define_intersection_lambda(name, func): - from dune.perftool.ufl.execution import Expression from ufl.classes import Expr if func is None: func = 0. if isinstance(func, (int, float)): return "auto {} = [&](const auto& x){{ return {}; }};".format(name, float(func)) - elif isinstance(func, Expression): - if expression.is_global: - return "auto {} = [&](const auto& x){{ {} }};".format(name, func.c_expr[0]) - else: - return "auto {} = [&](const auto& is, const auto& x){{ {} }};".format(name, func.c_expr[0]) elif isinstance(func, Expr): # Set up a visitor with global_context(integral_type="exterior_facet", formdata=get_formdata(), driver=True): diff --git a/python/dune/perftool/pdelab/driver/interpolate.py b/python/dune/perftool/pdelab/driver/interpolate.py index a3751c84785943fea71e91c06b984ba62687c607..4278294aab8a7655f4a2e852fe46b07cdb440f1d 100644 --- a/python/dune/perftool/pdelab/driver/interpolate.py +++ b/python/dune/perftool/pdelab/driver/interpolate.py @@ -95,15 +95,9 @@ def name_boundary_lambda(boundary): @preamble def define_boundary_lambda(name, boundary): - 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) - 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=get_formdata(), driver=True): diff --git a/python/dune/perftool/ufl/execution.py b/python/dune/perftool/ufl/execution.py index 0f7534515bf5aa990cd1fb12bfe05f27950ee50d..3c3583391f46c1f03a98d12dc839cd377675e5bb 100644 --- a/python/dune/perftool/ufl/execution.py +++ b/python/dune/perftool/ufl/execution.py @@ -30,8 +30,6 @@ class Coefficient(ufl.Coefficient): def split(obj): - if isinstance(obj, Expression): - return split_expression(obj) return ufl.split_functions.split2(obj) @@ -45,84 +43,3 @@ def TestFunctions(element): def TrialFunctions(element): return split(TrialFunction(element)) - - -class Expression(Coefficient): - def __init__(self, cppcode=None, element=None, cell="triangle", degree=None, is_global=True, on_intersection=False, cellwise_constant=False): - assert cppcode - - if isinstance(cppcode, str): - cppcode = (cppcode,) - - def wrap_return(e): - if "return" not in e: - return "return {};".format(e) - else: - return e - - cppcode = tuple(wrap_return(e) for e in cppcode) - - if cellwise_constant: - if element: - assert element.degree() == 0 - elif degree is not None: - assert degree == 0 - else: - element = FiniteElement("DG", cell, 0) - - if degree is None: - degree = 1 - - if element is None: - if len(cppcode) == 1: - element = FiniteElement("DG", cell, degree) - else: - element = VectorElement("DG", cell, degree, len(cppcode)) - - def element_length(elem): - if isinstance(elem, ufl.FiniteElement): - return 1 - else: - return elem.value_shape()[0] - - assert element_length(element) == len(cppcode) - - self.c_expr = cppcode - self.is_global = is_global - self.on_intersection = on_intersection - self.element = element - - # Initialize a coefficient with a dummy finite element map. - Coefficient.__init__(self, element) - - def __mul__(self, other): - # Allow the combination of Expressions - if isinstance(other, Expression): - def get_sub(elem): - if isinstance(elem, MixedElement) and not isinstance(elem, (VectorElement, TensorElement)): - return elem.sub_elements() - else: - return [elem] - - newelem = MixedElement(*(get_sub(self.element) + get_sub(other.element))) - return Expression(cppcode=self.c_expr + other.c_expr, element=newelem) - else: - return Coefficient.__mul__(self, other) - - # TODO the subdomain_data code needs a uflid, not idea how to get it here - # The standard way through class decorator fails here... - def ufl_id(self): - return 0 - - -def split_expression(expr): - assert isinstance(expr, Expression) - - def element_slice(expression, sub): - offset = sum(subel.value_size() for subel in expression.element.sub_elements()[:sub]) - return expression.c_expr[offset:offset + expr.element.sub_elements()[sub].value_size()] - - return tuple(Expression(cppcode=element_slice(expr, i), - element=expr.element.sub_elements()[i]) - for i in range(expr.element.num_sub_elements()) - ) diff --git a/python/dune/perftool/ufl/visitor.py b/python/dune/perftool/ufl/visitor.py index 5a5c74fa521bae4b08c38badc19d08fdd8bb4ea7..3692d7548bac09df0f998bbd01f0d4684ee0c0a9 100644 --- a/python/dune/perftool/ufl/visitor.py +++ b/python/dune/perftool/ufl/visitor.py @@ -8,7 +8,6 @@ from dune.perftool.ufl.flatoperators import get_operands from dune.perftool.ufl.modified_terminals import (ModifiedTerminalTracker, Restriction, ) -from dune.perftool.ufl.execution import Expression from dune.perftool.tools import maybe_wrap_subscript from loopy import Reduction @@ -122,6 +121,7 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker): # Check if this is a parameter function else: + raise NotImplementedError("Handling non-symbolic parameter functions is currently reevaluated!") # We expect all coefficients to be of type Expression! assert isinstance(o, Expression) diff --git a/test/heatequation/heatequation.ufl b/test/heatequation/heatequation.ufl index de1b61ec6723fd64b87a4550449aef528a306a94..9fe2e20bff8bead3b823bf86c9cb4b4372f155c9 100644 --- a/test/heatequation/heatequation.ufl +++ b/test/heatequation/heatequation.ufl @@ -1,5 +1,10 @@ -f = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());") -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", "triangle", 1) u = TrialFunction(V) diff --git a/test/heatequation/heatequation_dg.ufl b/test/heatequation/heatequation_dg.ufl index d6eb5141d3e7743870ffb3e1dc5a07b456bb00c8..8dfd74302e57c53f3672d4a0e25592dd0c4cd34c 100644 --- a/test/heatequation/heatequation_dg.ufl +++ b/test/heatequation/heatequation_dg.ufl @@ -1,7 +1,11 @@ -f = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());") -g = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return std::exp(-1.*c.two_norm2());", on_intersection=True) - 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("DG", cell, 1) u = TrialFunction(V) diff --git a/test/nonlinear/nonlinear.ufl b/test/nonlinear/nonlinear.ufl index b71ed1da3fa37a4ed68d46766b1dc1066acf8ef4..f81b287e504909b1381a4ac71f3a84af5e34a01c 100644 --- a/test/nonlinear/nonlinear.ufl +++ b/test/nonlinear/nonlinear.ufl @@ -1,7 +1,10 @@ -f = Expression("return -4;") -g = Expression("return x.two_norm2();") +cell = "triangle" +x = SpatialCoordinate(cell) -V = FiniteElement("CG", "triangle", 1) +f = -4. +g = x[0]*x[0] + x[1]*x[1] + +V = FiniteElement("CG", cell, 1) u = TrialFunction(V) v = TestFunction(V) diff --git a/test/poisson/CMakeLists.txt b/test/poisson/CMakeLists.txt index 7dbeb9aaac1ecc389eb9991150efebcb636fbb81..68c635bc79eef0818780ca09d36c808c1eb01954 100644 --- a/test/poisson/CMakeLists.txt +++ b/test/poisson/CMakeLists.txt @@ -43,12 +43,6 @@ dune_add_formcompiler_system_test(UFLFILE poisson_dg.ufl INIFILE poisson_dg_matrix_free.mini ) -# 7. Poisson Test Case: test constant expressions -dune_add_formcompiler_system_test(UFLFILE poisson_cellwise_constant.ufl - BASENAME poisson_cellwise_constant - INIFILE poisson_cellwise_constant.mini - ) - # 8. Poisson with operator counting dune_add_formcompiler_system_test(UFLFILE opcount_poisson_dg.ufl BASENAME opcount_poisson_dg_symdiff diff --git a/test/poisson/opcount_poisson_dg.ufl b/test/poisson/opcount_poisson_dg.ufl index 348ebff253126f87129c0addc5392816da80b308..d731f77c3a9c180a630af3a161e0deaa896bacf5 100644 --- a/test/poisson/opcount_poisson_dg.ufl +++ b/test/poisson/opcount_poisson_dg.ufl @@ -1,9 +1,10 @@ -dim = 2 degree = 1 cell = "quadrilateral" -f = Expression("return -2.0*x.size();", cell=cell) -g = Expression("return x.two_norm2();", on_intersection=True, cell=cell) +x = SpatialCoordinate(cell) + +f = -4. +g = x[0]*x[0] + x[1]*x[1] V = FiniteElement("DG", cell, degree) diff --git a/test/poisson/poisson_cellwise_constant.mini b/test/poisson/poisson_cellwise_constant.mini deleted file mode 100644 index d3730be03345b7844283a56f2ce0aeda8d1f2612..0000000000000000000000000000000000000000 --- a/test/poisson/poisson_cellwise_constant.mini +++ /dev/null @@ -1,14 +0,0 @@ -__name = poisson_cellwise_constant - -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 diff --git a/test/poisson/poisson_cellwise_constant.ufl b/test/poisson/poisson_cellwise_constant.ufl deleted file mode 100644 index 28ee63c28f918c8b5e7be43002b0b3c3c4af9308..0000000000000000000000000000000000000000 --- a/test/poisson/poisson_cellwise_constant.ufl +++ /dev/null @@ -1,12 +0,0 @@ -a = Expression("return 1.0;", cellwise_constant=True) -f = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());") -g = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return std::exp(-1.*c.two_norm2());") - -V = FiniteElement("CG", "triangle", 1) -u = TrialFunction(V) -v = TestFunction(V) - -forms = [(a*inner(grad(u), grad(v)) - f*v)*dx] -exact_solution = g -is_dirichlet = 1 -dirichlet_expression = g diff --git a/test/stokes/stokes_stress.ufl b/test/stokes/stokes_stress.ufl index 2efd1420b25d1f1c1837ce2fb590707a7e2dc3f4..186d614f3a846110b34c53f9fa078dd675b68e37 100644 --- a/test/stokes/stokes_stress.ufl +++ b/test/stokes/stokes_stress.ufl @@ -1,9 +1,8 @@ -v_bctype = Expression("if (x[0] < 1. - 1e-8) return 1; else return 0;", on_intersection=True) -g_v = Expression(("4*x[1]*(1.-x[1])", "0.0")) -g_p = Expression("8*(1.-x[0])") -g_S = Expression(("0.0", "0.0", "-8*x[1] + 4.", "0.0")) - cell = triangle + +x = SpatialCoordinate(cell) +v_bctype = conditional(x[0] < 1. - 1e-8, 1, 0) + P2 = VectorElement("Lagrange", cell) P1 = FiniteElement("Lagrange", cell, 1) P2_stress = TensorElement("DG", cell, 1) @@ -17,5 +16,5 @@ r = (inner(grad(v), S) + inner(grad(u) - S, T) - div(v)*p - q*div(u))*dx forms = [r] is_dirichlet = v_bctype, v_bctype, 0, 0, 0, 0, 0 -dirichlet_expression = g_v, None, None, None, None, None, None -exact_solution = g_v, g_p, g_S \ No newline at end of file +dirichlet_expression = 4*x[1]*(1.-x[1]), 0.0, None, None, None, None, None +exact_solution = 4*x[1]*(1.-x[1]), 0.0, 8*(1.-x[0]), 0.0, 0.0, -1.*8*x[1] + 4., 0.0 \ No newline at end of file diff --git a/test/stokes/stokes_stress_sym.ufl b/test/stokes/stokes_stress_sym.ufl index cff1dee52ffa81881619e7b9580ad76c18538990..f5dc520c07202ee05b4294b153d90f3dcd5b4ab0 100644 --- a/test/stokes/stokes_stress_sym.ufl +++ b/test/stokes/stokes_stress_sym.ufl @@ -1,9 +1,8 @@ -v_bctype = Expression("if (x[0] < 1. - 1e-8) return 1; else return 0;", on_intersection=True) -g_v = Expression(("4*x[1]*(1.-x[1])", "0.0")) -g_p = Expression("8*(1.-x[0])") -g_S = Expression(("0.0", "-8*x[1] + 4.", "0.0")) - cell = triangle + +x = SpatialCoordinate(cell) +v_bctype = conditional(x[0] < 1. - 1e-8, 1, 0) + P2 = VectorElement("Lagrange", cell, 2) P1 = FiniteElement("Lagrange", cell, 1) P2_stress = TensorElement("DG", cell, 1, symmetry={(0, 1): (1, 0)}) @@ -23,5 +22,5 @@ r = (inner(grad(v), S) + inner(2*sym(grad(u)) - S, T) - div(v)*p - q*div(u))*dx forms = [r] is_dirichlet = v_bctype, v_bctype, 0, 0, 0, 0 -dirichlet_expression = g_v, None, None, None, None, None -exact_solution = g_v, g_p, g_S \ No newline at end of file +dirichlet_expression = 4*x[1]*(1.-x[1]), 0.0, None, None, None, None +exact_solution = 4*x[1]*(1.-x[1]), 0.0, 8*(1.-x[0]), 0.0, 0.0, -1.*8*x[1] + 4. \ No newline at end of file diff --git a/test/stokes/stokes_sym.ufl b/test/stokes/stokes_sym.ufl index 0fe8fea45591acd9709d110d43ffb99dfdb20bf7..1ae13db697976a046b461096a530ef315a2d7417 100644 --- a/test/stokes/stokes_sym.ufl +++ b/test/stokes/stokes_sym.ufl @@ -1,7 +1,9 @@ -v_bctype = Expression("if (x[0] < 1. - 1e-8) return 1; else return 0;", on_intersection=True) -g_v = Expression(("4*x[1]*(1.-x[1])", "0.0")) - cell = triangle + +x = SpatialCoordinate(cell) +v_bctype = conditional(x[0] < 1. - 1e-8, 1, 0) +g_v = as_vector((4.*x[1]*(1.-x[1]), 0.0)) + P2 = VectorElement("Lagrange", cell, 2) P1 = FiniteElement("Lagrange", cell, 1) TH = P2 * P1 @@ -17,4 +19,4 @@ r = (inner(2*sym(grad(u)), grad(v)) - div(v)*p - q*div(u))*dx - inner(grad(u).T* forms = [r] is_dirichlet = v_bctype, v_bctype, 0 dirichlet_expression = g_v, None -exact_solution = g_v, 8.*(1.-x[0]) +exact_solution = g_v, 8.*(1.-x[0]) \ No newline at end of file