From 56c9b94df9ad0bcd3a6071097e710f18aef1ead1 Mon Sep 17 00:00:00 2001 From: Reuben Hill Date: Fri, 9 Oct 2020 15:50:26 +0100 Subject: [PATCH 1/3] Remove coords argument in compile_expression_dual_evaluation First step towards dealing with #232. Use fake coordinates coefficient when required instead of data carrying .coordinates firedrake Function. --- tsfc/driver.py | 24 +++++++++++++----------- 1 file changed, 13 insertions(+), 11 deletions(-) diff --git a/tsfc/driver.py b/tsfc/driver.py index b8160a85..95ff9bd6 100644 --- a/tsfc/driver.py +++ b/tsfc/driver.py @@ -267,7 +267,7 @@ def name_multiindex(multiindex, name): return builder.construct_kernel(kernel_name, impero_c, index_names, quad_rule) -def compile_expression_dual_evaluation(expression, to_element, coordinates, *, +def compile_expression_dual_evaluation(expression, to_element, *, domain=None, interface=None, parameters=None, coffee=False): """Compile a UFL expression to be evaluated against a compile-time known reference element's dual basis. @@ -276,8 +276,7 @@ def compile_expression_dual_evaluation(expression, to_element, coordinates, *, :arg expression: UFL expression :arg to_element: A FInAT element for the target space - :arg coordinates: the coordinate function - :arg domain: optional UFL domain the expression is defined on (useful when expression contains no domain). + :arg domain: optional UFL domain the expression is defined on (required when expression contains no domain). :arg interface: backend module for the kernel interface :arg parameters: parameters object :arg coffee: compile coffee kernel instead of loopy kernel @@ -328,17 +327,20 @@ def compile_expression_dual_evaluation(expression, to_element, coordinates, *, argument_multiindices = tuple(builder.create_element(arg.ufl_element()).get_indices() for arg in arguments) - # Replace coordinates (if any) - domain = expression.ufl_domain() - if domain: - assert coordinates.ufl_domain() == domain - builder.domain_coordinate[domain] = coordinates - builder.set_cell_sizes(domain) + # Replace coordinates (if any) unless otherwise specified by kwarg + if domain is None: + domain = expression.ufl_domain() + assert domain is not None + + # Create a fake coordinate coefficient for a domain. + coords_coefficient = ufl.Coefficient(ufl.FunctionSpace(domain, domain.ufl_coordinate_element())) + builder.domain_coordinate[domain] = coords_coefficient + builder.set_cell_sizes(domain) # Collect required coefficients coefficients = extract_coefficients(expression) if has_type(expression, GeometricQuantity) or any(fem.needs_coordinate_mapping(c.ufl_element()) for c in coefficients): - coefficients = [coordinates] + coefficients + coefficients = [coords_coefficient] + coefficients builder.set_coefficients(coefficients) # Split mixed coefficients @@ -346,7 +348,7 @@ def compile_expression_dual_evaluation(expression, to_element, coordinates, *, # Translate to GEM kernel_cfg = dict(interface=builder, - ufl_cell=coordinates.ufl_domain().ufl_cell(), + ufl_cell=domain.ufl_cell(), argument_multiindices=argument_multiindices, index_cache={}, scalar_type=parameters["scalar_type"]) From 649a4519482ff584dae8409d879d06966d5c9add Mon Sep 17 00:00:00 2001 From: Reuben Nixon-Hill Date: Mon, 19 Oct 2020 17:42:07 +0100 Subject: [PATCH 2/3] test compile_expression_dual_evaluation with symbolics Tests for non-dependence of compile_expression_dual_evaluation on numerics by checking it works with purely symbolic expressions. --- tests/test_dual_evaluation.py | 57 +++++++++++++++++++++++++++++++++++ 1 file changed, 57 insertions(+) create mode 100644 tests/test_dual_evaluation.py diff --git a/tests/test_dual_evaluation.py b/tests/test_dual_evaluation.py new file mode 100644 index 00000000..da6bf7e6 --- /dev/null +++ b/tests/test_dual_evaluation.py @@ -0,0 +1,57 @@ +import pytest +import ufl +from tsfc.finatinterface import create_element +from tsfc import compile_expression_dual_evaluation + + +def test_ufl_only_simple(): + mesh = ufl.Mesh(ufl.VectorElement("P", ufl.triangle, 1)) + V = ufl.FunctionSpace(mesh, ufl.FiniteElement("P", ufl.triangle, 2)) + v = ufl.Coefficient(V) + expr = ufl.inner(v, v) + W = V + to_element = create_element(W.ufl_element()) + ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + + +def test_ufl_only_spatialcoordinate(): + mesh = ufl.Mesh(ufl.VectorElement("P", ufl.triangle, 1)) + V = ufl.FunctionSpace(mesh, ufl.FiniteElement("P", ufl.triangle, 2)) + x, y = ufl.SpatialCoordinate(mesh) + expr = x*y - y**2 + x + W = V + to_element = create_element(W.ufl_element()) + ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + + +def test_ufl_only_from_contravariant_piola(): + mesh = ufl.Mesh(ufl.VectorElement("P", ufl.triangle, 1)) + V = ufl.FunctionSpace(mesh, ufl.FiniteElement("RT", ufl.triangle, 1)) + v = ufl.Coefficient(V) + expr = ufl.inner(v, v) + W = ufl.FunctionSpace(mesh, ufl.FiniteElement("P", ufl.triangle, 2)) + to_element = create_element(W.ufl_element()) + ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + + +def test_ufl_only_to_contravariant_piola(): + mesh = ufl.Mesh(ufl.VectorElement("P", ufl.triangle, 1)) + V = ufl.FunctionSpace(mesh, ufl.FiniteElement("P", ufl.triangle, 2)) + v = ufl.Coefficient(V) + expr = ufl.as_vector([v, v]) + W = ufl.FunctionSpace(mesh, ufl.FiniteElement("RT", ufl.triangle, 1)) + to_element = create_element(W.ufl_element()) + ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + + +def test_ufl_only_shape_mismatch(): + mesh = ufl.Mesh(ufl.VectorElement("P", ufl.triangle, 1)) + V = ufl.FunctionSpace(mesh, ufl.FiniteElement("RT", ufl.triangle, 1)) + v = ufl.Coefficient(V) + expr = ufl.inner(v, v) + assert expr.ufl_shape == () + W = V + to_element = create_element(W.ufl_element()) + assert to_element.value_shape == (2,) + with pytest.raises(ValueError): + ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) From ed89971538271dd5bdc8d2e02189a3155ca82cbe Mon Sep 17 00:00:00 2001 From: Reuben Nixon-Hill Date: Wed, 18 Nov 2020 11:50:09 +0000 Subject: [PATCH 3/3] Make clear if using manufactured coordinates coefficient So that caller can replace if necessary --- tests/test_dual_evaluation.py | 14 +++++++++----- tsfc/driver.py | 13 +++++++------ tsfc/kernel_interface/firedrake_loopy.py | 9 ++++++--- 3 files changed, 22 insertions(+), 14 deletions(-) diff --git a/tests/test_dual_evaluation.py b/tests/test_dual_evaluation.py index da6bf7e6..8cabba70 100644 --- a/tests/test_dual_evaluation.py +++ b/tests/test_dual_evaluation.py @@ -11,7 +11,8 @@ def test_ufl_only_simple(): expr = ufl.inner(v, v) W = V to_element = create_element(W.ufl_element()) - ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + ast, oriented, needs_cell_sizes, coefficients, first_coeff_fake_coords, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + assert first_coeff_fake_coords is False def test_ufl_only_spatialcoordinate(): @@ -21,7 +22,8 @@ def test_ufl_only_spatialcoordinate(): expr = x*y - y**2 + x W = V to_element = create_element(W.ufl_element()) - ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + ast, oriented, needs_cell_sizes, coefficients, first_coeff_fake_coords, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + assert first_coeff_fake_coords is True def test_ufl_only_from_contravariant_piola(): @@ -31,7 +33,8 @@ def test_ufl_only_from_contravariant_piola(): expr = ufl.inner(v, v) W = ufl.FunctionSpace(mesh, ufl.FiniteElement("P", ufl.triangle, 2)) to_element = create_element(W.ufl_element()) - ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + ast, oriented, needs_cell_sizes, coefficients, first_coeff_fake_coords, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + assert first_coeff_fake_coords is True def test_ufl_only_to_contravariant_piola(): @@ -41,7 +44,8 @@ def test_ufl_only_to_contravariant_piola(): expr = ufl.as_vector([v, v]) W = ufl.FunctionSpace(mesh, ufl.FiniteElement("RT", ufl.triangle, 1)) to_element = create_element(W.ufl_element()) - ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + ast, oriented, needs_cell_sizes, coefficients, first_coeff_fake_coords, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + assert first_coeff_fake_coords is True def test_ufl_only_shape_mismatch(): @@ -54,4 +58,4 @@ def test_ufl_only_shape_mismatch(): to_element = create_element(W.ufl_element()) assert to_element.value_shape == (2,) with pytest.raises(ValueError): - ast, oriented, needs_cell_sizes, coefficients, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) + ast, oriented, needs_cell_sizes, coefficients, first_coeff_fake_coords, _ = compile_expression_dual_evaluation(expr, to_element, coffee=False) diff --git a/tsfc/driver.py b/tsfc/driver.py index 95ff9bd6..d6db7120 100644 --- a/tsfc/driver.py +++ b/tsfc/driver.py @@ -332,15 +332,16 @@ def compile_expression_dual_evaluation(expression, to_element, *, domain = expression.ufl_domain() assert domain is not None - # Create a fake coordinate coefficient for a domain. - coords_coefficient = ufl.Coefficient(ufl.FunctionSpace(domain, domain.ufl_coordinate_element())) - builder.domain_coordinate[domain] = coords_coefficient - builder.set_cell_sizes(domain) - # Collect required coefficients + first_coefficient_fake_coords = False coefficients = extract_coefficients(expression) if has_type(expression, GeometricQuantity) or any(fem.needs_coordinate_mapping(c.ufl_element()) for c in coefficients): + # Create a fake coordinate coefficient for a domain. + coords_coefficient = ufl.Coefficient(ufl.FunctionSpace(domain, domain.ufl_coordinate_element())) + builder.domain_coordinate[domain] = coords_coefficient + builder.set_cell_sizes(domain) coefficients = [coords_coefficient] + coefficients + first_coefficient_fake_coords = True builder.set_coefficients(coefficients) # Split mixed coefficients @@ -433,7 +434,7 @@ def compile_expression_dual_evaluation(expression, to_element, *, # Handle kernel interface requirements builder.register_requirements([ir]) # Build kernel tuple - return builder.construct_kernel(return_arg, impero_c, index_names) + return builder.construct_kernel(return_arg, impero_c, index_names, first_coefficient_fake_coords) def lower_integral_type(fiat_cell, integral_type): diff --git a/tsfc/kernel_interface/firedrake_loopy.py b/tsfc/kernel_interface/firedrake_loopy.py index ee5cafa4..eb90cf72 100644 --- a/tsfc/kernel_interface/firedrake_loopy.py +++ b/tsfc/kernel_interface/firedrake_loopy.py @@ -17,7 +17,7 @@ # Expression kernel description type -ExpressionKernel = namedtuple('ExpressionKernel', ['ast', 'oriented', 'needs_cell_sizes', 'coefficients', 'tabulations']) +ExpressionKernel = namedtuple('ExpressionKernel', ['ast', 'oriented', 'needs_cell_sizes', 'coefficients', 'first_coefficient_fake_coords', 'tabulations']) def make_builder(*args, **kwargs): @@ -153,12 +153,14 @@ def register_requirements(self, ir): provided by the kernel interface.""" self.oriented, self.cell_sizes, self.tabulations = check_requirements(ir) - def construct_kernel(self, return_arg, impero_c, index_names): + def construct_kernel(self, return_arg, impero_c, index_names, first_coefficient_fake_coords): """Constructs an :class:`ExpressionKernel`. :arg return_arg: loopy.GlobalArg for the return value :arg impero_c: gem.ImperoC object that represents the kernel :arg index_names: pre-assigned index names + :arg first_coefficient_fake_coords: If true, the kernel's first + coefficient is a constructed UFL coordinate field :returns: :class:`ExpressionKernel` object """ args = [return_arg] @@ -173,7 +175,8 @@ def construct_kernel(self, return_arg, impero_c, index_names): loopy_kernel = generate_loopy(impero_c, args, self.scalar_type, "expression_kernel", index_names) return ExpressionKernel(loopy_kernel, self.oriented, self.cell_sizes, - self.coefficients, self.tabulations) + self.coefficients, first_coefficient_fake_coords, + self.tabulations) class KernelBuilder(KernelBuilderBase):