Adds some comments

dune/perftool/localbasiscache.hh

@@ -31,7 +31,7 @@ class CacheReturnProxy
   const std::vector<T>* v;
 };
 
-//TODO FS hier entfernen und nur LocalBasisType benutzen
+
 template<typename LocalBasisType >
 class LocalBasisCacheWithoutReferences
 {

python/dune/perftool/blockstructured/__init__.py

@@ -24,7 +24,7 @@ class BlockStructuredInterface(PDELabInterface):
     # Local function space related generator functions
     #
 
-    # TODO current way to squeeze subelem iname in, not really understandable
+    # TODO current way to squeeze subelem iname in, not really ideal
     def lfs_inames(self, element, restriction, number=None, context=''):
         return lfs_inames(element, restriction, number, context) + sub_element_inames()
 

python/dune/perftool/blockstructured/argument.py

@@ -19,4 +19,5 @@ def pymbolic_coefficient(container, lfs, element, index):
     if not isinstance(lfs, prim.Expression):
         lfs = prim.Variable(lfs)
 
+    # use higher order FEM index instead of Q1 index
     return prim.Call(CoefficientAccess(container), (lfs, micro_index_to_macro_index(element, index),))

python/dune/perftool/blockstructured/basis.py

@@ -20,6 +20,7 @@ from dune.perftool.pdelab.quadrature import pymbolic_quadrature_position_in_cell
 from dune.perftool.pdelab.spaces import type_gfs
 
 
+# define FE basis explicitly in localoperator
 @backend(interface="typedef_localbasis", name="blockstructured")
 @class_member(classtag="operator")
 def typedef_localbasis(element, name):

python/dune/perftool/blockstructured/geometry.py

@@ -12,25 +12,31 @@ from dune.perftool.blockstructured.tools import sub_element_inames
 import pymbolic.primitives as prim
 
 
+# scale determinant according to the order of the blockstructure
 def pymbolic_jacobian_determinant():
     return prim.Quotient(prim.Variable(name_jacobian_determinant()),
                          prim.Power(get_option("number_of_blocks"), local_dimension()))
 
 
+# scale Jacobian according to the order of the blockstructure
 def pymbolic_jacobian_inverse_transposed(i, j, restriction):
     return prim.Product((get_option("number_of_blocks"),
                          prim.Subscript(prim.Variable(name_jacobian_inverse_transposed(restriction)), (j, i))))
 
 
+# scale determinant according to the order of the blockstructure
 def pymbolic_facet_jacobian_determinant():
     return prim.Quotient(prim.Variable(name_facet_jacobian_determinant()),
                          prim.Power(get_option("number_of_blocks"), local_dimension()))
 
 
+# translate a point in the micro element into macro coordinates
 def define_point_in_macro(name, point_in_micro):
     dim = local_dimension()
     temporary_variable(name, shape_impl=('fv',), shape=(dim,))
 
+    # point_macro = (point_micro + index_micro) / number_of_blocks
+    # where index_micro = tensor index of the micro element
     subelem_inames = sub_element_inames()
     for i in range(dim):
         if isinstance(point_in_micro, prim.Subscript):

python/dune/perftool/blockstructured/tools.py

@@ -15,6 +15,8 @@ from dune.perftool.options import get_option
 import pymbolic.primitives as prim
 
 
+# add inames over the micro elements in tensor representation,
+# i.e. each element has (i_1,i_2,...,i_d) indices
 @iname
 def sub_element_inames():
     name = "subel"
@@ -27,20 +29,30 @@ def sub_element_inames():
     return inames
 
 
+# define inames for boundary integration
+# In the case of boundary integrationsince we have only d-1 loops,
+# but we need always d inames to compute the macro index.
+# Therefore we must find out which iname is constant.
+# Since loo.py cannot handle if-else blocks very well, this whole
+# computation seems a bit cumbersome.
 def sub_facet_inames():
     subelem_inames = sub_element_inames()
 
     center = pymbolic_in_cell_coordinates(prim.Variable(name_localcenter()), Restriction.NEGATIVE)
 
+    # check if iname[index] must be constant or not
     def predicate(index):
         return prim.Comparison(prim.Subscript(center, index), "==", 0.5)
 
     def conditional_instruction(index):
+        # instruction for not constant iname
+        # special case for the third iname
         instruction(assignee=prim.Variable(inames[index]),
                     expression=prim.Variable(subelem_inames[1 if index == 2 else 0]),
                     within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
                     predicates=frozenset([predicate(index)])
                     )
+        # instruction for constant iname
         instruction(assignee=prim.Variable(inames[index]),
                     expression=prim.Product(((k - 1), prim.Subscript(center, (index,)))),
                     within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
@@ -59,6 +71,7 @@ def sub_facet_inames():
     else:
         inames = inames + ("y",)
         temporary_variable(inames[1])
+        # one additional condition is needed in 3d for the second iname
         instruction(assignee=prim.Variable(inames[1]),
                     expression=prim.Variable(subelem_inames[0]),
                     within_inames=frozenset(subelem_inames).union(frozenset(quadrature_inames())),
@@ -82,18 +95,22 @@ def sub_facet_inames():
     return inames
 
 
+# compute sequential index for given tensor index, the same as index in base-k to base-10
 def tensor_index_to_sequential_index(indices, k):
     return prim.Sum(tuple(index * k ** i for i, index in enumerate(indices)))
 
 
+# compute tensor index for given sequential index, only for Q1
 def to_tensor_index(iname, k):
     return tuple(prim.RightShift(prim.BitwiseAnd((prim.Variable(iname), 2**n)), n) for n in range(world_dimension()))
 
 
+# compute tensor index for given sequential index, the same as index in base-10 to base-k
 def to_tensor_index_v2(iname, k):
     return tuple(prim.Remainder(prim.Variable(iname) / (k**i), k) for i in range(world_dimension()))
 
 
+# compute index for higher order FEM for a given Q1 index
 def micro_index_to_macro_index(element, iname):
     it = get_global_context_value("integral_type")
     if it == "cell":