First working example: spatial coordinate on volumes

python/dune/perftool/sumfact/__init__.py

@@ -54,7 +54,11 @@ class SumFactInterface(PDELabInterface):
 
     def pymbolic_quadrature_weight(self):
         return quadrature_weight()
+#
+#    def pymbolic_spatial_coordinate(self):
+#        from dune.perftool.pdelab.geometry import to_global
+#        return to_global(pymbolic_quadrature_position())
 
     def pymbolic_spatial_coordinate(self):
-        from dune.perftool.pdelab.geometry import to_global
-        return to_global(pymbolic_quadrature_position())
+        from dune.perftool.sumfact.geometry import pymbolic_spatial_coordinate
+        return pymbolic_spatial_coordinate()
\ No newline at end of file

python/dune/perftool/sumfact/geometry.py

+""" Sum factorized geometry evaluations """
+
+from dune.perftool.generation import (domain,
+                                      get_backend,
+                                      get_counted_variable,
+                                      iname,
+                                      instruction,
+                                      kernel_cached,
+                                      temporary_variable,
+                                      )
+from dune.perftool.loopy.buffer import get_buffer_temporary
+from dune.perftool.pdelab.geometry import world_dimension
+from dune.perftool.sumfact.symbolic import SumfactKernelInputBase
+from dune.perftool.sumfact.vectorization import attach_vectorization_info
+
+from pytools import ImmutableRecord
+
+import pymbolic.primitives as prim
+
+
+@iname
+def corner_iname():
+    name = get_counted_variable("corneriname")
+    domain(name, 2 ** world_dimension())
+    return name
+
+
+class GeoCornersInput(SumfactKernelInputBase, ImmutableRecord):
+    def __init__(self, dir):
+        ImmutableRecord.__init__(self, dir=dir)
+
+    def realize(self, sf, index, insn_dep):
+        name = get_buffer_temporary(sf.buffer,
+                                    shape=(2 ** world_dimension(), sf.vector_width),
+                                    name="input_{}".format(sf.buffer)
+                                    )
+
+        ciname = corner_iname()
+
+        from dune.perftool.pdelab.geometry import name_geometry
+        geo = name_geometry()
+
+        # NB: We need to realize this as a C instruction, because the corner
+        #     method does return a non-scalar, which does not fit into the current
+        #     loopy philosophy for function calls. This problem will be solved once
+        #     #11 is resolved.
+        code = "{}[{}] = {}.corner({})[{}];".format(name,
+                                                    ciname,
+                                                    geo,
+                                                    ciname,
+                                                    self.dir,
+                                                    )
+
+        instruction(code=code,
+                    within_inames=frozenset({ciname}),
+                    assignees=(name,),
+                    tags=frozenset({"sumfact_stage{}".format(sf.stage)}),
+                    )
+
+
+@kernel_cached
+def pymbolic_spatial_coordinate():
+    # Construct the matrix sequence for this sum factorization
+    from dune.perftool.sumfact.tabulation import construct_basis_matrix_sequence
+    matrix_sequence = construct_basis_matrix_sequence(basis_size=2)
+
+    expressions = []
+    insn_dep = frozenset()
+    for i in range(world_dimension()):
+        inp = GeoCornersInput(i)
+
+        from dune.perftool.sumfact.symbolic import SumfactKernel
+        sf = SumfactKernel(matrix_sequence=matrix_sequence,
+                           input=inp,
+                           )
+
+        vsf = attach_vectorization_info(sf)
+
+        # Add a sum factorization kernel that implements the evaluation of
+        # the basis functions at quadrature points (stage 1)
+        from dune.perftool.sumfact.realization import realize_sum_factorization_kernel
+        var, insn_dep = realize_sum_factorization_kernel(vsf.copy(insn_dep=vsf.insn_dep.union(insn_dep)))
+
+        expressions.append(prim.Subscript(var, vsf.quadrature_index(sf)))
+
+    # Return an indexable temporary with the results!
+    name = "pos_global"
+    temporary_variable(name, shape=(world_dimension(),))
+    for i, expr in enumerate(expressions):
+        assignee = prim.Subscript(prim.Variable(name), (i,))
+        instruction(assignee=assignee,
+                    expression=expr,
+                    within_inames=frozenset(get_backend("quad_inames")()),
+                    within_inames_is_final=True,
+                    depends_on=insn_dep,
+                    )
+
+    return prim.Variable(name)

python/dune/perftool/sumfact/symbolic.py

@@ -118,7 +118,6 @@ class SumfactKernel(SumfactKernelBase, ImmutableRecord, prim.Variable):
 
         if stage == 1:
             assert isinstance(input, SumfactKernelInputBase)
-            restriction = input.restriction
 
         if stage == 3:
             assert isinstance(restriction, tuple)