Commits

Martin Alnæs  committed 10136cc

Swap ordering of test/trial loops so the trial loop is innermost as this index runs faster in A[].

  • Participants
  • Parent commits 0147ef5

Comments (0)

Files changed (6)

File sandbox/indexmapping/indexmapping.py

-from ufl.common import product
-
-class Range(object):
-    def __init__(self, begin, end):
-        self.begin = begin
-        self.end = end
-
-class IndexMapping(object):
-    def __init__(self, ranges):
-        self.names = sorted(ranges.keys())
-
-        self.ranges = {}
-        self.begin = {}
-        self.end = {}
-        self.size = {}
-
-        for name in self.names:
-            r = ranges[name]
-
-            if isinstance(r, int):
-                begin, end = 0, r
-            elif isinstance(r, str):
-                begin, end = "0", r
-            elif isinstance(r, Range):
-                begin, end = r.begin, r.end
-            else:
-                begin, end = r
-
-            if begin in (0, "0"):
-                size = end
-            elif isinstance(end, str) or isinstance(begin, str):
-                size = "({0} - {1})".format(end, begin)
-            else:
-                size = end - begin
-
-            self.ranges[name] = (begin, end)
-            self.begin[name] = begin
-            self.end[name] = end
-            self.size[name] = size
-
-    def __str__(self):
-        l = ['{!r}: ({!r}, {!r})'.format(name, self.begin[name], self.end[name]) for name in self.names]
-        return '{{ {} }}'.format(', '.join(l))
-
-    def __repr__(self):
-        return 'IndexMapping({!s})'.format(self)
-
-def as_tuple(a):
-    return a if isinstance(a, tuple) else (a,)
-
-def any_str(*args):
-    return any(isinstance(a, str) for a in args)
-
-def dim_mul(i, j):
-    if any_str(i, j):
-        return '{!r} * {!r}'.format(i, j)
-    else:
-        return i*j
-
-def mul_dims(dims):
-    if not dims:
-        return 1
-    ints = [d for d in dims if isinstance(d, int)]
-    strs = [d for d in dims if isinstance(d, str) and d != "1"]
-    if ints and ints != [1]:
-        strs = ['{!r}'.format(product(ints))] + strs
-    if strs:
-        return ' * '.join(strs)
-    return '1'
-
-class AxisMapping(object):
-    def __init__(self, index_mapping, axes):
-        self.index_mapping = index_mapping
-        self.axis_index_names = [as_tuple(a) for a in axes]
-
-        self.num_axes = len(self.axis_index_names)
-        self.num_dims = [len(a) for a in self.axis_index_names]
-
-        # The size of each axis equals the product of the associated index dimensions
-        self.axis_size = [mul_dims([self.index_mapping.size[name] for name in names])
-                          for names in self.axis_index_names]
-
-        # The stride of each axis equals the product of the following axis sizes
-        self.axis_stride = [mul_dims(self.axis_size[i+1:]) for i in range(self.num_axes)]
-
-        # For each axis, the internal strides of each dimension within that axis
-        self.dim_stride = [tuple(mul_dims([self.index_mapping.size[name] for name in self.axis_index_names[i][j+1:]])
-                                 for j in range(self.num_dims[i]))
-                           for i in range(self.num_axes)]
-
-    def format_decl(self):
-        return ''.join('[{}]'.format(size) for size in self.axis_size)
-
-    def format_access(self, **kwargs):
-        return ''.join('[{}]'.format(' + '.join(mul_dims([self.dim_stride[i][j], kwargs[self.axis_index_names[i][j]]])
-                                                for j in range(self.num_dims[i])))
-                       for i in range(self.num_axes))
-                          
-
-# TODO: Change this into unit tests
-def test():
-    ranges = {"i": 4, "j": (10,12), "n": "N", "m": (1, "M") }
-    im = IndexMapping(ranges)
-    print str(eval(repr(im)))
-
-    axes = ["i", "j", "n", "m"]
-    am = AxisMapping(im, axes)
-    print
-    print am.axis_index_names
-    print am.axis_size
-    print am.axis_stride
-    print am.dim_stride
-    print am.format_decl()
-    print am.format_access(i='ii', j='jj', n='nn', m='mm')
-
-    axes = [("i", "j", "n", "m")]
-    am = AxisMapping(im, axes)
-    print
-    print am.axis_index_names
-    print am.axis_size
-    print am.axis_stride
-    print am.dim_stride
-    print am.format_decl()
-    print am.format_access(i='ii', j='jj', n='nn', m='mm')
-
-test()

File site-packages/uflacs/analysis/graph_ssa.py

     Partition 0: Piecewise constant on each cell (including Real and DG0 coefficients)
     Partition 1: Depends on x
     Partition 2: Depends on x and coefficients
