Pass the contribution into stage 3 correctly

python/dune/perftool/sumfact/sumfact.py

@@ -22,6 +22,7 @@ from dune.perftool.loopy.buffer import (get_buffer_temporary,
                                         )
 from dune.perftool.sumfact.quadrature import nest_quadrature_loops
 from dune.perftool.pdelab.localoperator import determine_accumulation_space
+from dune.perftool.pdelab.restriction import restricted_name
 from dune.perftool.pdelab.spaces import name_lfs
 from dune.perftool.sumfact.amatrix import (AMatrix,
                                            quadrature_points_per_direction,
@@ -72,6 +73,11 @@ def setup_theta(element, restriction, component, a_matrices):
                        )
 
 
+def name_test_function_contribution(test):
+    grad = "grad_" if test.reference_grad else ""
+    return restricted_name("contrib_{}phi".format(grad), test.restriction)
+
+
 @backend(interface="accum_insn", name="sumfact")
 def generate_accumulation_instruction(visitor, accterm, measure, subdomain_id):
     pymbolic_expr = visitor(accterm.term)
@@ -83,29 +89,32 @@ def generate_accumulation_instruction(visitor, accterm, measure, subdomain_id):
     formdata = get_global_context_value('formdata')
     dim = formdata.geometric_dimension
 
-    # Get a matrix that holds the result of the pymbolic_expr at all quadrature points
-    name = "blubb"
-    temporary_variable(name, shape=(quadrature_points_per_direction(),) * dim, managed=True)
-    insn_dep = instruction(assignee=Subscript(Variable(name), tuple(Variable(i) for i in quadrature_inames())),
-                           expression=pymbolic_expr,
-                           forced_iname_deps=frozenset(quadrature_inames() + visitor.inames),
-                           forced_iname_deps_is_final=True,
-                           )
-
+    # Get the a matrices needed for this accumulation term
     theta_transposed = name_theta_transposed()
     rows = basis_functions_per_direction()
     cols = quadrature_points_per_direction()
     a_matrix = AMatrix(theta_transposed, rows, cols)
     a_matrices = (a_matrix, a_matrix)
 
-    # Do stage 1
-    initialize_buffer("reffub",
-                      base_storage_size=product(max(mat.rows, mat.cols) for mat in a_matrices),
-                      num=2
-                      )
+    # Get a matrix that holds the result of the pymbolic_expr at all quadrature points
+    buffer = name_test_function_contribution(accterm.argument)
+    temp = initialize_buffer(buffer,
+                             base_storage_size=product(max(mat.rows, mat.cols) for mat in a_matrices),
+                             num=2
+                             ).get_temporary(shape=(quadrature_points_per_direction(),) * dim)
+
+    # Issue an instruction in the quadrature loop that fills the buffer
+    # with the evaluation of the contribution at all quadrature points
+    insn_dep = instruction(assignee=Subscript(Variable(temp), tuple(Variable(i) for i in quadrature_inames())),
+                           expression=pymbolic_expr,
+                           forced_iname_deps=frozenset(quadrature_inames() + visitor.inames),
+                           forced_iname_deps_is_final=True,
+                           )
 
+    # Add a sum factorization kernel that implements the multiplication
+    # with the test function (step 3)
     result, insn_dep = sum_factorization_kernel(a_matrices,
-                                                "reffub",
+                                                buffer,
                                                 insn_dep=frozenset({insn_dep}),
                                                 additional_inames=frozenset(visitor.inames),
                                                 )