diff --git a/python/dune/perftool/cgen/clazz.py b/python/dune/perftool/cgen/clazz.py index 8d3c9985da2d10b436742e7c0b941fa0c6d5f6ce..76a86ad0748a1f3d56f6133d0b36f269d32875fb 100644 --- a/python/dune/perftool/cgen/clazz.py +++ b/python/dune/perftool/cgen/clazz.py @@ -39,10 +39,12 @@ class ClassMember(Generable): if isinstance(member, str): from cgen import Line self.member = Line(member) + else: + from collections import Iterable + assert all(isinstance(m, str) for m in self.member) - # We only consider a Generable or a list thereof as member - from collections import Iterable - assert isinstance(self.member, Generable) or (isinstance(self.member, Iterable) and all(isinstance(m, Generable) for m in self.member)) + from cgen import LiteralLines + self.member = LiteralLines("\n" + "\n".join(self.member)) def generate(self): yield "\n\n" diff --git a/python/dune/perftool/compile.py b/python/dune/perftool/compile.py index edd87d4db4d9f1e5364b33219495c7982bf571cb..450e0b1ed1d399eb19206c1f16884a28801a1040 100644 --- a/python/dune/perftool/compile.py +++ b/python/dune/perftool/compile.py @@ -48,7 +48,7 @@ def read_ufl(uflfile): form = transform_form(form, reindexing) # form = transform_form(form, split_arguments) - return form + return form, data.object_names def generate_driver(form, filename): @@ -78,7 +78,7 @@ def generate_driver(form, filename): # This function is the entrypoint of the ufl2pdelab executable def compile_form(): from dune.perftool.options import get_option - form = read_ufl(get_option("uflfile")) + form, namedata = read_ufl(get_option("uflfile")) from dune.perftool.generation import cache_context if get_option("driver_file"): @@ -87,7 +87,7 @@ def compile_form(): if get_option("operator_file"): from dune.perftool.pdelab.localoperator import generate_localoperator_kernels - kernels = generate_localoperator_kernels(form) + kernels = generate_localoperator_kernels(form, namedata) # TODO insert sophisticated analysis/feedback loops here if get_option("interactive"): diff --git a/python/dune/perftool/generation/__init__.py b/python/dune/perftool/generation/__init__.py index a51143e6040c83338910a4352385a88662d31fbb..a65fa989cc055d62aa88a6a24590270c37ab4fed 100644 --- a/python/dune/perftool/generation/__init__.py +++ b/python/dune/perftool/generation/__init__.py @@ -29,4 +29,5 @@ from dune.perftool.generation.loopy import (domain, from dune.perftool.generation.context import (cache_context, generic_context, get_generic_context_value, + namedata_context, ) diff --git a/python/dune/perftool/generation/context.py b/python/dune/perftool/generation/context.py index 13099e01966bce4ac5308c90c638b4f49c8c8366..10abd4f4e7459c53de0cab0ac5e6fac4020f3db6 100644 --- a/python/dune/perftool/generation/context.py +++ b/python/dune/perftool/generation/context.py @@ -53,3 +53,6 @@ def generic_context(key, value): def get_generic_context_value(key): return _generic_context_cache[key] + +import functools +namedata_context=functools.partial(generic_context, "namedata") diff --git a/python/dune/perftool/loopy/transformer.py b/python/dune/perftool/loopy/transformer.py index 3bf90701472525ffc935405977d94445d6808b15..96477b8ec6b8c901bf977f632761fe03d87ca60c 100644 --- a/python/dune/perftool/loopy/transformer.py +++ b/python/dune/perftool/loopy/transformer.py @@ -38,8 +38,22 @@ class UFL2LoopyVisitor(ModifiedTerminalTracker, UFL2PymbolicMapper): def coefficient(self, o): # All trial functions should already be handled assert o.count() != 0 - # TODO implement non-trialfunction coefficients - raise NotImplementedError + + # We expect all coefficients to be of type Expression! + from dune.perftool.ufl.execution import Expression + assert isinstance(o, Expression) + + # Determine the name of the parameter function + from dune.perftool.generation import get_generic_context_value + name = get_generic_context_value("namedata")[id(o)] + + # Trigger the generation of code for this thing in the parameter class + from dune.perftool.pdelab.parameter import parameter_function + parameter_function(name, o) + + # And return a symbol + from pymbolic.primitives import Variable + return Variable(name) def index_sum(self, o): from dune.perftool.ufl.shape import determine_shape diff --git a/python/dune/perftool/pdelab/basis.py b/python/dune/perftool/pdelab/basis.py index 06b0a797b2a418da8ef975574de40ae016930042..f78e8ee56bac75ce442d255d3fdc7318a7fc48a9 100644 --- a/python/dune/perftool/pdelab/basis.py +++ b/python/dune/perftool/pdelab/basis.py @@ -181,7 +181,7 @@ def evaluate_basis(element, name): temporary_variable(name, shape=(name_lfs_bound(element),), shape_impl=('vec',)) cache = name_localbasis_cache(element) lfs = name_lfs(element) - qp = name_quadrature_point() + qp = name_quadrature_position() instruction(inames=(quadrature_iname(), ), code='auto& {} = {}.evaluateFunction({}, {}.finiteElement().localBasis());'.format(name, diff --git a/python/dune/perftool/pdelab/geometry.py b/python/dune/perftool/pdelab/geometry.py index e22c87bdddd5e568db7f5534e9851ec7147e1313..1cf85c4e16f4e2781a6d073a148a75fb90887daa 100644 --- a/python/dune/perftool/pdelab/geometry.py +++ b/python/dune/perftool/pdelab/geometry.py @@ -7,9 +7,23 @@ from dune.perftool.pdelab.quadrature import (name_quadrature_position, ) +@symbol +def name_elementgeometry(): + return 'eg' + + +@symbol +def name_element(): + eg = name_elementgeometry() + return "{}.element()".format(eg) + + @preamble def define_geometry(name): - return "auto {} = eg.geometry();".format(name) + eg = name_elementgeometry() + return "auto {} = {}.geometry();".format(name, + eg + ) @symbol @@ -57,7 +71,8 @@ def define_jacobian_inverse_transposed(name): return quadrature_preamble("{} = {}.jacobianInverseTransposed({});".format(name, geo, pos, - ) + ), + assignees=frozenset({name}), ) diff --git a/python/dune/perftool/pdelab/localoperator.py b/python/dune/perftool/pdelab/localoperator.py index d3c61f83d5261ce69ce2a29004d51d2bc9e353b2..56f8fe232cf34f6ad4b20c39a53fc576bfbc1c9f 100644 --- a/python/dune/perftool/pdelab/localoperator.py +++ b/python/dune/perftool/pdelab/localoperator.py @@ -178,7 +178,10 @@ class AssemblyMethod(ClassMember): def __init__(self, signature, kernel): from loopy import generate_code from cgen import LiteralLines, Block - content = [LiteralLines('\n' + '\n'.join(signature)), Block([LiteralLines('\n' + generate_code(kernel)[0])])] + content = signature + content.append('{') + content.extend(' ' + l for l in generate_code(kernel)[0].split('\n')) + content.append('}') ClassMember.__init__(self, content) @@ -201,7 +204,7 @@ def cgen_class_from_cache(tag, members=[]): return Class(basename, base_classes=base_classes, members=members + pm, constructors=[constructor], tparam_decls=tparams) -def generate_localoperator_kernels(form): +def generate_localoperator_kernels(form, namedata): # For the moment, I do assume that there is but one integral of each type. This might differ # if you use different quadrature orders for different terms. assert len(form.integrals()) == len(set(i.integral_type() for i in form.integrals())) @@ -220,15 +223,20 @@ def generate_localoperator_kernels(form): localoperator_basename() lop_template_ansatz_gfs() lop_template_test_gfs() + from dune.perftool.pdelab.parameter import parameterclass_basename + parameterclass_basename() base_class('Dune::PDELab::LocalOperatorDefaultFlags', classtag="operator") # Have a data structure collect the generated kernels operator_kernels = {} + from dune.perftool.generation import namedata_context + # Generate the necessary residual methods for integral in form.integrals(): - kernel = generate_kernel(integrand=integral.integrand(), measure=integral.integral_type()) + with namedata_context(namedata): + kernel = generate_kernel(integrand=integral.integrand(), measure=integral.integral_type()) operator_kernels[(integral.integral_type(), 'residual')] = kernel # Generate the necessary jacobian methods @@ -241,7 +249,8 @@ def generate_localoperator_kernels(form): jacform = expand_derivatives(derivative(form, form.coefficients()[0])) for integral in jacform.integrals(): - kernel = generate_kernel(integrand=integral.integrand(), measure=integral.integral_type()) + with namedata_context(namedata): + kernel = generate_kernel(integrand=integral.integrand(), measure=integral.integral_type()) operator_kernels[(integral.integral_type(), 'jacobian')] = kernel # TODO: JacobianApply for matrix-free computations. diff --git a/python/dune/perftool/pymbolic/inameset.py b/python/dune/perftool/pymbolic/inameset.py index ba1a9d7ee94f35f3fb10a10c15445a45e1468770..f87c58d614687a297253d5f101a051336dab0fd0 100644 --- a/python/dune/perftool/pymbolic/inameset.py +++ b/python/dune/perftool/pymbolic/inameset.py @@ -9,6 +9,12 @@ class INameMapper(CombineMapper): else: return frozenset([str(i) for i in e.index]) + def map_constant(self, e): + return frozenset({}) + + def map_algebraic_leaf(self, e): + return frozenset({}) + def combine(self, values): return frozenset().union(*values) diff --git a/python/dune/perftool/ufl/execution.py b/python/dune/perftool/ufl/execution.py index 9e91a2680e9f6b16158520ef6bbfe8e4ea977294..a2cd7a143290f0ff461056ef24122e3d0de21ad0 100644 --- a/python/dune/perftool/ufl/execution.py +++ b/python/dune/perftool/ufl/execution.py @@ -32,3 +32,13 @@ def Coefficients(element): def TrialFunctions(element): return split(TrialFunction(element)) + + +class Expression(ufl.Coefficient): + def __init__(self, expr, is_global=True): + assert isinstance(expr, str) + self.c_expr = expr + self.is_global = is_global + + # Initialize a coefficient with a dummy finite element map. + ufl.Coefficient.__init__(self, FiniteElement("DG", "triangle", 0)) diff --git a/python/dune/perftool/ufl/transformations/extract_accumulation_terms.py b/python/dune/perftool/ufl/transformations/extract_accumulation_terms.py index 3dcb8b428781d5cba5dd099d206963a283c912a8..2b8505c56d77e790ad875bc454bd995c1423f2b8 100644 --- a/python/dune/perftool/ufl/transformations/extract_accumulation_terms.py +++ b/python/dune/perftool/ufl/transformations/extract_accumulation_terms.py @@ -16,14 +16,13 @@ import itertools class _ReplacementDict(dict): - def __init__(self, *args): + def __init__(self, good=[], bad=[]): dict.__init__(self) - for a in args: + for a in bad: + self[a] = Zero() + for a in good: self[a] = a - def __getitem__(self, key): - return self.get(key, Zero()) - @ufl_transformation(name="accterms2", extraction_lambda=lambda l: l) def split_into_accumulation_terms(expr): @@ -35,13 +34,23 @@ def split_into_accumulation_terms(expr): if len(filter(lambda ma: ma.argexpr.count() == 1, mod_args)) == 0: for arg in mod_args: # Do the replacement on the expression - accumulation_terms.append(replace_expression(expr, replacemap=_ReplacementDict(arg.expr))) + accumulation_terms.append(replace_expression(expr, + replacemap=_ReplacementDict(good=(arg.expr,), + bad=[ma.expr for ma in mod_args], + ) + ) + ) # and now the case of a rank 2 form: else: for arg1, arg2 in itertools.product(filter(lambda ma: ma.argexpr.count() == 0, mod_args), filter(lambda ma: ma.argexpr.count() == 1, mod_args) ): - accumulation_terms.append(replace_expression(expr, replacemap=_ReplacementDict(arg1.expr, arg2.expr))) + accumulation_terms.append(replace_expression(expr, + replacemap=_ReplacementDict(good=(arg1.expr, arg2.expr), + bad=[ma.expr for ma in mod_args], + ) + ) + ) # and return the result return accumulation_terms diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 7c7839cc3cf6acdfa31f5175f483281be56bcee9..8c530894a2465c7e6de3fa5b352e631c3b4b1dc9 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -5,3 +5,12 @@ add_generated_executable(UFLFILE ../examples/laplace.ufl dune_add_system_test(TARGET laplace INIFILE laplace.mini) + + +add_generated_executable(UFLFILE ../examples/poisson.ufl + TARGET poisson + FORM_COMPILER_ARGS --numerical-jacobian + ) + +dune_add_system_test(TARGET poisson + INIFILE laplace.mini)