-    Partitions [3,3+rank): depend on argument with count 2+rank-partition
+    Partitions [3,3+rank): depend on argument with count partition-3
     """
     rank = len(argument_mapping) # TODO: Take rank as input instead?
+    # TODO: Use named constants for the partition numbers here
 
     modifiers = (Grad, Restricted, Indexed)
     if isinstance(expr, modifiers):
         arg = argument_mapping.get(expr)
         uflacs_assert(arg is not None, "Expecting to have mapping for all arguments.")
         ac = arg.count()
-        assert ac >= 0
-        assert ac < rank
-        p = 2 + rank - ac # Magical formula for mapping argument count to partition
-        assert p >= 3
-        assert p < 3+rank
+        assert 0 <= ac < rank
+        poffset = 3
+        p = poffset + ac
         return p
 
     elif isinstance(expr, Coefficient):

File site-packages/uflacs/backends/ffc/ffc_statement_formatter.py

 from uflacs.utils.log import uflacs_assert, warning, error
 
 from uflacs.codeutils.format_code_structure import ForRange
+from uflacs.codeutils.indexmapping import IndexMapping, AxisMapping
 
 from uflacs.geometry.default_names import names
 from uflacs.geometry.generate_geometry_snippets import (
         self._domain_type = ir["domain_type"]
         self._cell = ir["cell"]
         self._entitytype = ir["entitytype"]
-        self._argument_space_dimensions = ir["prim_idims"]
+        self._argument_space_dimensions = ir["prim_idims"] # FIXME: Dimensions with dS?
 
         self._ir = ir
 
+        # Practical variables
+        rank = len(dependency_handler.mapped_arguments)
+        self._idofs = ["%s%d" % (names.ia, i) for i in range(rank)] # FIXME: Make reusable function for idof
+
         # TODO: Support multiple quadrature loops, affects several places...
         quadrature_rules = ir["quadrature_weights"]
         if quadrature_rules:
         return code
 
     def define_argument_for_loop(self, argument_count):
-        idof = "%s%d" % (names.ia, argument_count) # FIXME: Make reusable function
+        idof = self._idofs[argument_count]
         isize = self._argument_space_dimensions[argument_count]
         return langfmt.for_loop("int", idof, "0", isize)
 
 
         v = dh.mapped_arguments[argument_count]
         assert v.count() == argument_count
-        idof = "%s%d" % (names.ia, argument_count) # FIXME: Make reusable function
+        idof = self._idofs[argument_count]
 
         reqdata = dh.required.get(v)
         if reqdata:
         return code
 
     def output_variable_names(self, num_variables):
-        assert num_variables == 1, "Confused about this design now..."
-        dh = self._dependency_handler
-        rank = len(dh.mapped_arguments)
-
-        idofs = ["%s%d" % (names.ia, i) for i in range(rank)] # FIXME: Make reusable function for idof
-
-        if rank == 0:
-            ii = "0"
-        elif rank == 1:
-            ii = idofs[0]
-        elif rank == 2:
-            # TODO: Validate this. And make the fastest dimension the inner loop in the compiler.
-            ii = langfmt.sum(langfmt.product(idofs[0], self._argument_space_dimensions[1]), idofs[1])
-            #ii = langfmt.sum(langfmt.product(idofs[1], self._argument_space_dimensions[0]), idofs[0])
+        # For integrands, the number of variables is 1,
+        # but that single variable is actually an array entry
+        # indexed with one or more loop variables.
+        uflacs_assert(num_variables == 1, "Only single integrands implemented for ffc.")
+
+        if self._idofs:
+            im = IndexMapping({ idof: idim for idof, idim in zip(self._idofs, self._argument_space_dimensions) })
+            am = AxisMapping(im, [self._idofs])
+            ii = am.format_access()
         else:
-            error # TODO: Generalize to rank r tensors
-        return [langfmt.array_access(names.A, ii)]
+            ii = "0"
+
+        #return [langfmt.array_access(names.A, ii)] # Adds [] TODO: Consolidate fragmented code generation utils
+        return [names.A + ii]

File site-packages/uflacs/codeutils/indexmapping.py

+from ufl.common import product
+
+class Range(object):
+    def __init__(self, begin, end):
+        self.begin = begin
+        self.end = end
+
+class IndexMapping(object):
+    def __init__(self, ranges):
+        self.names = sorted(ranges.keys())
+
+        self.ranges = {}
+        self.begin = {}
+        self.end = {}
+        self.size = {}
+
+        for name in self.names:
+            r = ranges[name]
+
+            if isinstance(r, int):
+                begin, end = 0, r
+            elif isinstance(r, str):
+                begin, end = "0", r
+            elif isinstance(r, Range):
+                begin, end = r.begin, r.end
+            else:
+                begin, end = r
+
+            if begin in (0, "0"):
+                size = end
+            elif isinstance(end, str) or isinstance(begin, str):
+                size = "({0} - {1})".format(end, begin)
+            else:
+                size = end - begin
+
+            self.ranges[name] = (begin, end)
+            self.begin[name] = begin
+            self.end[name] = end
+            self.size[name] = size
+
+    def __str__(self):
+        l = ['{!r}: ({!r}, {!r})'.format(name, self.begin[name], self.end[name]) for name in self.names]
+        return '{{ {} }}'.format(', '.join(l))
+
+    def __repr__(self):
+        return 'IndexMapping({!s})'.format(self)
+
+def as_tuple(a):
+    if isinstance(a, tuple):
+        return a
+    elif isinstance(a, list):
+        return tuple(a)
+    else:
+        return (a,)
+
+def any_str(*args):
+    return any(isinstance(a, str) for a in args)
+
+def dim_mul(i, j):
+    if any_str(i, j):
+        return '{!r} * {!r}'.format(i, j)
+    else:
+        return i*j
+
+def mul_dims(dims):
+    if not dims:
+        return 1
+    ints = [d for d in dims if isinstance(d, int)]
+    strs = [d for d in dims if isinstance(d, str) and d != "1"]
+    if ints and ints != [1]:
+        strs = ['{!r}'.format(product(ints))] + strs
+    if strs:
+        return ' * '.join(strs)
+    return '1'
+
+class AxisMapping(object):
+    def __init__(self, index_mapping, axes):
+        self.index_mapping = index_mapping
+        self.axis_index_names = [as_tuple(a) for a in axes]
+
+        self.num_axes = len(self.axis_index_names)
+        self.num_dims = [len(a) for a in self.axis_index_names]
+
+        # The size of each axis equals the product of the associated index dimensions
+        self.axis_size = [mul_dims([self.index_mapping.size[name] for name in names])
+                          for names in self.axis_index_names]
+
+        # The stride of each axis equals the product of the following axis sizes
+        self.axis_stride = [mul_dims(self.axis_size[i+1:]) for i in range(self.num_axes)]
+
+        # For each axis, the internal strides of each dimension within that axis
+        self.dim_stride = [tuple(mul_dims([self.index_mapping.size[name] for name in self.axis_index_names[i][j+1:]])
+                                 for j in range(self.num_dims[i]))
+                           for i in range(self.num_axes)]
+
+    def format_decl(self):
+        return ''.join('[{}]'.format(size) for size in self.axis_size)
+
+    def format_access(self, **kwargs):
+        return ''.join('[{}]'.format(' + '.join(mul_dims([self.dim_stride[i][j],
+                                                          kwargs.get(self.axis_index_names[i][j],self.axis_index_names[i][j])
+                                                          ])
+                                                for j in range(self.num_dims[i])))
+                       for i in range(self.num_axes))
+
+
+# TODO: Change this into unit tests
+def _test():
+    ranges = {"i": 4, "j": (10,12), "n": "N", "m": (1, "M") }
+    im = IndexMapping(ranges)
+    print str(eval(repr(im)))
+
+    axes = ["i", "j", "n", "m"]
+    am = AxisMapping(im, axes)
+    print
+    print am.axis_index_names
+    print am.axis_size
+    print am.axis_stride
+    print am.dim_stride
+    print am.format_decl()
+    print am.format_access(i='ii', j='jj', n='nn', m='mm')
+
+    axes = [("i", "j", "n", "m")]
+    am = AxisMapping(im, axes)
+    print
+    print am.axis_index_names
+    print am.axis_size
+    print am.axis_stride
+    print am.dim_stride
+    print am.format_decl()
+    print am.format_access(i='ii', j='jj', n='nn', m='mm')
+
+#_test()

File site-packages/uflacs/generation/compiler.py

         if p in partition_codes:
             del partition_codes[p]
 
-    # TODO: This partition numbering is maybe too "magic"?
+    # TODO: This partition numbering is maybe too "magic"? At least make named constants!
     # --- Partition 0: independent of x and arguments.
     p = 0
     if 1:
 
     # --- Partitions 3...3+rank-1: argument function dependency
     rank = len(dh.mapped_arguments)
-    for p in range(3,3+rank):
-        ac = 2+rank-p # Magic formula for mapping partition to argument count
+    poffset = 3
+    for ac in range(rank):
+        p = poffset + ac
         update_loops(tfmt.define_argument_for_loop(ac),
                      tfmt.define_argument_loop_vars(ac),
                      p)
 
     # --- Final partition: final assignments
-    p = 3 + rank
+    p = poffset + rank
     if 1:
         assign_to_variables = tfmt.output_variable_names(len(final_variable_names))
         scaling_factor = tfmt.accumulation_scaling_factor()

File tests/test.sh

 rm -f generated/*
 
 rm -f python_tests.log
+#python runpytests.py 2>&1 | tee -a python_tests.log
 python runpytests.py &> python_tests.log
 if [ $? -ne 0 ]; then
     RESULT=1
 fi
 
 rm -f cpp_tests.log
+#make 2>&1 | tee -a cpp_tests.log
 make &> cpp_tests.log
 if [ $? -ne 0 ]; then
     RESULT=2